University of Southern California Dissertations and Theses

Control of spacecraft with flexible structures using pulse-modulated thrusters

(USC Thesis Other)

Control of spacecraft with flexible structures using pulse-modulated thrusters

PDF

Download a page range

Download transcript

Copy asset link

Request this asset

Transcript (if available)

Content CONTROL OF SPACECRAFT WITH FLEXIBLE STRUCTURES USING PULSE-MODULATED THRUSTERS by Shalini Reddy A Dissertation Presented to the FACULTY OF THE USC GRADUATE SCHOOL UNIVERSITY OF SOUTHERN CALIFORNIA In Partial Fulllment of the Requirements for the Degree DOCTOR OF PHILOSOPHY (MECHANICAL ENGINEERING) August 2012 Copyright 2012 Shalini Reddy Dedication To my niece ii Acknowledgments I would like to thank my research advisor, Prof. Henryk Flashner of the University of Southern California's Department of Aerospace and Mechanical Engineering. His contin- ued patience, guidance, and understanding during the course of our research is immensely appreciated, and I am sincerely grateful to him for helping through this chapter in my life. I also want to thank the members of my Dissertation Committee for their guidance and support: Prof. Larry Redekopp, Prof. Bingen Yang, and Prof. Georey Shi ett of the Aerospace and Mechanical Engineering Department and Prof. Igor Kukavica of the Mathematics Department. I am especially thankful for the support and encouragement from my family and friends during my years as a graduate student. iii Table of Contents Dedication ii Acknowledgments iii List of Tables vi List of Figures vii Abstract x Chapter 1: Introduction 1 1.1 Review of Previous Work . . . . . . . . . . . . . . . . . . . . . . . . . . . 3 1.2 Dissertation Outline . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 6 Chapter 2: Mathematical Preliminaries 8 2.1 Matrix Theory . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 8 2.2 Cross Product . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 10 2.3 Quaternions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 11 2.4 Dynamics . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 12 2.4.1 Lagrangian Mechanics . . . . . . . . . . . . . . . . . . . . . . . . . 12 2.4.2 Boltzmann-Hamel Equation . . . . . . . . . . . . . . . . . . . . . . 13 2.5 Stability Concepts . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 16 2.5.1 Passivity . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 16 2.5.2 Lyapunov Stability . . . . . . . . . . . . . . . . . . . . . . . . . . . 22 2.5.3 Boundedness . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 24 Chapter 3: Modeling of Closed-Loop Systems in a Three Dimensional Space 27 3.1 Modeling of Spacecraft with Flexible Structures . . . . . . . . . . . . . . . 27 3.1.1 Modeling of Flexible Structures using Modal Analysis . . . . . . . 30 3.1.2 Position and Velocity . . . . . . . . . . . . . . . . . . . . . . . . . 33 3.1.3 Kinetic Energy of the Flexible Structure . . . . . . . . . . . . . . . 36 3.1.4 Elastic Potential Energy of the Flexible Structure . . . . . . . . . 38 3.1.5 Kinetic Energy of the Rigid Central Body . . . . . . . . . . . . . . 39 3.1.6 Gravitational Potential Energy of the Spacecraft . . . . . . . . . . 40 3.1.7 Energy of the Spacecraft . . . . . . . . . . . . . . . . . . . . . . . . 41 3.1.8 Derivation of the Equations of Motion . . . . . . . . . . . . . . . . 42 iv 3.2 Implementation of Commanded Control . . . . . . . . . . . . . . . . . . . 50 3.2.1 Attitude and Guidance Control Law . . . . . . . . . . . . . . . . . 50 3.2.2 Thruster Selection and Jet Firing Logic . . . . . . . . . . . . . . . 51 3.3 Thruster Dynamics . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 56 3.4 Implementation of Vibration Suppression . . . . . . . . . . . . . . . . . . 56 Chapter 4: Trajectory Tracking Control in a Three Dimensional Space 61 4.0.1 Controller States and Terms . . . . . . . . . . . . . . . . . . . . . . 61 4.0.2 Stability of Continuous Trajectory Tracking . . . . . . . . . . . . . 66 Chapter 5: Simulation 77 5.1 Description of Simulation . . . . . . . . . . . . . . . . . . . . . . . . . . . 77 5.2 Spacecraft and Control Parameters . . . . . . . . . . . . . . . . . . . . . . 82 5.2.1 Rigid Central Body . . . . . . . . . . . . . . . . . . . . . . . . . . 82 5.2.2 Flexible Structures . . . . . . . . . . . . . . . . . . . . . . . . . . . 82 5.2.3 Control Parameters . . . . . . . . . . . . . . . . . . . . . . . . . . 83 5.3 Simulation Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 84 5.3.1 Vertical Trajectory . . . . . . . . . . . . . . . . . . . . . . . . . . . 84 5.3.2 Vertical Trajectory with Angular Acceleration . . . . . . . . . . . . 96 5.3.3 Diagonal Trajectory with Commanded Angular Acceleration about all Axes . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 108 5.3.4 Vertical Trajectory with Resting End Conditions . . . . . . . . . . 125 Chapter 6: Concluding Remarks and Future Work 135 References 138 Appendix A: Quaternion 146 A.1 Quaternion Denition . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 146 A.2 Quaternion Algebra . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 147 A.3 Unit Quaternion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 149 A.4 Quaternion Rotation Operator . . . . . . . . . . . . . . . . . . . . . . . . 150 A.5 Identities with Arbitrary Vectors . . . . . . . . . . . . . . . . . . . . . . . 154 A.6 Quaternion Time Derivative . . . . . . . . . . . . . . . . . . . . . . . . . . 154 Appendix B: Kronecker Product and Matrix Derivative 156 B.1 Kronecker Product . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 156 B.2 Matrix Derivative . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 157 Appendix C: Modeling of Closed-Loop Systems in a Two Dimensional Space 161 C.1 Derivation of Equations of Motion using Lagrangian Mechanics . . . . . . 163 C.2 Implementation of Commanded Control . . . . . . . . . . . . . . . . . . . 174 Appendix D: Trajectory Tracking in a Two Dimensional Space 179 v List of Tables 5.1 Parameters of Vertical Trajectory . . . . . . . . . . . . . . . . . . . . . . . 85 5.2 Parameters of Vertical Trajectory with Constant Angular Acceleration . . 96 5.3 Parameters of Diagonal Trajectory with Constant Angular Accelerations about all Axes . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 108 5.4 Parameters of Vertical Trajectory with End Conditions . . . . . . . . . . . 125 vi List of Figures 2.1 Feedback Interconnection of Nonlinear Systems . . . . . . . . . . . . . . . 19 3.1 Model of Spacecraft in 3-Dimensional Reference Frame . . . . . . . . . . . 29 3.2 Cantilever Beam . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 31 3.3 Firing Logic for Thruster . . . . . . . . . . . . . . . . . . . . . . . . . . . 52 3.4 Attitude and Guidance Control System of Spacecraft . . . . . . . . . . . . 53 3.5 Mars Lander Cruise Conguration [49] . . . . . . . . . . . . . . . . . . . . 54 3.6 Gas Jet/ Thruster Firing Prole . . . . . . . . . . . . . . . . . . . . . . . 56 5.1 Matlab Simulink Model used for Simulation . . . . . . . . . . . . . . . . . 81 5.2 Commanded and Actual Displacement in Z-Direction . . . . . . . . . . . . 86 5.3 Commanded and Actual Velocity in Z-Direction Without Error Compen- sation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 87 5.4 Commanded and Actual Velocity in Z-Direction . . . . . . . . . . . . . . . 88 5.5 Commanded and Actual Displacement in Z-Direction Without Error Com- pensation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 89 5.6 Commanded and Actual Rotation Angle of Spacecraft . . . . . . . . . . . 90 5.7 Commanded and Actual Rotation Angle of Spacecraft Without Error Com- pensation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 91 5.8 End De ection of Flexible Structure Position on Central Body's x-axis . 92 vii 5.9 End De ection of Flexible Structure Position on Central Body's x-axis Without Error Compensation . . . . . . . . . . . . . . . . . . . . . . . . . 93 5.10 End De ection of Flexible Structure Position on Central Body's x-axis Without Vibration Suppression . . . . . . . . . . . . . . . . . . . . . . . . 94 5.11 Commanded and Applied Forces in Z-Directions . . . . . . . . . . . . . . 95 5.12 Commanded and Actual Displacement in Z-Direction . . . . . . . . . . . . 98 5.13 Commanded and Actual Velocity in Z-Direction Without Error Compen- sation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 99 5.14 Commanded and Actual Velocity in Z-Direction . . . . . . . . . . . . . . . 100 5.15 Commanded and Actual Displacement in Z-Direction Without Error Com- pensation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 101 5.16 Commanded and Actual Rotation Angle of Spacecraft . . . . . . . . . . . 102 5.17 Commanded and Actual Rotation Angle of Spacecraft Without Error Com- pensation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 103 5.18 End De ection of Flexible Structure Position on Central Body's x-axis . 104 5.19 End De ection of Flexible Structure Position on Central Body's x-axis Without Error Compensation . . . . . . . . . . . . . . . . . . . . . . . . . 105 5.20 End De ection of Flexible Structure Position on Central Body's x-axis Without Vibration Suppression . . . . . . . . . . . . . . . . . . . . . . . . 106 5.21 Commanded and Applied Forces in Z-Directions . . . . . . . . . . . . . . 107 5.22 Commanded and Actual Displacement in X-Direction . . . . . . . . . . . 110 5.23 Commanded and Actual Displacement in Y-Direction . . . . . . . . . . . 111 5.24 Commanded and Actual Displacement in Z-Direction . . . . . . . . . . . . 112 5.25 Commanded and Actual Velocity in X-Direction . . . . . . . . . . . . . . 113 5.26 Commanded and Actual Velocity in Y-Direction . . . . . . . . . . . . . . 114 5.27 Commanded and Actual Velocity in Z-Direction . . . . . . . . . . . . . . . 115 viii 5.28 Commanded and Actual Rotation Angle of Spacecraft . . . . . . . . . . . 116 5.29 End De ection of Flexible Structure Position on Central Body's x-axis . 117 5.30 End De ection of Flexible Structure Position on Central Body's x-axis Without Vibration Suppression . . . . . . . . . . . . . . . . . . . . . . . . 118 5.31 Commanded and Applied Forces in X-Directions . . . . . . . . . . . . . . 119 5.32 Commanded and Applied Forces in Y-Directions . . . . . . . . . . . . . . 120 5.33 Commanded and Applied Forces in Z-Directions . . . . . . . . . . . . . . 121 5.34 Commanded and Applied Torques about X-Axis . . . . . . . . . . . . . . 122 5.35 Commanded and Applied Torques in Y-Axis . . . . . . . . . . . . . . . . . 123 5.36 Commanded and Applied Torques in Z-Axiz . . . . . . . . . . . . . . . . . 124 5.37 Commanded and Actual Displacement in Z-Direction . . . . . . . . . . . . 127 5.38 Commanded and Actual Velocity in Z-Direction Without Error Compen- sation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 128 5.39 Commanded and Actual Velocity in Z-Direction . . . . . . . . . . . . . . . 129 5.40 Commanded and Actual Displacement in Z-Direction Without Error Com- pensation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 130 5.41 End De ection of Flexible Structure Position on Central Body's x-axis . 131 5.42 End De ection of Flexible Structure Position on Central Body's x-axis Without Error Compensation . . . . . . . . . . . . . . . . . . . . . . . . . 132 5.43 End De ection of Flexible Structure Position on Central Body's x-axis Without Vibration Suppression . . . . . . . . . . . . . . . . . . . . . . . . 133 5.44 Commanded and Applied Forces in Z-Directions . . . . . . . . . . . . . . 134 6.1 Directions of Vibration in Three Dimensions . . . . . . . . . . . . . . . . . 136 C.1 Translation and Rotation of Spacecraft in a Planar Reference Frame . . . 162 C.2 Beam De ection . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 164 ix Abstract A methodology for the design of robust control systems for exible spacecraft undergo- ing a large angle maneuver using pulse-modulated thrusters is proposed. The need for developing such a methodology stems from the need to reduce the weight of the mod- ern spacecraft and to satisfy strict performance specications resulting in vibrational modes within the control system bandwidth. Moreover, the discontinuous operation of pulse-modulation may excite large exible motion that can lead to deterioration of perfor- mance, limit cycles, and even instability. Consequently, the current practice of designing the control system using a rigid body model and evaluating performance by simulation using a exible model is inadequate for the modern spacecraft. In the proposed method- ology structural exibility is included in the control design model thus resulting in control systems with robust stability and improved performance. This dissertation has two main parts: formulation of exible spacecraft dynamics performing large angle motion and control system design. The formulation of equations of motion is performed using quasi-coordinates that allow for parameterization of the attitude in terms of Euler-Rodrigues parameters that avoid singularities. Moreover, the formulation is in terms of variables measurable by the control system, such as angular velocities and is in a form that is convenient for establishing closed-loop robust stability. x The control law structure proposed here results in guaranteed closed-loop performance in the presence of parameter uncertainties, measurement noise, and command implemen- tation errors. A theorem that guarantees the bounds of trajectory following errors of a spacecraft performing a general translational and rotational motion in three-dimensional space is formulated and proven using Lyapunov stability theory. Finally, a model of a exible spacecraft is used to demonstrate the performance of a control system developed using the proposed methodology are presented. A number of dierent maneuvers are studied by simulation showing satisfactory performance in the presence of uncertainties and implementation errors. xi Chapter 1 Introduction With an increase in complexity of spacecrafts being launched today, a comprehensive methodology for modeling and control design for analysis that could be applicable to a variety of missions would be benecial when it is important to determine the impact of the controller on the system. In an eort to reduce mission costs, engineers try to reduce the mass of a spacecraft which may result in exibility in some of the spacecraft's structures. Additionally, the spacecraft may contain structures that are inherently exible such as antennas, booms, and solar panels. There is minimal research in control theory that accounts for the coupling that exists between both the translational and rotational degrees of freedom with the exible structures of spacecraft. Most research in controls of a spacecraft with exible structures only accounts for the coupling between the rotational degrees of freedom with the exible structures, which may cause limit cycles in the exible structures and attitude that may take the spacecraft o course. In the derivation of the equations of motion for a spacecraft with exible structures, we have included the translational, rotational and vibrational degrees of freedom, thereby allowing us to account more accurately for the coupling in the control design and stability 1 analysis. The rotational motion is represented using the quaternions, which are dependent coordinates and subject the system to a nonholonomic constraint that may be accounted for using Lagrange multipliers. To obtain a reduced-order model of the system, the equations of motions are derived in terms of quasi-coordinates which are independent, so the Lagrange multipliers are no longer needed. Presented here is an approach for designing a controller for the translational and ro- tational degrees of freedom during trajectory tracking, which accounts for the coupling that exists between all degrees of freedom of the spacecraft. Furthermore, the controller is implemented using pulse-modulated (on-o) thrusters. This theoretical approach to control design incorporates Lyapunov theory, passivity concepts, and active vibration suppression during simultaneous guidance and attitude control. The vibration suppres- sion is necessary for the exible structures because of the disturbances introduced during trajectory tracking due to the coupling, and a low inherent damping in the exible struc- ture may not be enough to compensate for such disturbances. Additionally, the controller is designed using the nonlinear model of the spacecraft with exible structures. The com- manded controller is calculated using sensor measurements, the commanded trajectory, and the system's nominal parameters. It is then used to determine control implementa- tion based on some thruster selection and jet ring logic. It also determines the forces to be applied by the actuators for damping vibration. However, given that the control implementation is modulated and discontinuous, the commanded controller cannot be implemented as desired when attempting to simultaneously control the translational and rotational degrees of freedom using the same pulse-modulated thrusters. Additionally, the discontinuous control implementation introduces bounded disturbances to the system 2 which excites vibration in the exible structures. To account for control implementation uncertainty and any other known bounded uncertainties such as those due to the envi- ronment, materials used in building the spacecraft, sensor measurement inputs to the commanded controller, and variations in the system's parameters, an error compensa- tion term is added to the commanded control input which is designed using Lyapunov Redesign concepts. To compensate for disturbances that excite vibration in the exi- ble structures, collocated actuators and sensors are placed on the exible components to increase damping. In addition to the error compensation term, there is an inverse dynamics control term for tracking and a proportional-derivative feedback control term. We have established stability of a spacecraft comprised of a rigid central body with six degrees of freedom and exible structures that are susceptible to latitudinal vibration in two directions. 1.1 Review of Previous Work Various design methods for attitude control of rigid spacecraft have been investigated us- ing the dierent parameters for attitude representation, such as Euler angles and quater- nions. Wen et al. [53] have proved global asymptotic stability using Lyapunov stability theory with quaternions for attitude tracking using proportional-derivative feedback and feed-forward compensation. Wie et al. [54, 55] have proved stability of three-axis large angle feedback maneuvers with control implementation using pulse modulated reaction jets. Sliding mode controllers, which aim to constrain the states to a manifold or sur- face, have also been used for attitude control. Singh and Iyer [44] have derived a sliding 3 mode controller such that the controlled outputs of the closed-loop system track a given trajectory in the presence of bounded disturbances. Jing and Xu [23] have used the unit error quaternion for nonlinear attitude tracking where the control is implemented using an on-o controller where the equations of motion are in terms of the error quaternion and its derivative. There has also been research done in simultaneous control of the rotational and trans- lational degrees of freedom for rigid spacecraft with on-o thrusters. Thurman [49] has developed a theoretical control design methodology for a rigid spacecraft that motivated this research topic. In his approach, the coupling terms are incorporated with the distur- bances, and the inverse dynamics and proportional-derivative control terms are designed for the linearized system. To account for uncertainty from the model assumptions, sensor measurements, and control implementation, an error compensation term is designed using concepts presented by Corless [9]. The error compensation term accounts for bounded disturbances in a linear system and guarantees exponential convergence of the error to a ball around the origin. Curti et al. [10] have developed a theoretical methodology using on-o thrusters for guidance and attitude control using a Lyapunov-based thruster selec- tion to determine the minimum number of thrusters needed at each time step to maintain stability while tracking a linear reference model. The controller is tested via simulation of a 6-DOF model and experimentally using a 3-DOF model. Motivation for research in control of exible spacecraft include wanting to reduce launch costs, thereby resulting in a reduction of the spacecraft's mass, which in turn may cause exibility is some of the spacecraft's structures. Or, the spacecraft may have structures that are inherently exible, and designing a controller which incorporates the 4 dynamics of the exible structure into the model allows us to account for the impact of the controller on the entire system because it is a hybrid of states for position, attitude and innite dimensional modal amplitudes which are all coupled [3]. Most research in controls of spacecrafts with exible structures only accounts for the coupling that exists between attitude control and active vibration suppression. Using passivity concepts for control design exploit the known physical structure of the system [35]. Researchers have used passivity concepts to design a feedback controller for point-to-point attitude control of exible multi-body system that is implemented con- tinuously to the nonlinear system [30], [25]. When not accounting for any uncertainties in the system during point-to-point trajectory tracking, the closed-loop system is glob- ally asymptotically stable; due to the inherent damping in the exible structures, active vibration suppression is not necessary. Additionally, Kelkar and Joshi [29] have stated without proof that it would only be possible to show boundedness for continuous trajec- tory tracking. However, due to coupling between the rigid and exible degrees of freedom, it is necessary to implement some form of vibration suppression during continuous tra- jectory tracking due to the resulting disturbance that is introduced in the states for exibility. Erdong et al. [15] have used passivity concepts to design a feedback controller to account for arbitrary disturbances entering the system through torque implementation to show the closed-loop system is globally asymptotically stable. Furthermore, the system remains globally asymptotically stable even when the control law is implemented using on-o type thrusters. Sliding mode techniques have also been used for attitude control of a exible spacecraft [3]. Hu [19] has proved asymptotic stability using a sliding mode controller for attitude control in conjunction with sensors and actuators for vibration 5 suppression. Optimal control for a single-axis maneuver is implemented by a continuous torque to a linearized model of a exible spacecraft [52]. Fulton [24] researched optimal attitude control of small exible spacecrafts. Under the assumption that vibration in the exible appendages is bounded when implementing active vibration suppression, a sliding mode controller can be implemented such that the tracking error is asymptotically stable [22]. Due to coupling of the rigid and exible dynamics, trajectory tracking causes disturbances impacting the exible elements thereby requiring active vibration suppression. Benhabib et al. [5] have proven input- output stability of the modal amplitudes when collocated actuators and sensors are used for vibration suppression. Experimental research has also be done using collocated and non-collocated piezoelectric sensors and actuators for vibrations suppression in [12], [47], [20], and [17]. 1.2 Dissertation Outline In the following sections, the mathematical model is derived and a controller is designed for attitude and guidance of a spacecraft with exible structures. The controller is de- signed to compensate for known bounded disturbances including those due to discon- tinuous control implementation such as on-o thrusters. In Chapter 2, we present the mathematical backgrounded needed for deriving the equations of motion for the system and nonlinear control stability concepts. In Chapter 3, we derive a reduced-order math- ematical model of the system using quasi-coordinates and the Kronecker product. In Chapter 4, we show the proof of stability for tracking a trajectory continuously in the 6 presence of disturbances with a nonlinear controller implemented discontinuously. Chap- ter 5 includes a simulation of the spacecraft tracking a various trajectories. Finally, we discuss the results presented and topics of future research. 7 Chapter 2 Mathematical Preliminaries 2.1 Matrix Theory The following denitions and theorems pertaining to matrices are from Bernstein's Matrix Mathematics [6]. Denition 2.1.1. The symbols R and C denote the real and complex number elds containing scalar elements, respectively. The symbols R n and C n denote the real and complex vector elds whose elements are n-dimensional column vectors. The symbols R mn and C mn denote the real and complex matrix elds whose elements are mn matrices. The symbolF is used to denoteR orC, when they elements may be either real or complex. Similarly,F n is used forR n orC n andF mn is used forR mn orC mn . Notation 2.1.1. LetH n ,N n , andP n denote the set of Hermitian, positive semidenite, and positive denite matrices inF nn , respectively. The following proposition contains properties of positive (semi)denite matrices. Proposition 2.1.1. Let matrices A, B, C, D2H n and S2F mn . Then the following statements hold: 8 (i) If A 0 and B> 0, then A + B> 0. (ii) A 2 0. (iii) A 2 > 0 if and only if det A6= 0. (iv) If A B and B< C, then A< C. (v) If A< B and B C, then A< C. (vi) If A B and C D, then A + C B + D. (vii) If A B and C< D, then A + C< B + D. (viii) If A B, then SAS SBS . (ix) If A< B and rank S =m, then SAS < SBS . (x) If SAS SBS and rank S =n, then A B. (xi) If SAS < SBS and rank S =n, then m =n and A< B. The following proposition is used to determine the rank of a product of matrices. Proposition 2.1.2. Let A2F nm and let C2F kn be left invertible and B2F ml be right invertible. Then rank A = rank CA = rank AB: (2.1) Here we have the necessary conditions for showing a matrix is positive (semi)denite. Proposition 2.1.3. Let A, 2 6 6 4 A 11 A 12 A 12 A 22 3 7 7 5 2H n+m . The the following statements are equivalent: 9 (i) A 0. (ii) A 11 0, A 12 = A 11 A 11 1 A 1 2 and A 12 A 11 1 A 12 A 22 . (iii) A 22 0, A 12 = A 12 A 22 A 22 1 and A 12 A 22 1 A 12 A 11 . The following statements are also equivalent: (i) A> 0. (ii) A 11 > 0 and A 12 A 11 1 A 12 < A 22 . (iii) A 22 > 0 and A 12 A 22 1 A 12 < A 11 . 2.2 Cross Product The symbolJK denotes the cross product matrix of a vector, i.e. for u; v2R 3 , uv =JuKv = 2 6 6 6 6 6 6 4 0 u 3 u 2 u 3 0 u 1 u 2 u 1 0 3 7 7 7 7 7 7 5 2 6 6 6 6 6 6 4 v 1 v 2 v 3 3 7 7 7 7 7 7 5 (2.2) Properties of cross-products: for a, b, c2R 3 ab =ba (2.3) a (b +c) = (ab) + (ac) (2.4) a (bc) =b (ca) =c (ab) (2.5) ab +cd = (ac) (bd) +ad +cb (2.6) 10 a (bc) =b (ac)c (ab) (2.7) (ab) (ac) = (a (bc))a (2.8) 2.3 Quaternions Let q I b 2R 4 be the unit quaternion describing the rotation from the inertial coordinate frame to the body coordinate frame with a kinematical relationship to ~ ! b I b 2 R 3 , the angular velocity of the body with respect to the inertial coordinate frame expressed in the body coordinate frame. For simplicity, ~ ! b I b will be denoted by ~ ! and q I b will be denoted as q = 2 6 6 4 q 0 ~ q 3 7 7 5 = 2 6 6 4 cos 2 ~ e sin 2 3 7 7 5 = 2 6 6 6 6 6 6 6 6 6 6 4 q 0 q 1 q 2 q 3 3 7 7 7 7 7 7 7 7 7 7 5 : (2.9) Note that~ q; ~ e2R 3 andj~ ej = 1. Physically, the angle represents the angular rotation of the spacecraft about the unit vector~ e in the body coordinate frame. The unit quaternion is subject to the following constraints q T q =q 2 0 +q 2 1 +q 2 2 +q 2 3 = 1 (2.10) _ q T q +q T _ q = 0! _ q T q = 0 (2.11) 11 q T q + _ q T _ q = 0 (2.12) Additionally, the complex conjugate is equal to the inverse of the quaternion, that is q =q 1 = 2 6 6 4 q 0 ~ q 3 7 7 5 (2.13) where q is the complex conjugate and q 1 denotes the quaternion inverse. The derivative of the quaternion _ q = 1 2 q! = 1 2 2 6 6 4 ~ ! T ~ q q 0 ~ ! +~ !~ q 3 7 7 5 = 1 2 L T ~ ! (2.14) illustrates the kinematical relationship between the quaternion and the angular velocity where denotes quaternion multiplication, ! = 0 ~ ! T T and the matrix L is dened in Equation A.18. 2.4 Dynamics 2.4.1 Lagrangian Mechanics Let the states of any dynamical system be specied by the vectorq = q 1 q 2 q n T 2 R n , whereq i fori = 1, 2,:::,n are the independent generalize coordinates of the system. LetT =T ~ q; _ ~ q denote the kinetic energy andV =V (~ q) denote the potential energy of a dynamical system. The Lagrangian is dened as L =TV: (2.15) 12 To determine the equations of motion using Lagrangian mechanics, the Lagrangian is substituted into the Euler-Lagrange equation d dt @L @ _ q i @L @q i =Q 0 i (2.16) for i = 1, 2, :::, n where Q i 's are the generalized forces not included in the potential energy, to obtain the equations of motion for the system. Euler-Lagrange's equation can be rewritten as d dt @T @ _ q i @T @q i + @V @q i =Q 0 i ) d dt @T @ _ q i @T @q i =Q 0 i @V @q i =Q i : (2.17) The latter form of the Euler-Lagrange equation will be used to derive the equations of motion. 2.4.2 Boltzmann-Hamel Equation When using the Euler-Lagrange equation to derive the equations of motion, the La- grangian is in terms of n-generalized coordinates q i and their derivatives _ q i . When a nonholonomic constraint is placed on the the generalized coordinates, like the unit norm constraint on a quaternion used to describe rotational motion, Lagrange multipliers are necessary. The Boltzmann-Hamel equation is an alternative to the Euler-Lagrange equa- tion for deriving the equations of motion in terms of independent quasi-coordinates where there exists a linear relation between the dependent generalized coordinates and the quasi- coordinates [8,37]. 13 To derive the Boltzmann-Hamel equation, consider d'Alembert-Lagrange equation in matrix format [8] f _ qg T d dt @T @ _ q @T @q fQg = 0; (2.18) which is an expression of the principle of virtual work, where q = q 1 q 2 q n T is a vector of the virtual variation of the generalized coordinates q = q 1 q 2 q n T ; the kinetic energy T = 1 2 _ ~ q T M _ ~ q (2.19) for M positive semidenite with the partial derivatives; @T @ _ q = 2 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 4 @T @ _ q 1 @T @ _ q 2 . . . @T @ _ q n 3 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 5 ; and @T @q = 2 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 4 @T @q 1 @T @q 2 . . . @T @q n 3 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 5 : (2.20) The linear relationship between the derivatives of the generalized coordinates q and the quasi-coordinates ! has the form ! = T _ q (2.21) 14 and _ q = !: (2.22) where both and are matrix functions of q and have the property T = T = I (2.23) where I is the nn identity matrix. Substituting Equation 2.22 into the kinetic energy yields T = 1 2 ! T T M !, 1 2 ! T M! (2.24) where M = T M is a positive denite matrix. Using Equation (2.21), the virtual displacement can be written as fqg = f!g (2.25) Substituting Equation 2.24 and 2.25 into Equation 2.18 gives f!g T d dt @T @ _ q @T @q fQg = 0: (2.26) The quasi-coordinates ! are independent; thus in order to satisfy the principle of virtual work for any virtual displacement !, it follows that T d dt @T @ _ q @T @q fQg = 0: (2.27) 15 As a result of the principle of virtual work, the equations of motion can be derived by pre-multiplying the Euler-Lagrange equation by T . 2.5 Stability Concepts 2.5.1 Passivity The following denitions and theorems related to passivity are from Kelkar and Joshi's Control of Nonlinear Multibody Flexible Space Structures [29]. Denition 2.5.1. LetX be the set of measurable real-valued n-vector functions x(t) of time t, that is x :R + 7!R n . Denition 2.5.2. DeneL n p =L n p [0;1). be the set of all x(t)2X such that Z 1 0 jx i (t)j p dt<1 (2.28) for any p2Z + . For p = 2, the inner product of x and y inL n p is dened as hx;yi = Z 1 0 x T (t)y(t) dt (2.29) Denition 2.5.3. Letx (t)2X . The truncation operator P T , is dened as the function P T (x (t)) = 8 > > < > > : x(t) for tT; 0 for t>T: (2.30) 16 We will denote x T (t) = P T (x (t)). Furthermore, for a xed integer p2 Z + , the set L n pe =L n pe [0;1) denotes the set of all functions x(t) whose truncations x T (t) are inL n p for all nite T .L n pe is called the extension of the spaceL n p . The truncated inner producthx;yi T is dened as: hx;yi T =hx T ;y T i = Z T 0 x T (t)y(t) dt (2.31) Denition 2.5.4. A system is said to be internally passive if there exists a nonnegative storage function E [] such that hu;yi T E [x (T )]E [x (0)]8u2L m 2e ;8T 0: (2.32) Denition 2.5.5. A system is said to be internally general passive if there exists a nonnegative number and and a non-negative storage function E [] such that hu;yi T E [x (T )]E [x (0)] +kuk 2 T +kyk 2 T 8u2L m 2e : (2.33) Denition 2.5.6. A system is said to be internally input strictly passive if there exists an and a non-negative storage function E [] such that hu;yi T E [x (T )]E [x (0)] +kuk 2 T 8u2L m 2e ;8T 0: (2.34) 17 Denition 2.5.7. A system is said to be internally general passive if there exists a nonnegative number and and a non-negative storage function E [] such that hu;yi T E [x (T )]E [x (0)] +kyk 2 T 8u2L m 2e ;8T 0: (2.35) The nonlinear system of a exible spacecraft can be described by _ x = f (x) +g (x)u y = h (x) +N (x)u (2.36) where x2 R n , u2U, f (x)2 R n , h (x)2 R l , g (x)2 R nm with f, g, h, N smooth, f (0) = 0, and h (0) = 0. Theorem 2.5.1. The nonlinear system (2.36) is internally passive if and only if there exists a non-negative function E (x)2C 1 , E (0), and functions l (x)2R k , and W (x)2 R km such that r T E (x)f (x) = l T (x)l (x) g T (x)rE (x) = h (x) 2W T (x)l (x) W T (x)W (x) = 1 2 N (x) +N T (x) : (2.37) Denition 2.5.8. The gain of the system :U!Y is dened as: p B () = inf u2U T0 K such that kuk Tp Kkuk Tp +B for some finiteB: (2.38) 18 Figure 2.1: Feedback Interconnection of Nonlinear Systems Denition 2.5.9. The system :U!Y is said to beL p -stable ifu2L m p [0;1))y2 L l p [0;1). is said to be nite-gainL p -stable if nite p B () and B exist. Two com- monly used denitions of stability are bounded-input-bounded-output (BIBO) stability which corresponds toL 1 -stability andL 2 -stability, which implies that every nite-energy input produces a nite-energy output. Theorem 2.5.2. Passivity Theorem Consider the feedback system in Figure 2.1 above. Suppose both G and H are general passive; that is, they satisfy Z T 0 y T i u i i + i Z T 0 u T i u i dt + i Z T 0 y T i y i dt (2.39) for all T 0; for all u i 2U; and for some i , i 0, i 0 for i = 1, 2. Then the feedback system is nite-gain stable if 1 + 2 > 0 and 2 + 1 > 0: (2.40) 19 Corollary 2.5.3. The feedback system in Figure 2.1 is nite-gain stable if both G and H are input strictly passive, or both G and H are output strictly passive. Theorem 2.5.4. Suppose both G and H in Figure 2.1 are zero-state observable and internally general passive; that is, there exist positive denite storage functions E (x i ) for i = 1, 2 such that Z T 0 y T i u i E i (x i (T ))E i (x i (0)) + i Z T 0 u T i u i dt + i Z T 0 y T i y i dt (2.41) for all T 0; for all u i 2U for i = 1, 2 along the trajectories of the system. Then the origin (0; 0) is asymptotically stable if G and H are zero-state observable and 1 + 2 > 0 and 2 + 1 > 0: (2.42) Corollary 2.5.5. The G and H are passive, and satisfy any one of the following condi- tions,then the origin of the closed-loop system is asymptotically stable: both G and H are input strictly passive and zero-state observable; both G and H are output strictly passive and zero-state observable; either G or H is zero-state observable and either G is output strictly passive or H is input strictly passive; (GH) is zero-state observable and either G is input strictly passive or H is output strictly passive. 20 For nite-dimensional linear, time-invariant systems, passivity is equivalent to "pos- itive realness" of the transfer functions [29]. The following denitions for "positive re- alness" are from Benhabib et al. "Stability of Large Space Structure Control Systems Using Positivity Concepts." Denition 2.5.10. A square transfer matrix Z(s) is called "positive real" if 1. Z (s) has real elements for real s. 2. Z (s) has elements which are analytic for< [s]> 0. 3. Z (s) +Z (s) is non-negative denite for< [s] > 0, where denotes the complex conjugate transpose. Denition 2.5.11. A square transfer matrix Z(s) is called "strictly positive real" if 1. Z (s) has real elements for real s. 2. Z (s) has elements which are analytic for< [s] 0. 3. Z (j!) +Z (j!) is positive denite for all real !. Theorem 2.5.6. If in the feedback system in Figure 2.1, assume r 2 = 0 and G =G (s) and H = H (s) where G (s) and H (s) are both square transfer matrices. The system is asymptotically stable in the input/output sense it at least one of transfer matrices is strictly positive real and to other is positive real. 21 2.5.2 Lyapunov Stability The following denitions and theorems regarding Lyapunov Stability and Boundedness are from Khalil's Nonlinear Systems [31]. Denition 2.5.12. For the autonomous system _ x =f (x); (2.43) let x = 0 be an equilibrium point that is stable, if for each > 0, there exists = ()> 0 such that kx (0)k<)kx (t)k<8t 0 (2.44) unstable, if it is not stable asymptotically stable, if it is stable and can be chosen such that kx (0)k<) lim t!1 x (t) = 0 (2.45) Denition 2.5.13. LetDR n be a domain containingx = 0. The functionV :D!R is positive denite if V (0) = 0 and V (x)> 0 for all x2D 0. The function V :D!R is positive semidenite if V (x) 0 for all x2D. Denition 2.5.14. The V (x) is negative denite or negative semidenite ifV (x) is positive denite or positive semidenite, respectively. 22 Theorem 2.5.7. Let x = 0 be an equilibrium point for (2.43) and DR n be a domain containing x = 0. Let V : D ! R n be a continuously dierentiable positive denite function such that _ V (x) 0 8 x2D; (2.46) that is, _ V (x) is negative semidenite. Then, x = 0 is stable. Moreover, if _ V (x)< 0 8 x2D 0; (2.47) that is, _ V (x) is negative denite; then x = 0 is asymptotically stable. Theorem 2.5.8. Let x = 0 be an equilibrium point for (2.43). Let V : R n ! R be a continuously dierentiable positive denite function such that kxk!1 ) V (x)!1 (2.48) _ V (x)< 0; 8 x6= 0 (2.49) then the equilibrium point is globally asymptotically stable. Theorem 2.5.9. LaSalle's Theorem Let x = 0 be an equilibrium point for (2.43). Let V :R n !R be a continuously dieren- tiable, radially unbounded, positive denite function such that _ V (x) 0 for all x2R n . Let S = n x2R n j _ V (x) = 0 o and suppose that no solution can stay identically in S, other than the trivial solution x(t) 0. Then the origin is globally asymptotically stable. 23 2.5.3 Boundedness Denition 2.5.15. A continuous function : [0;a)! [0;1) is said to belong to class K if it is strictly increasing and (0) = 0. It is said to belong to classK 1 if =1 and (r)!1 as r!1. Denition 2.5.16. A continuous function : [0;) [0;1)! [0;1) is said to belong to classKL if, for each xed s the mapping (r;s) belongs to classK with respect to r and for each xed r, the mapping (r;s) is decreasing with respect to s and (r;s)! 0 as s!1. Denition 2.5.17. The solutions to (2.43) are uniformly bounded if there exists a positive constant c, independent of t 0 0, and for every 2 (0;c), there is =(a)> 0, independent of t 0 such that kx (t 0 )ka)kx (t)k; 8 tt 0 (2.50) globally uniformly bounded if (2.50) holds for arbitrarily large a; uniformly ultimately bounded with ultimate bound b if there exist positive constants b and c, independent of t 0 0, and for every a2 (0;c), there is T = T (a;b) 0, independent of t 0 such that kx (t 0 )ka)kx (t)kb; 8 tt 0 +T (2.51) globally uniformly ultimately bounded if (2.51) holds for arbitrarily large a. 24 Consider the following system _ x =f (t;x) +G (t;x) [u; (t;x;u)] (2.52) where x2 R n is the state and u2 R p is the control input. The functions f, G, and are dened for (t;x;u)2 [0;1)DR p , where DR n is a domain that contains the origin. The function f and G are known precisely, while is a bounded but unknown function that combines the various uncertain terms. Theorem 2.5.10. (Ultimate Boundedness) Let D R n be a domain that contains the origin and V :D!R be a continuously dierentiable function such that 1. 1 (kxk)V (x) 2 (kxk) 2. @V @x f (x)W 3 (x) 8 kxk>> 0 for all x2D where 1 and 2 are class K functions and W 3 (x) is a continuous positive denite function. Take r> 0 such that B r D and suppose that < 1 2 ( 1 (r)) (2.53) Then there exists a class KL function and for every initial state x (t 0 ), satisfying kx (t 0 )k, there isT > 0 (dependent onx (t 0 ) and) such that the solution of _ x =f (x) satises 1. kx (t)k (kx (t 0 )k;tt 0 ) 8 t 0 <t<t 0 +T 2. kx (t)kr = 1 1 ( 2 ()) 8 tt 0 +T 25 Thus, the solution x (t) is uniformly bounded8t t 0 and uniformly ultimately bounded with the ultimate bound r = 1 1 ( 2 ()). Moreover, if D =R n and 1 belongs to class K 1 , then these conditions hold for any initial state x (t 0 ) with no restrictions on how large can be. 26 Chapter 3 Modeling of Closed-Loop Systems in a Three Dimensional Space 3.1 Modeling of Spacecraft with Flexible Structures The rst step in control design of a spacecraft or any other system is to develop a mathematical representation of the system. The exible structures of a spacecraft are distributed-parameter (continuous) systems; that is, they are systems with an innite dimensional state space. Modal analysis can be used to calculate the natural frequencies, modal damping, and mode shapes which are determined by the material properties and boundary conditions of the exible structures and are necessary for deriving the equations of motion for the exible structures. The assumed modes methods is one methodology used to derive the equations of motions using Lagrangian mechanics for the exible struc- tures by discretizing the kinetic energy, potential energy and virtual work. A spacecraft moving in a three dimensional reference plane has six degrees of freedom, and exible structures may be susceptible to longitudinal, lateral, and/or torsional vibration, thereby adding degrees of freedom to an already innitely dimensional state space. In deriving 27 the equations of the spacecraft, we have assumed that the spacecraft has six degrees of freedom for rigid body motion in a three dimensional reference frame. Attached to the central rigid body are N exible structures with latitudinal vibration in two directions where then-th exible structure for i = 1; 2;:::;N hask n vibrational modes included in the model. In this section we begin with a derivation of the mode shapes of the exible structures. We then use the mode shapes to derive the equations of motion for the entire system. Lastly, we discuss the uncertainty introduced via control inputs in the system. 28 Figure 3.1: Model of Spacecraft in 3-Dimensional Reference Frame 29 3.1.1 Modeling of Flexible Structures using Modal Analysis The mode shapes of the beam are necessary for determining the mass and stiness matri- ces for the system. We have assumed that the exible structures are only susceptible to lateral vibration. Thus the beam can be modeled by the Euler-Bernoulli Beam Equation, a fourth-order partial dierential equation (PDE) from which we derive the mode shapes. The Euler Bernoulli Equation is dened as @ 2 @x 2 EI(x) @ 2 y @x 2 +(x) @ 2 y @t 2 =p(x;t) (3.1) where y (x;t) is the lateral displacement of the beam; I(x) is Moment of Inertia of the beam's cross section; E is the Young's Modulus of the beam; y(x;t) is vertical displacement of the beam; (x) is the mass of beam per unit length; p(x;t) is the vertical force on the beam per unit length. Assuming that E, I, and (x) are constant, the beam equation can be reduced to EI @ 4 y @x 4 + @ 2 y @t 2 = 0 (3.2) A cantilever-free beam with length L n is subjected to the following boundary conditions 30 1. y(0;t) = 0 2. y x (0;t) = 0 3. y xx (L n ;t) = 0 4. y xxx (L n ;t) = 0 Figure 3.2: Cantilever Beam Using the technique of separation of variable, assumingy(x;t) =(x)(t) and substituting this into the beam equation gives EI (4) (x)(t) +(x) 00 (t) = 0) EI (4) (x) (x) = 00 (t) q(t) =! 2 k (3.3) where the constant ! k is the natural frequency of the beam. This allows us to separate the PDE of the beam into the following two ODEs: (4) (x) 4 (x) = 0 where 4 = ! 2 k EI (3.4) 00 (t) +! 2 k (t) = 0 (3.5) 31 Equation (3.4) is a fourth-order dierential equation whose solution has the form (x) = P K k=1 k (x) where k (x) =C 1 (cos k x + cosh k x) +C 2 (cos k x cosh k x) +C 3 (sin k x + sinh k x) +C 4 (sin k x sinh k x): (3.6) The coecients C i for i = 1; 2; 3; 4 are determined using the boundary conditions of a cantilever-free beam. The boundary conditions of the clamped end (x = 0) give, 1: k (0) = 2C 1 = 0!C 1 = 0 2: 0 k (0) = 2 k C 3 = 0!C 3 = 0 The boundary conditions of the free end give 3: 00 k (L n ) =C 2 2 k ( cos k L n cosh k L n ) +C 4 2 k ( sin k L n sinh k L n ) = 0 )C 4 =C 2 cos k L n + cosh k L n sin k L n + sinh k L n 4: 000 k (L n ) =C 2 3 k (sin k L n sinh k L n ) +C 4 3 k ( cos k L n cosh k L n ) = 0 ) (sin k L n sinh k L n ) + cos k L n + cosh k L n sin k L n + sinh k L n (cos k L n + cosh k L n ) = 0 ) sin 2 k x sinh 2 k L n + cos 2 k L n + 2 cos k L n cosh k L n + cosh 2 k L n = 0 ) 1 + 1 cos k L n cosh k L n = 0 (3.7) 32 The frequency equation (3.7) is used to solve for the k . The equation for the nth mode shape reduces to k (x) = C k sin k L n + sinh k L n h (sin k L n + sinh k L n )(cos k x cosh k x) (cos k L n + cosh k L n )(sin k x sinh k x) i (3.8) In deriving the equations of motion, the modes shapes are used to determine the mass matrix, and the modal amplitudes add degrees of freedom to the system. Additionally, the coecients C k of the modes shapes are determined by solving k (L n ) = 1 for C k for all k; it follows that C k = sin k L n + sinh k L n sinh k L n cos k L n sin k L n cosh k L n = 1 2 : (3.9) For a cantilevered beam with a square cross section, the mode shapes for latitudinal vibrational motion in the y and z directions are equivalent. 3.1.2 Position and Velocity The model in Figure 3.1 is used to derive the equations of motion. The position vector of k-th element on the n-th beam, denoted by n k in Figure 3.1, relative to the central body's coordinate systemfbg expressed in the central body's coordinate system ~ rn k=b b = ~ Rn =b b + T n b ~ r n k (3.10) 33 where ~ Rn =b b is the position of the joint where the exible structure is attached to the central body with respect to and expressed in the coordinate frame xed to the central body at the central body's center of gravity; T b n is the rotation matrix from the coordinate system xed to then-the exible structure to the body coordinate system; and~ r n k is the position of then k -th element relative to the coordinate system xed to then-th appendage at the joint. The position of n k relative to the inertial coordinate systemfIg expressed in the inertial coordinate system ~ r n k = I I = ~ r b= I I + T b I ~ Rn =b + T n b ~ r n k (3.11) where ~ r b= I I is the position of the central body with respect to and expressed in the inertial coordinate system; and T b I is the rotation matrix from the central body's coordi- nate system to the inertial coordinate system. The position of n k relative to the inertial coordinate system expressed in the central body's coordinate system ~ r n k = I b = T I b ~ r n k = I I = T I b n ~ r b= I I + T b I ~ Rn =b + T n b ~ r n k = b b o = T I b ~ rb =I I + ~ Rn =b b + T n b ~ r n k = b b (3.12) For simplicity, let ~ r, ~ rb =I I (3.13) _ ~ r, _ ~ rb =I I : (3.14) 34 The velocity ofn k relative to the inertial reference frame expressed in the body reference frame _ ~ rn k=I b = _ T I b ~ r + T I b _ ~ r + T n b _ ~ r n k = 2 _ LG T ~ r + T I b _ ~ r + T n b _ ~ r n k = 2L + r _ q + T I b v + T n b _ ~ r n k = T I b 2L + r T n b 2 6 6 6 6 6 6 4 ~ v _ q _ ~ r n k 3 7 7 7 7 7 7 5 = T I b 2L + r T n b n T 2 6 6 6 6 6 6 4 ~ v _ q _ ~ n 3 7 7 7 7 7 7 5 = T I b 2L + r 0 3k 1 :k n1 T n b n T 0 3k n+1 :k N 2 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 4 ~ v _ q _ ~ 1 . . . _ ~ n . . . _ ~ N 3 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 7 5 , N n _ p (3.15) 35 where the G and L matrix is dened in Equation A.17 and A.18, respectively; + r is dened in Equation A.34; and 0 3k i :k j = 0 3 j P n=i kn (3.16) Additionally, let 0 k i :k j 3 = 0 j P n=i kn3 (3.17) 0 k i :k j = 0 j P n=i kn j P n=i kn ; (3.18) k = N X n=1 k n ; (3.19) that is, k is the total number of modal amplitudes included in the model. 3.1.3 Kinetic Energy of the Flexible Structure The kinetic energy of the n-th exible structure T n = 1 2 Z mn ~ vn k=I T b ~ vn k=I b dm n (3.20) = 1 2 Z mn _ p T N n T N n _ p dm n (3.21) 36 Consider the term N n T N n = 2 6 6 6 6 6 6 6 6 6 6 6 6 6 6 4 T b I 2 + r T L T 0 k 1 :k n1 3 n T b n 0 k n+1 :k N 3 3 7 7 7 7 7 7 7 7 7 7 7 7 7 7 5 T I b 2L + r 0 3k 1 :k n1 T n b n T 0 3k n+1 :k N (3.22) The kinetic energy of the n-th structure can than be rewritten as T n = 1 2 Z mn _ p T N n T N n _ p dm n = 1 2 _ p T Z mn N n T N n dm n _ p = 1 2 _ p T 2 6 6 4 M n 11 M n 12 M T n 12 M n 22 3 7 7 5 _ p = 1 2 _ p T M n _ p (3.23) where M n 11 = 2 6 6 4 m n I 3 2m n T b I L + r 2m n + rL T T I b 4m n + r T L T L + r 3 7 7 5 (3.24) M n 12 = 0 3k 1 :k n1 T b I T n b R mn n T dm n 0 3k n+1 :k N (3.25) 37 M n 22 = 2 6 6 6 6 6 6 4 0 k 1 :k n1 R mn n n T dm n 0 k n+1 :k N 3 7 7 7 7 7 7 5 (3.26) 3.1.4 Elastic Potential Energy of the Flexible Structure The potential energy of the n-th structure due to exibility V ns = 1 2 T n K n n (3.27) where n = T ny T nz T has length n which is the sum of the number of vibrational modes in the y and z direction to be included in the model. The stiness matrix K n is derived as follows. First, consider the integral for the elastic potential energy in the y-direction due to latitudinal vibration: Z Ln 0 EI d 2 z n dx 2 n 2 dx n = Z Ln 0 EI T ny (t) 00 ny (x n ) 00 ny (x n ) T ny (t) dx n = T ny (t) Z Ln 0 EI(t) 00 ny (x n ) 00 ny (x n ) T dx n ny (t) , T ny K ny ny : (3.28) Assuming that the exible structure has a square cross section, the elastic potential energy is the same in the y and z-direction, that is Z Ln 0 EI d 2 z n dx 2 n 2 dx n = Z Ln 0 EI d 2 y n dx 2 n 2 dx n : (3.29) 38 Therefore assuming only latitudinal vibration in the y and z-direction, K nz = K ny , and the stiness matrix for the n-th exible structure K n = 2 6 6 4 K ny K nz 3 7 7 5 : (3.30) 3.1.5 Kinetic Energy of the Rigid Central Body The translational kinetic energy of the central body T btrans = 1 2 m b _ ~ r T _ ~ r: (3.31) where m b is the mass of the central body. The rotational kinetic energy of the central body T brot = 1 2 ! b I T b (J b ) b ! b I b = 1 2 4 _ q T L T (J b ) b L _ q = 1 2 _ q T J b _ q (3.32) where (J b ) b is the inertia tensor of the central body, and J b = 4L T (J b ) b L. 39 3.1.6 Gravitational Potential Energy of the Spacecraft The gravitational potential energy of the spacecraft V g =m b ~ g ~ rbcg=I I + N X n=1 m n ~ g h ~ rbcg=I I + T b I ~ rncg=b b i (3.33) where m n is the mass of the n-th structure; ~ g is the gravitational acceleration vector expressed in the inertial coordinate system; ~ rbcg=I I is the position vector of the central body's cg relative to the inertial coordinate system expressed in the inertial coordinate system; and the position vector of the n-th structure's cg relative to the central body's coordinate system expressed in the central body's coordinate system is ~ rncg=b b and is taken as a constant vector since we assume small vibration in the exible structures. Let~ r n = ~ rncg=b b . The gravitational potential energy can be rewritten as V g =m b ~ g ~ rbcg=I I + N X n=1 m n ~ g h ~ rbcg=I I + T b I ~ rncg=b b i =m b (~ g) I I 3 0 34 0 3k p + N X n=1 m n (~ g) I h ~ rbcg=I I + GL T ~ r n i =m b (~ g) I I 3 0 34 0 3k p + N X n=1 m n (~ g) I h ~ rbcg=I I + G r n q i =m b (~ g) I I 3 0 34 0 3k p + N X n=1 m n (~ g) I I 3 G r n 0 3k p =m b (~ g) T I E b p + N X n=1 m n (~ g) T I E n p: (3.34) 40 3.1.7 Energy of the Spacecraft The kinetic energy of the spacecraft T = T btrans +T brot + N X n=1 T n = 1 2 m b ~ v T ~ v + 1 2 _ q T J b _ q + N X n=1 1 2 _ p T M n t _ p = 1 2 _ p T Mp (3.35) = 1 2 _ p T 2 6 6 4 M 11 M 11 M T 12 M 22 3 7 7 5 p (3.36) where M 11 = 2 6 6 4 (m b +M n ) I 3 2M N T b I L + r 2M N + r T L T T I b 4M N + r T L T L + r + J b 3 7 7 5 ; (3.37) M 12 = 2 6 6 4 T b I T T 1 b R m 1 1 dm 1 T b I T T N b R m N N dm N 2 + rL T T b 1 R m 1 1 dm 1 2 + rL T T b N R m N N dm N 3 7 7 5 ; (3.38) M 22 = 2 6 6 6 6 6 6 6 4 R m 1 1 T 1 dm 1 . . . R m N N T N dm N 3 7 7 7 7 7 7 7 5 ; (3.39) M N = N X n=1 m n : (3.40) 41 The potential energy of the spacecraft V =V g +V s =V bg + N X n=1 +V ng + N X n=1 (V ns ) (3.41) =m b (~ g) T I E b p + N X n=1 m n (~ g) T I E n p + T n K n n : (3.42) 3.1.8 Derivation of the Equations of Motion The kinetic energy of the system is given by T = 1 2 _ p M _ p: (3.43) Using the kinematical relationship from Equation 2.14 between the quaternion and an- gular velocity, the velocity vector _ p can be rewritten as _ p = 2 6 6 6 6 6 6 4 _ r _ q _ 3 7 7 7 7 7 7 5 = 2 6 6 6 6 6 6 4 I 3 1 2 L T I k 3 7 7 7 7 7 7 5 2 6 6 6 6 6 6 4 _ r ! _ 3 7 7 7 7 7 7 5 , v: (3.44) wherev is the quasi-coordinate velocity vector and v can also be written in terms of _ p as v = 2 6 6 6 6 6 6 4 I 3 2L I k 3 7 7 7 7 7 7 5 _ p, T _ p (3.45) 42 Substituting _ p = into the kinetic energy gives T = 1 2 v T T Mv = 1 2 v T Mv (3.46) where the mass matrix M = T M is positive denite. This expression for the kinetic energy is then substituted into Equation 2.24, that is T d dt @T @ _ q @T @q fQg (3.47) to derive the equations of motion. We will calculate the partial derivatives of the kinetic energy using the product rule dened using the Kronecker product. Then we manipulate the terms following the methodology presented by Lewis et al. [16] and properties dened in Appendices B and A to determine the Coriolis/centripetal/coupling matrix. Consider rst, the term @T @ _ p = @ @ _ p 1 2 v T M v = 1 2 I s v T M @v @ _ p + @ @ _ p v T M v = 1 2 I s v T M I s T @ _ p @ _ p + I s v T @M @ _ p v + Mv = 1 2 I s v T M T @ _ p @ _ p + I s v T @M @ _ p v + Mv = 1 2 I s v T M T @ _ p @ _ p + Mv = Mv (3.48) 43 Now consider d dt @T @ _ p = _ Mv + _ Mv + M _ v = M _ v + _ M + _ M v (3.49) Finally consider @T @p = 1 2 @ @p v T Mv = 1 2 I s v T M @v @p + I s v T @M @p v + @v T @p Mv = 1 2 I s v T M I s _ T @p @p + I s v T @M @p v _ Mv = 1 2 I s v T M _ T @p @p + I s v T @M @p v _ Mv = 1 2 _ Mv + I s v T @M @p v _ Mv = 1 2 I s v T @M @p v _ Mv (3.50) The equations of motion, as a result of substituting equations (3.48), (3.49), (3.50) into equation (2.27), are T d dt @T @ _ p @T @p = T M _ v + _ Mv + _ Mv 1 2 I s v T @M @p v _ Mv = T M _ v + T _ Mv + T _ Mv 1 2 T I s v T @M @p v + T _ Mv = M _ v + 2 T _ Mv + _ Mv 1 2 T I s v T @M @p v (3.51) 44 Now we are left to manipulate the Coriolis/centripetal/coupling term C (!; _ p;p)v 2 T _ Mv + _ Mv 1 2 T I s v T @M @p v (3.52) so that the following matrix 1 2 _ M C (3.53) is a skew-symmetric matrix. Let U := T I s v T @M @p : (3.54) Then U T v = @M @p T (I s v) z = @M @p T (I s v) _ p = @M @p T ( _ p I m )v = _ M T v = _ Mv (3.55) Now consider the term T _ Mv = 2 6 6 6 6 6 6 4 0 q 1 2 ! y 0 3 7 7 7 7 7 7 5 (3.56) 45 where T _ = 2 6 6 6 6 6 6 4 0 L _ L T 0 3 7 7 7 7 7 7 5 = 2 6 6 6 6 6 6 4 0 q 1 2 ! y 0 3 7 7 7 7 7 7 5 (3.57) , Mv = 2 6 6 6 6 6 6 4 1 2 3 3 7 7 7 7 7 7 5 (3.58) Then T _ = 2 6 6 6 6 6 6 4 0 q 1 2 ! y 2 0 3 7 7 7 7 7 7 5 = 2 6 6 6 6 6 6 4 0 1 2 J 2 K! 0 3 7 7 7 7 7 7 5 = 2 6 6 6 6 6 6 4 0 1 2 J 2 K 0 3 7 7 7 7 7 7 5 2 6 6 6 6 6 6 4 _ r ! _ 3 7 7 7 7 7 7 5 ,Av (3.59) whereA is a skew-symmetric matrix. 46 Thereby allowing us to rewrite the Coriolis/centripetal/coupling term as C (!; _ p;p)v 2 T _ Mv + _ Mv 1 2 T I s v T @M @p v = 2 T _ Mv + _ M n v 1 2 Uv = 2Av + 1 2 h 2 _ Mv Uv i = 2Av + 1 2 h _ Mv + U T v Uv i = 1 2 h _ M + 4A + U T U i v (3.60) Property 3.1.1. The matrix N (r; _ r) = 1 2 _ M (r)C (r; _ r) (3.61) is skew-symmetric where _ M (r) is the derivative of the mass matrix, and therefore has the property z T N (r; _ r)z = 0 (3.62) for any vector z2R n . For any 2R m , consider the term 1 2 T M T C (!; _ p;p) = T 1 2 _ M C = T 1 2 _ M 1 2 h _ M + 4A + U T U i = T 4A + U T U = 0 47 because the matrix 4A + U T U is skew-symmetric. Finally, T d dt @T @ _ p @T @p = M (p) _ v + C (p;v)v: (3.63) Now consider Q =Q 0 + @V @p (3.64) whereQ 0 is vector of the generalized forces and @V @p is a vector of the potential forces. The generalized forces are those that are applied externally to the spacecraft. The potential forces due to gravity and elasticity are determined as follows. The gravitational forces impacting the system can be determined by taking the partial derivative of the gravitational potential energy with respect to the position state vector p using the Kronecker product, that is @V g @p = I s m b ~ g T E b @p @p + N X n=1 I s m n ~ g T E n @p @p +m n I s ~ g T @E n @p p = I s m b ~ g T I E b vec (I s ) + N X n=1 I s m n ~ g T E n vec (I s ) +m n I s ~ g T @E n @p p ,G (p) (3.65) From the potential energy of the spacecraft, we can derive the stiness and gravity term, that is @V @q = Kp +G(p): (3.66) 48 These terms then need to be transformed into quasi-coordinates. Since the stiness matrix K only multiplies the modal amplitudes, with respect to the quasi-coordinates, the stiness term can be rewritten as K 2 6 6 4 0 61 3 7 7 5 = 2 6 6 4 0 61 ~ K 3 7 7 5 : (3.67) where K = T K, and ~ K is a block diagonal matrix of the stiness matrices for each exible structure. The stiness matrix K and the mass matrix M is then used to generate the damping matrix D. The gravity term simply becomes G(p) = T G(p): (3.68) The equations of motion of the spacecraft can then be written as M _ v + C (p;v)v + Dv + K 2 6 6 4 0 61 3 7 7 5 +G(p) = B (u applied ): (3.69) The forces and torques applied to the model are not the same as those commanded by the controller. We assume that the applied forces are an output of some thruster selection and jet ring logic. As a result the applied forces u applied =u c +u t (3.70) 49 whereu t is the uncertainty introduced to the system because the commanded controller u c can not be applied at desired. For control design and stability analysis, we use the following equations of motion M _ v + C (p;v)v + Dv + K 2 6 6 4 0 61 3 7 7 5 +G(p) = B (u c +u t ): (3.71) 3.2 Implementation of Commanded Control Here we examine the control system as it is applied to the spacecraft to determine the uncertainty introduced to the system via the control input. Figure 3.4 is a general rep- resentation of the control system for a spacecraft. The errors in position, attitude, and modal amplitude are sent to the control law, which determines the necessary forces and torques for tracking and forces to compensate for the vibration in order to maintain sta- bility during trajectory tracking. The thruster selection and jet ring logic algorithm determines which thrusters are turned on, the magnitude of the thrust to be applied, and the duration of the thrust ring. These commands are then sent to the thrusters. 3.2.1 Attitude and Guidance Control Law The guidance and attitude control laws determine the required forces and torques, re- spectively, necessary for tracking a given trajectory. The input to the control law is the error between the desired and measured values of the attitude, angular velocity, position, and velocity. Due to singularities that arise from using the Euler angles, quaternions are 50 often used for attitude tracking. Additionally, the inertial measurement unit measures the rates of acceleration and attitude changes which are then used to calculate the atti- tude, angular velocity, position and velocity. The measurements are susceptible to noise which accumulates with each new measurement, therefore the control law may need an additional term to compensate for the uncertainty in the measurements. 3.2.2 Thruster Selection and Jet Firing Logic The required forces and torques for tracking determine the combination of thrusters and level of thrust required of each of the thrusters that will best satisfy the force and torque requirements. Development of the thruster selection and jet ring logic is not a trivial task due to the various possible arrangements of the thrusters and the number of thrusters utilized. Most likely the force and torque requirements will not be satised exactly as desired thereby introducing additional uncertainty to the system. Ideally, the logic used for selecting the thrusters and determining the thrust levels would try to minimize the uncertainty. The ring logic prole in Figure 3.3 ignites the thrusters when the commanded thrust is greater than half of some specied force A and shuts down the thruster when the commanded thrust is less than half of A. Figure 3.5 illustrates the conguration of the thrusters on the Mars Lander. 51 Figure 3.3: Firing Logic for Thruster 52 Figure 3.4: Attitude and Guidance Control System of Spacecraft 53 Figure 3.5: Mars Lander Cruise Conguration [49] 54 In order to quantify the uncertainty due to the pulse-modulated thrusters, we need to determine a mapping from the commanded controller to the thrusters. The commanded controller species forces in and torques about the x, y, and z-direction in the central body's coordinate frame, that is u c b = 2 6 6 4 f c b c b 3 7 7 5 : (3.72) We assume that there are thrusters that can apply forces in the x, y, andz-direction placed in six locations around the central body. Assuming the central body is a sphere, we use spherical coordinates to specify the locations of the thrusters and determine the position vector ~ r b t i for i = 1; 2;:::; 6 of the i-th thruster relative to the body expressed in the body coordinate frame. The forces and torques applied to the spacecraft can be calculated as 2 6 6 4 f a a 3 7 7 5 = 2 6 6 4 I 3 I 3 q ~ r b t 1 y q ~ r b t 6 y 3 7 7 5 2 6 6 6 6 6 6 4 t 1 . . . t 6 3 7 7 7 7 7 7 5 (3.73) = T b t 2 6 6 6 6 6 6 4 t 1 . . . t 6 3 7 7 7 7 7 7 5 : (3.74) The pseudo-inverse of the matrix T b t is used to determine the thruster commands which in turn determine the thrust that will be applied based on the ring logic, that is t c = T b t + u c b (3.75) 55 where the superscript + denotes the pseudo-inverse and t c is the vector of thruster com- mands. 3.3 Thruster Dynamics When commands from the thruster selection and jet ring logic are sent to the thrusters, there exists a delay in the ignition and termination due to the mechanics of the thruster, that is the thrusters are not instantly turn on or o. Time is needed to build up to the desired thrust level or to shut down the thruster. Figure 3.6 illustrates the thrust prole of a gas jet or thruster ring [49]. Figure 3.6: Gas Jet/ Thruster Firing Prole 3.4 Implementation of Vibration Suppression Consider the input applied by the actuators placed along the exible structure. The equations of motion are functions of the modal amplitudes which can not be measured 56 by the sensors or controlled by actuators. However, it is possible to measure the displace- ment and/or velocity of nodes along the exible structure and place actuators to damp vibration. To determine the impact of the actuators placed along the exible appendages on the system, we use the concept of virtual work to determine the input matrix for the actuators, which we'll denote ~ B f . Assume there are k actuators placed along the nth exible structure, where x n i is the position of the i-th actuator relative to the point of attachment between the n-th exible structure and central body and y n i denotes the de- formation of the exible structure at the location of the actuator. Each actuator applies a force u n i for i = 1; 2;:::;k and the vector of forces applied to the nth appendage is u n = 2 6 6 6 6 6 6 6 6 6 6 4 u n 1 u n 2 . . . u n k 3 7 7 7 7 7 7 7 7 7 7 5 (3.76) and u f = 2 6 6 6 6 6 6 6 6 6 6 4 ~ u 1 ~ u 2 . . . ~ u N 3 7 7 7 7 7 7 7 7 7 7 5 (3.77) 57 is a vector of all inputs from all the actuators on all N exible structures. The virtual work of the forces applied by the actuators on the nth appendage is W n = k X i=1 u n i (@y n i +@ [(x n i +R b ) sin n ]) = k X i=1 u n i " M X m=1 nm (x n i )@ nm (t) # + k X i=1 u n i (x n i +R b ) cos n @ (3.78) whereM is the number of modes we want to control. It follows that the input matrix for the vibration suppression forces applied to the n-th exible structure is ~ B f = 2 6 6 6 6 6 6 6 6 6 6 6 6 6 6 4 (x 1 +R b ) cos n (x 2 +R b ) cos n (x k +R n ) cos n n 1 (x 1 ) n 1 (x 2 ) n 1 (x k ) n 2 (x 1 ) n 2 (x 2 ) n 2 (x k ) . . . . . . . . . . . . n M (x 1 ) n M (x 2 ) n M (x k ) 3 7 7 7 7 7 7 7 7 7 7 7 7 7 7 5 (3.79) ~ B f u f = 2 6 6 4 N P n=1 k P i=1 (x n i +R b ) cos n u n i u f 3 7 7 5 (3.80) To account for errors in implementation by the actuators or some modulation as with the thrusters, let u fa = u fc +u f (3.81) 58 whereu fa is the vector of forces applied by the actuators;u fc is the vector of commanded forces sent to the actuators; and u f denotes the known bounded uncertainty in imple- mentation. Then ~ B f u fa = 2 6 6 4 fa ~ u fa 3 7 7 5 (3.82) where fa = N X n=1 k X i=1 (x n i +R b ) cos n (u n i +u n i ) (3.83) ~ u fa = (u f +u) (3.84) From this analysis we see that the actuators also apply a torque on the system, fa . The total torque applied a = fa + ta + (3.85) where is additional torque applied by a control moment gyro or any other attitude control device. The uncertainty in torque implementation is then dened as = a c (3.86) 59 where c is the commanded torque. It follows that the input to the system can be dened as follows ~ Bu I a = 2 6 6 6 6 6 6 4 f I c +f a c + u fc + u f 3 7 7 7 7 7 7 5 = ~ B u I c +u (3.87) 60 Chapter 4 Trajectory Tracking Control in a Three Dimensional Space The controller designed for spacecraft consists of three components. Here we discuss each of the components of the commanded controller, demonstrate how they are implemented, and how implementation allows disturbances to enter the system through the control terms. 4.0.1 Controller States and Terms We begin with a review of the states used in determining the dierent components of the controller. Recall that the states used for the model are p = 2 6 6 6 6 6 6 4 ~ r q 3 7 7 7 7 7 7 5 ; v = 2 6 6 6 6 6 6 4 _ ~ r ! _ 3 7 7 7 7 7 7 5 ; v = 2 6 6 6 6 6 6 4 ~ r _ ! 3 7 7 7 7 7 7 5 : (4.1) 61 The commanded trajectory is denoted with the subscript d, that is p d , v d , and _ v d . The following vectors comprised of errors in positions and velocities are needed for the con- troller e = 2 6 6 6 6 6 6 4 e ~ r e ~ q e 3 7 7 7 7 7 7 5 ; e v = 2 6 6 6 6 6 6 4 _ e ~ r ~ ! eq _ e 3 7 7 7 7 7 7 5 : (4.2) The error vector for the translational motion is dened as e ~ r =~ r~ r d (4.3) with derivatives _ e ~ r = _ ~ r _ ~ r d and e ~ r = ~ r ~ r d : (4.4) The error quaternion is dened as e q =q 1 d q = 2 6 6 4 e q 0 e ~ q 3 7 7 5 = 2 6 6 4 q d 0 q 0 +~ q d ~ q q d 0 ~ qq 0 ~ q d ~ q d ~ q 3 7 7 5 ; (4.5) where denotes the quaternion product, and has the derivative _ e q = 1 2 e q ! eq = 2 6 6 4 1 2 e T ~ q ~ ! eq 1 2 e q 0 ~ ! eq +e ~ q ~ ! eq 3 7 7 5 (4.6) 62 The term ! eq = 2 6 6 4 0 ~ ! eq 3 7 7 5 has a kinematical relationship to the error quaternion, where ~ ! eq =~ ! (~ ! d ) b ; (4.7) ~ ! is the angular velocity of the central body with respect to the inertial coordinate frame expressed in the actual body coordinate frame; and the term (! d ) b = 2 6 6 4 0 (~ ! d ) b 3 7 7 5 =e 1 q ! d e q (4.8) can be thought of as the commanded angular velocity expressed in the actual body coordinate frame. The derivative of ~ ! eq is dened as _ ~ ! eq = _ ~ ! _ ~ ! d b : (4.9) Note that we want to suppress any vibration in the exible structures, thus d = _ d = d = 0. Furthermore e =; _ e = _ ; e = (4.10) We begin with the inverse dynamics control term, which is included for tracking a continuous or point-to-point trajectory. The inverse dynamics control term is dened as u 0 = M(p) + C(p;v) _ +G(p) (4.11) 63 where M(p) and C(p;v) are the mass and coupling matrices, respectively; G(p) is the gravity vector; the states for the control term are _ = 2 6 6 6 6 6 6 4 _ ~ r d (~ ! d ) b _ 3 7 7 7 7 7 7 5 K p 2 6 6 6 6 6 6 4 e ~ r e ~ q e 3 7 7 7 7 7 7 5 : (4.12) The gain matrix K p = blkdiag K pr ; K p! ; K p where K pr = pr I 3 ; K p! = p! I 3 ; and K p = p I k with positive scalars pr , p! , and p . Due to uncertainty in data measurements and the parameters characterizing the space- craft, the control terms can not be implemented as desired. The measured states describ- ing the position and velocity are denoted ~ p and ~ v, respectively. The mass matrix, coupling matrix, and gravity matrix determined using estimates of the spacecraft's parameters are denoted ~ M (p), ~ C (p;v), and ~ G, respectively. The commanded inverse dynamics control term, which is determined using the measured data and the estimations of the spacecraft's parameters, is denoted u 0 = ~ M (~ p) ~ + ~ C (~ p; ~ v) _ ~ + ~ G (~ p) (4.13) Assuming there are some known bounds on the data measurement and parameter uncer- tainties, the commanded inverse dynamics control term can be written as u 0 = (M (p) +M) + + (C +C) _ + _ +G (p) +G (4.14) = M(p) + C(p;v) _ +G(p) +u 0 (4.15) 64 The feedback control term is dened as u f =K d + 1 2 K d 2 6 6 6 6 6 6 4 0 31 e ~ q 0 k1 3 7 7 7 7 7 7 5 (4.16) where the gain matrix K d = blkdiag K dr ; K d! ; K d where K dr = dr I 3 ; K d! = d! I 3 ; and K d = d I k with positive scalars dr , d! , and d . The state for the feedback controller =v _ =e v + K p e = 2 6 6 6 6 6 6 4 r ! 3 7 7 7 7 7 7 5 ; (4.17) _ = _ v = _ e v + K p _ e = 2 6 6 6 6 6 6 4 _ r _ ! _ 3 7 7 7 7 7 7 5 ; (4.18) and , e q 0 1; it follows that2 0. Note that both = 0 and =2 are equilibrium points, but due to the unity constraint on the error quaternion e q , e ~ q = 0 in both cases. 65 When the feedback controller is calculated using measurement data, we have u f =K d ~ e v + 1 2 K p ~ e + K d 2 6 6 6 6 6 6 4 0 31 ~ ~ e ~ q 0 k1 3 7 7 7 7 7 7 5 (4.19) =K d (e v +e v + K p (e +e)) + 1 2 K d 2 6 6 6 6 6 6 4 0 31 ( +) e ~ q +e ~ q 0 k1 3 7 7 7 7 7 7 5 (4.20) =K d (e v + K p e) + 1 2 K d 2 6 6 6 6 6 6 4 0 31 e ~ q 0 k1 3 7 7 7 7 7 7 5 +u f (4.21) 4.0.2 Stability of Continuous Trajectory Tracking Proposition 4.0.1. Assume that the number of the collocated actuators and sensors for vibration measurement and suppression is equal to the number of the vibrational we want to control. Let the commanded controller u c consist of three components: the inverse dy- namics termu 0 , the proportional-derivative feedback termu f , and the error compensation term u . That is u c =u 0 +u f +u = 2 6 6 6 6 6 6 4 ~ f c ~ c f c 3 7 7 7 7 7 7 5 (4.22) 66 where ~ f c denotes the vector of commanded forces;~ c denotes the commanded torque vector; f c denotes the vector of forces applied to the included modes of the exible structures. The inverse dynamics control term, which is used for tracking, is dened as u 0 = M (p) + C (p;v) _ +G(p) +u 0 : (4.23) The feedback control term u f =K d + 1 2 K d 2 6 6 6 6 6 6 4 0 31 e ~ q 0 k1 3 7 7 7 7 7 7 5 +u f : (4.24) The error compensation term, which compensates for the uncertainty due implemen- tation of the controller and measurement error, is dened as u =K n (;) (4.25) where K is a positive-denite diagonal matrix such that max (K )>kuk (4.26) where u is the system's bounded uncertainty. The function n (;) = 8 > > > > < > > > > : 0 if kk< kk if kk; (4.27) 67 and has the following properties 1:kkn (;) =kn (;)k (4.28) 2:kn (;)k> 1 kk : (4.29) (4.30) If the commanded input is dened as u c =u 0 +u f +u ; (4.31) then the solutions to the system in Equation 3.71, M (p) _ v + C (p;v)v + D _ v + K 2 6 6 4 0 61 3 7 7 5 +G (p) = B (u c +u t ) (4.32) _ q 0 = 1 2 ~ q T ~ ! (4.33) _ ~ q = 1 2 (q 0 ~ ! +~ q~ !) (4.34) are uniformly bounded. Proof. If the number of collocated actuators and sensors is equal to the number of vibra- tional modes, the input and output matrix is invertible and all states are controllable, thereby satisfying the matching conditions. Additionally, the input matrix B is invert- ible, so during stability analysis we can assume that the input matrix is the identity 68 matrix. We begin by substituting the inverse dynamics and feedback control term into the equations of motion, that is M (p) _ v + C (p; _ v)v + Dv + K 2 6 6 4 0 61 3 7 7 5 +G(p) (4.35) = M (p) + C (p;v) _ +G(p) +u 0 K d + 1 2 K d 2 6 6 6 6 6 6 4 0 31 e ~ q 0 k1 3 7 7 7 7 7 7 5 +u f +u +u t Then M (p) _ + C (p;v) + Dv + K 2 6 6 4 0 61 3 7 7 5 =K d + 1 2 K d 2 6 6 6 6 6 6 4 0 31 e ~ q 0 k1 3 7 7 7 7 7 7 5 +u +u 0 +u f +u t M (p) _ =C (p;v) Dv K 2 6 6 4 0 61 3 7 7 5 K d + 1 2 K d 2 6 6 6 6 6 6 4 0 31 e ~ q 0 k1 3 7 7 7 7 7 7 5 +u +u (4.36) 69 Recall e = e T ~ r e T ~ q e T T and e v = _ e T ~ r ~ ! T eq _ e T T . Dene the Lyapunov function candidate to be V = 1 2 T M (p) + 1 2 e T (K + K p D)e + 1 2 d! 2 (4.37) = 1 2 T M (p) + 1 2 e T ~ K + K p ~ D e + 1 2 d! 2 (4.38) = 1 2 T e T 2 6 6 6 6 6 6 4 M (p) ~ K + K p ~ D d! 3 7 7 7 7 7 7 5 2 6 6 6 6 6 6 4 e 3 7 7 7 7 7 7 5 (4.39) , 1 2 T M (4.40) where the matrixM is symmetric positive denite matrix and T = T e T , there- fore V 0. The Lyapunov function has the following lower and upper bounds 1 2 min (M)kk 2 V 1 2 max (M)kk 2 (4.41) Dierentiating the Lyapunov function gives _ V = T M(p) _ + 1 2 T _ M(p) +e T (K + K p D) _ e + d! _ : (4.42) 70 Substituting the equations of motion and _ = 1 2 e T ~ q ~ ! eq gives _ V = T 8 > > > > > > < > > > > > > : C(p;v) De v Ke K d + 1 2 K d 2 6 6 6 6 6 6 4 0 31 e ~ q 0 k1 3 7 7 7 7 7 7 5 +u +u 9 > > > > > > = > > > > > > ; + 1 2 T _ M(p) +e T (K + K p D) _ e 1 2 d! e T ~ q ~ ! eq : (4.43) Using the skew-symmetric Property 3.1.1 and substituting =e v + K p e gives _ V = e T v +e T K p T De v e T v +e T K p T Ke T K d + 1 2 e T v +e T K p T K d 2 6 6 6 6 6 6 4 0 31 e ~ q 0 k1 3 7 7 7 7 7 7 5 +e T (K + K p D)e v 1 2 d! e T ~ q ~ ! eq + T (u +u): (4.44) Note that the damping matrix D and stiness matrix K only multiply the modal am- plitudes, thereby allowing us to replace _ e with e v in the above equation. Expanding the terms gives _ V =e T v De v e T K p De v e T v Kee T K T p Ke T K d + 1 2 d! ~ ! T eq e ~ q + 1 2 d! p! e T ~ q e ~ q +e T (K + K p D)e v 1 2 d! e T ~ q ~ ! eq + T (u +u) (4.45) =e T v De v e T K T p Ke T K d + 1 2 d! p! e T ~ q e ~ q + T (u +u) (4.46) = _ e T ~ D _ e e T K T p ~ Ke T K d + 1 2 d! p! e T ~ q e ~ q + T (u +u) (4.47) 71 Note that 0 2 1, therefore _ V _ e T ~ D _ e e T K T p ~ Ke T K d 1 4 d! p! 2 e T ~ q e q + T (u +u) T 2 6 6 6 6 6 6 4 K d K p T ~ K 1 4 d! p! e T ~ q e ~ q 3 7 7 7 7 7 7 5 + T (u +u) (4.48) , T K + T (u +u) (4.49) min (K)kk 2 + T (u +u) (4.50) For the nominal system, that is when the uncertainty u = 0, the derivative of the Lyapunov equation is _ V = _ e T ~ D _ e e T K T p ~ Ke T K d + 1 2 d! p! e T ~ q e ~ q : (4.51) If _ V = 0, then e , _ e and are all zero vectors and either or e ~ q are zero. If = 0 then e q 0 = 1 implying that e ~ q is a zero vector, that is they are both simultaneously zero. Therefore, e , _ , _ , _ e q 0 , and _ e ~ q are all zero. Now we want to show that for = 0, the errors e v and e are zero. Recall = 2 6 6 6 6 6 6 4 r ! 3 7 7 7 7 7 7 5 =e v + K p e = 2 6 6 6 6 6 6 4 _ e ~ r + pr e ~ r ! eq + p! e ~ q _ e + p! e 3 7 7 7 7 7 7 5 (4.52) 72 Consider rst the system r = _ e ~ r + pr e ~ r which can be rewritten as _ e ~ r = pr e ~ r + r (4.53) If the input r tends to zero, the solution to this system will tend asymptotically to zero. Now consider ! =! eq + p! e ~ q , which can be written as e ~ q = 1 p! ! ~ ! eq (4.54) = 1 p! ! 2 e q 0 _ e ~ q _ e q 0 e ~ q e ~ q _ e ~ q : (4.55) For ! , _ e ~ q and _ e q 0 all zero, e ~ q must be zero. Therefore, the nominal closed-loop system is globally asymptotically stable. Now we use the concepts of Lyapunov Redesign to dene the error compensation control term u . First, determine the upper bound on the uncertainty term. Let E, T (u +u m ) (4.56) where u m is the error compensation term derived using the sensor measurements and is implemented as u m , n ;m = n (; m ) = 8 > > > > < > > > > : 0 if k m k< m k m k if k m k: (4.57) 73 Recall that m = +, where is determined by the sensor measurement un- certainty. Let = maxkk be the maximum absolute uncertainty in due sensor measurement inaccuracies. If the thrusters are turned o, that is u m = 0, then E = T u = ( m +) T u k m kkuk +kkkuk + ,E off (4.58) If the thrusters are on, then E, T (u +u m ) = ( m ) T (u n ;m ) = T m T u kn ;m k m k m k = T m T u T m T kn ;m k m k m k = T m T u kn ;m k k m k 2 k m k + kn ;m k k m k T m k m kkuk +kkkuk kn ;m kk m k + kn ;m k k m k kkk m k =k m kkuk +kkkuk + kn ;m kkk kn ;m kk m k k m k + + kn ;m k kn ;m kk m k = k m k (1kn ;m k) + (1 +kn ;m k) k m k kn ;m k + 2 = + 2 ,E on : (4.59) 74 We see that the disturbances entering the system are greatest when the thrusters are on. Denote the upper bound on the disturbances as E =E on . If _ V min (K)kk T +E < 0 (4.60) then the system is stable when kk E min (K) 1 2 ,: (4.61) For the ultimate bound r, max (M) min (M) (4.62) the state eventually converges to a ballB r with radiusr fort>t 0 . Therefore according to Theorem 2.5.10 on ultimate boundedness, is uniformly bounded [31]. Given that the signal , e , and are bounded, we want to show that the errors e v ande are also bounded. Note that =e q 0 1, implying thate q 0 is bounded and because the error quaternion e q is a unit vector, we have that e ~ q is also bounded. First consider r = _ e ~ r + pr e ~ r , which can be rewritten as _ e ~ r = r pr e ~ r ; (4.63) which is a linear time-invariant feedback system with plantG(s) = 1 and feedbackH(s) = pr s . The plant is strictly positive real and the feedback is positive real functions, therefore 75 the system output _ e ~ r is bounded for bounded input r . Similarly, it can also be rewritten as e ~ r = 1 p! ( r _ e ~ r ); (4.64) which is also a linear time-invariant feedback system with plantG(s) = 1 p! and feedback H(s) =s, and both are strictly positive real functions. Thus for bounded input r , the output e ~ r is bounded. Now consider ! =~ ! eq + p! e ~ q . For this system, both ! ande ~ q , therefore~ ! eq is bounded. Similarly, for = _ e + p! e , the signals ande are bounded, therefore _ e is bounded. Due to limits on the magnitude of the forces that can be applied by the thruster and actuators, the solutions are not globally uniformly bounded. The solutions are bounded only when the commanded input u c to the thrusters and actuators are less than the maximum force u max they are able to generate, therefore the u max is chosen such that ku c ku max (4.65) 76 Chapter 5 Simulation 5.1 Description of Simulation An extensive simulation was designed to test the controller and consists of the following components: I. Commanded Trajectory - The inputs to the s-function designed for generat- ing the commanded trajectory are some user-specied translational and angular velocities. The s-function then integrates a second-order model to determine the commanded positions and velocities, which along with the accelerations are inputs to the functions for each of the control terms. The commanded trajectory is dened only for the states describing the rigid body motion since vibration is undesired. II. Inverse Dynamics - The s-function for the inverse dynamics term determines the force necessary for tracking the commanded trajectory. The inputs are the commanded and measured trajectory, either of which can be used to determine the mass matrix M, Coriolis/coupling C matrix and the gravity vector G. Both the 77 commanded and measured trajectories are used to determine _ and as dened in Equation (4.12). III. Feedback Control - The feedback s-function calculates theu f as dened in Equa- tion (4.16) using the commanded and measured trajectory. IV. Error Compensation - The error compensation s-function calculatesu as dened in Equation (4.25), where the user-species and . V. Control Implementation - The commanded controller, which is a sum of the inverse dynamics, feedback and error compensation term, is used to determine the force and torques to be applied to the body for trajectory tracking and the forces needed to damp the vibration. Since the commanded controller can not be applied as desired, we need to determine which thrusters are needed and the magnitude of the thrust to be applied to track the commanded trajectory. (a) Trajectory Tracking i. Thruster Selection - The thrusters are placed around the spacecraft such that an invertible mapping can be dened from the body forces and torques determined by the commanded controller to commands for the thrusters. ii. Jet Firing Logic - The s-functionsfun_thruster_logic uses the following logic to determine the level of thrust to be applied for each thruster. if abs(cmd_thrust(i)) < 0.25*max_thrust applied_thrust(i) = 0 78 elseif (abs(cmd_thrust(i)) >= 0.25*max_thrust and abs(cmd_thrust(i)) < 0.75*max_thrust) applied_thrust(i) = sign(cmd_thrust(i))*0.5*max_thrust else applied_thrust(i) = sign(cmd_thrust(i))*max_thrust end iii. Delay and Firing Duration - A memory block is used to account for the delay in igniting or terminating the thrusters. A zero-order hold block is used to account for the minimum time required for ring a thruster. As a thruster can not be ignited or terminated instantaneously, and it must be on or o for some period time. (b) Vibration Suppression - Similarly, n damping actuators are placed in the exible structure, wheren is also the number of modes included in the modal. An invertible mapping is determined using the locations of the actuators and the mode shapes that calculates the force to be applied by the actua- tor. Additionally, the same ring logic is used to determine the force that the actuator should apply. Note that the forces applied for vibration sup- pression can apply some torque to the spacecraft which is determined by sfun_reaction_torques. VI. Equations of Motion - In the s-function sfun_scmdl the equations of motion derived in Chapter 3 are implemented. The s-function calls functions that determine the mass matrix M, the Coriolis/coupling matrix C, the damping matrix D, the 79 stiness matrix K, and the gravity vector G. The equations of motion are then integrated to generate the trajectory of the spacecraft. The measured trajectory is generated by adding some noise to the signal to account for sensor measurements errors. Figure 5.1 is of the Simulink model illustrating how the components described above are connected. 80 Figure 5.1: Matlab Simulink Model used for Simulation 81 5.2 Spacecraft and Control Parameters 5.2.1 Rigid Central Body The central body is modeled as a solid sphere with mass m b = 400 kg and radius r b = 5 meter. The inertia tensor of the central body is dened as (J b ) b = 2 6 6 6 6 6 6 4 2 5 m b r 2 b 0 0 0 2 5 m b r 2 b 0 0 0 2 5 m b r 2 b 3 7 7 7 7 7 7 5 : 5.2.2 Flexible Structures Attached to the central body are two aluminum exible structures modeled as cantilever beams. The beams are identical with one located on the positive x-axis of the central body and the other located on the negative x-axis. The dimensions for the exible structure in Figure 3.2 for n = 1; 2 L n = 3 meters b n = 0:02 meters h n = 0:02 meters: It is assumed the center of gravity (cg) located at the center of the structure. If n is the density of exible structure's material, then the mass of the nth exible structure m n = n L n b n h n . 82 5.2.3 Control Parameters The commanded controller is dened as u c =u 0 +u f +u (5.1) where the inverse dynamics term u 0 = M (p) + C (p;v) _ +G (p) (5.2) for _ as dened in Equations 4.12 with the gain matrix K p = 2 6 6 6 6 6 6 6 6 6 6 4 1:1m b I 3 1:1J b I 3 1:1m 1 I 2 1:1m 2 I 2 3 7 7 7 7 7 7 7 7 7 7 5 ; (5.3) where m is the total mass of the spacecraft and m n is the mass of the n-the exible structure for n = 1, 2. The feedback term u f =K d + K d 2 6 6 6 6 6 6 4 0 31 e ~ q 0 k1 3 7 7 7 7 7 7 5 (5.4) 83 with the gain matrix K d = 2 6 6 6 6 6 6 6 6 6 6 4 1:2m b I 3 1:2J b I 3 1:2m 1 I 2 1:2m 2 I 2 3 7 7 7 7 7 7 7 7 7 7 5 : (5.5) The error compensation term u = n (;) (5.6) where = 0:1, F zmax is the maximum thrust applicable in the z-direction, = 0:275 6F zmax to account for a 25% error, at most, in implementation in the six thruster. 5.3 Simulation Results 5.3.1 Vertical Trajectory The objective is to track a vertical trajectory using the pulse-modulated thrusters. The commanded acceleration in the z-direction is denoted z c , which is integrated with zero initial conditions to determine the commanded positions and velocities. To get a better view of the impact the thrusters have on the exible structures, we have assumed there is no measurement uncertainty. In the following plots we'll show that both error com- pensation and vibration suppression is needed during trajectory tracking. The following table includes additional parameter specications for the simulation. 84 z c 1 m =sec 2 Maximum Thrust in x 50 N Maximum Thrust in y 50 N Maximum Thrust in z 900 N Thrust Duration 0.25 s Maximum Thrust applied by Actuators 10 N Actuator Duration 0.01 s Table 5.1: Parameters of Vertical Trajectory In the following plots, we see that the error compensation term is necessary as is the vibration suppression itself. Note that the subscript c denotes the term describes the commanded trajectory. Figure 5.2 is of the commanded and actual displacement along the z-axis. Figure 5.3 shows the displacement without the error compensation term. In Figures 5.4 and 5.5, we have the commanded and actual velocity in the z-direction with and without error compensation, respectively. In Figures 5.6 and 5.7, we have the com- manded and actual rotation angle about the z-direction, where the commanded rotation angle is zero, with and without error compensation, respectively. Figures 5.8-5.10 show the end de ection in one of the exible structures; it's evident that not only in the vibra- tion suppression necessary, but the error compensation term makes a large contribution in determining the damping forces. Lastly, Figure 5.11 depicts the commanded and applied force in the z-direction. 85 0 1 2 3 4 5 6 7 8 9 10 −10 0 10 20 30 40 50 time (s) z−displacement (m) z z c Figure 5.2: Commanded and Actual Displacement in Z-Direction 86 0 1 2 3 4 5 6 7 8 9 10 −10 0 10 20 30 40 50 60 time (s) z−displacement (m) z z c Figure 5.3: Commanded and Actual Velocity in Z-Direction Without Error Compensation 87 0 1 2 3 4 5 6 7 8 9 10 −2 0 2 4 6 8 10 time (s) velocity in z−direction (m/s) dz dz c Figure 5.4: Commanded and Actual Velocity in Z-Direction 88 0 1 2 3 4 5 6 7 8 9 10 −2 0 2 4 6 8 10 12 time (s) velocity in z−direction (m/s) dz dz c Figure 5.5: Commanded and Actual Displacement in Z-Direction Without Error Com- pensation 89 0 1 2 3 4 5 6 7 8 9 10 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 time (s) Rotation angle φ (deg) φ φ c Figure 5.6: Commanded and Actual Rotation Angle of Spacecraft 90 0 1 2 3 4 5 6 7 8 9 10 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 time (s) Rotation angle φ (deg) φ φ c Figure 5.7: Commanded and Actual Rotation Angle of Spacecraft Without Error Com- pensation 91 0 1 2 3 4 5 6 7 8 9 10 −60 −50 −40 −30 −20 −10 0 10 20 30 40 time (s) end deflection (mm) y−direction z−direction Figure 5.8: End De ection of Flexible Structure Position on Central Body's x-axis 92 0 1 2 3 4 5 6 7 8 9 10 −100 −80 −60 −40 −20 0 20 40 60 80 100 time (s) end deflection (mm) y−direction z−direction Figure 5.9: End De ection of Flexible Structure Position on Central Body's x-axis With- out Error Compensation 93 0 1 2 3 4 5 6 7 8 9 10 −100 −80 −60 −40 −20 0 20 40 60 80 100 time (s) end deflection (mm) y−direction z−direction Figure 5.10: End De ection of Flexible Structure Position on Central Body's x-axis Without Vibration Suppression 94 0 1 2 3 4 5 6 7 8 9 10 0 1000 2000 3000 4000 5000 6000 7000 8000 time (s) commanded applied Figure 5.11: Commanded and Applied Forces in Z-Directions 95 5.3.2 Vertical Trajectory with Angular Acceleration Here, the objective is to track a vertical trajectory while maintaining a constant angu- lar acceleration using the pulse-modulated thrusters. The commanded acceleration z c is integrated with zero initial conditions to determine the commanded position and veloc- ity. The commanded angular acceleration denoted _ ! c is also integrated with zero initial conditions to determine the commanded angular velocity and quaternion. Some uncer- tainty has been added to the measurements to account for sensor measurement error. The following table includes additional parameter specications for the simulation. z c 1 m =sec 2 _ ! c 0:1 rad =sec 2 Maximum Thrust in x 100 N Maximum Thrust in y 100 N Maximum Thrust in z 900 N Thrust Duration 0.25 s Maximum Thrust applied by Actuators 10 N Actuator Duration 0.01 s Table 5.2: Parameters of Vertical Trajectory with Constant Angular Acceleration Figure 5.12 is of the commanded and actual displacement along thez-axis. Figure 5.13 shows the displacement without the error compensation term. In Figures 5.14 and 5.15, we have the commanded and actual velocity in the z-direction with and without error compensation. In Figures 5.14 and 5.15, we have the commanded and actual rotation 96 angle in the z-direction with and without error compensation. Figures 5.18-5.20 show the end de ection in one of the exible structures; again it's evident that not only in the vibration suppression necessary, but the error compensation term has a large impact on determining the damping forces. Lastly, Figure 5.21 depicts the commanded and applied force in the z-direction. 97 0 1 2 3 4 5 6 7 8 9 10 −10 0 10 20 30 40 50 time (s) z−displacement (m) z z c Figure 5.12: Commanded and Actual Displacement in Z-Direction 98 0 1 2 3 4 5 6 7 8 9 10 −10 0 10 20 30 40 50 60 time (s) z−displacement (m) z z c Figure 5.13: Commanded and Actual Velocity in Z-Direction Without Error Compensa- tion 99 0 1 2 3 4 5 6 7 8 9 10 −2 0 2 4 6 8 10 12 time (s) velocity in z−direction (m/s) dz dz c Figure 5.14: Commanded and Actual Velocity in Z-Direction 100 0 1 2 3 4 5 6 7 8 9 10 −2 0 2 4 6 8 10 12 time (s) velocity in z−direction (m/s) dz dz c Figure 5.15: Commanded and Actual Displacement in Z-Direction Without Error Com- pensation 101 0 1 2 3 4 5 6 7 8 9 10 0 50 100 150 200 250 300 time (s) Rotation angle φ (deg) φ φ c Figure 5.16: Commanded and Actual Rotation Angle of Spacecraft 102 0 1 2 3 4 5 6 7 8 9 10 0 50 100 150 200 250 300 time (s) Rotation angle φ (deg) φ φ c Figure 5.17: Commanded and Actual Rotation Angle of Spacecraft Without Error Com- pensation 103 0 1 2 3 4 5 6 7 8 9 10 −60 −50 −40 −30 −20 −10 0 10 20 30 40 time (s) end deflection (mm) y−direction z−direction Figure 5.18: End De ection of Flexible Structure Position on Central Body's x-axis 104 0 1 2 3 4 5 6 7 8 9 10 −100 −80 −60 −40 −20 0 20 40 60 80 100 time (s) end deflection (mm) y−direction z−direction Figure 5.19: End De ection of Flexible Structure Position on Central Body's x-axis Without Error Compensation 105 0 1 2 3 4 5 6 7 8 9 10 −100 −50 0 50 100 150 time (s) end deflection (mm) y−direction z−direction Figure 5.20: End De ection of Flexible Structure Position on Central Body's x-axis Without Vibration Suppression 106 0 1 2 3 4 5 6 7 8 9 10 0 1000 2000 3000 4000 5000 6000 7000 8000 time (s) commanded applied Figure 5.21: Commanded and Applied Forces in Z-Directions 107 5.3.3 Diagonal Trajectory with Commanded Angular Acceleration about all Axes The objective is to track a trajectory with constant acceleration along the x, y, and z axes and constant angular acceleration about the x, y, and z axes. The accelerations are integrated with zero initial conditions to determine the corresponding commanded positions and velocities. The variables describing the commanded trajectory are denoted with the subscript c. Some uncertainty has again been added to the measurements to account for sensor measurement error. In the following plots we'll show that both error compensation and vibration suppression is needed during trajectory tracking. The following table includes additional parameter specications for the simulation. x c 1 m =sec 2 y c 1 m =sec 2 z c 1 m =sec 2 _ ! xc 0:1 rad =sec 2 _ ! yc 0:1 rad =sec 2 _ ! zc 0:1 rad =sec 2 Maximum Thrust in x 100 N Maximum Thrust in y 100 N Maximum Thrust in z 900 N Thrust Duration 0.5 s Maximum Thrust applied by Actuators 10 N Actuator Duration 0.01 s Table 5.3: Parameters of Diagonal Trajectory with Constant Angular Accelerations about all Axes Figures 5.24 - 5.24 are of the displacement in the x, y, and z-direction with error compensation. Figures 5.25 - 5.27 show the commanded and actual translational velocities with error compensation. Figure 5.28 is of the rotation angle. Figure 5.29 shows the end de ection in one of the exible structures. Figures 5.31 - 5.33 depicts the body forces 108 commanded by the controller and resultant force applied by the thrusters. Figures 5.34 - 5.36, similarly, show the torques. 109 0 1 2 3 4 5 6 7 8 9 10 0 5 10 15 20 25 30 35 40 45 50 time (s) x−displacement (m) x x c Figure 5.22: Commanded and Actual Displacement in X-Direction 110 0 1 2 3 4 5 6 7 8 9 10 0 5 10 15 20 25 30 35 40 45 50 time (s) y−displacement (m) y y c Figure 5.23: Commanded and Actual Displacement in Y-Direction 111 0 1 2 3 4 5 6 7 8 9 10 −10 0 10 20 30 40 50 time (s) z−displacement (m) z z c Figure 5.24: Commanded and Actual Displacement in Z-Direction 112 0 1 2 3 4 5 6 7 8 9 10 0 1 2 3 4 5 6 7 8 9 10 time (s) velocity in x−direction (m/s) dx dx c Figure 5.25: Commanded and Actual Velocity in X-Direction 113 0 1 2 3 4 5 6 7 8 9 10 0 1 2 3 4 5 6 7 8 9 10 time (s) velocity in y−direction (m/s) dy dy c Figure 5.26: Commanded and Actual Velocity in Y-Direction 114 0 1 2 3 4 5 6 7 8 9 10 −2 0 2 4 6 8 10 12 time (s) velocity in z−direction (m/s) dz dz c Figure 5.27: Commanded and Actual Velocity in Z-Direction 115 0 1 2 3 4 5 6 7 8 9 10 0 50 100 150 200 250 300 350 400 time (s) Rotation angle φ (deg) φ φ c Figure 5.28: Commanded and Actual Rotation Angle of Spacecraft 116 0 1 2 3 4 5 6 7 8 9 10 −60 −40 −20 0 20 40 60 80 time (s) end deflection (mm) y−direction z−direction Figure 5.29: End De ection of Flexible Structure Position on Central Body's x-axis 117 0 1 2 3 4 5 6 7 8 9 10 −300 −200 −100 0 100 200 300 time (s) end deflection (mm) y−direction z−direction Figure 5.30: End De ection of Flexible Structure Position on Central Body's x-axis Without Vibration Suppression 118 0 1 2 3 4 5 6 7 8 9 10 −1000 −500 0 500 1000 1500 2000 time (s) commanded applied Figure 5.31: Commanded and Applied Forces in X-Directions 119 0 1 2 3 4 5 6 7 8 9 10 −1000 −500 0 500 1000 1500 2000 2500 3000 time (s) commanded applied Figure 5.32: Commanded and Applied Forces in Y-Directions 120 0 1 2 3 4 5 6 7 8 9 10 0 1000 2000 3000 4000 5000 6000 7000 8000 time (s) commanded applied Figure 5.33: Commanded and Applied Forces in Z-Directions 121 0 1 2 3 4 5 6 7 8 9 10 −4000 −3000 −2000 −1000 0 1000 2000 3000 4000 5000 time (s) commanded applied Figure 5.34: Commanded and Applied Torques about X-Axis 122 0 1 2 3 4 5 6 7 8 9 10 −6000 −4000 −2000 0 2000 4000 6000 8000 10000 time (s) commanded applied Figure 5.35: Commanded and Applied Torques in Y-Axis 123 0 1 2 3 4 5 6 7 8 9 10 −2000 0 2000 4000 6000 8000 10000 time (s) commanded applied Figure 5.36: Commanded and Applied Torques in Z-Axiz 124 5.3.4 Vertical Trajectory with Resting End Conditions The objective is to track a vertical trajectory then stop at some some specied location z f after time t f . The accelerations are integrated with zero initial conditions to deter- mine the corresponding commanded positions and velocities. The variables describing the commanded trajectory are denoted with the subscript c. In the following plots we'll show that both error compensation and vibration suppression is needed during trajec- tory tracking. The following table includes additional parameter specications for the simulation. z c 2 2 + 6 3 t + 12 4 t 2 Maximum Thrust in x 50 N Maximum Thrust in y 50 N Maximum Thrust in z 900 N Thrust Duration 0.5 s Maximum Thrust applied by Actuators 10 N Actuator Duration 0.01 s Table 5.4: Parameters of Vertical Trajectory with End Conditions 125 The coecients 2 , 3 , 4 of the polynomial for the angular acceleration are deter- mined using zero initial conditions, specifying some nal position z f at time t f and zero for the nal velocity and acceleration. Figure 5.37 is of the commanded and actual displacement along the z-axis. Figure 5.38 shows the displacement without the error compensation term. In Figures 5.39 and 5.40, we have the commanded and actual velocity in the z-direction. In Figures 5.39 and 5.40, we have the commanded and actual rotation angle in the z-direction. In Figures 5.41 - 5.43 show the end de ection in one of the exible structures; it's evident that not only in the vibration suppression necessary, but the error compensation term has a large impact on determining the damping forces. Lastly, Figure 5.44 depicts the commanded and applied force in the z-direction. 126 0 1 2 3 4 5 6 7 8 9 10 −2 0 2 4 6 8 10 12 time (s) z−displacement (m) z z c Figure 5.37: Commanded and Actual Displacement in Z-Direction 127 0 1 2 3 4 5 6 7 8 9 10 −2 0 2 4 6 8 10 12 time (s) z−displacement (m) z z c Figure 5.38: Commanded and Actual Velocity in Z-Direction Without Error Compensa- tion 128 0 1 2 3 4 5 6 7 8 9 10 −2 −1 0 1 2 3 4 time (s) velocity in z−direction (m/s) dz dz c Figure 5.39: Commanded and Actual Velocity in Z-Direction 129 0 1 2 3 4 5 6 7 8 9 10 −2 −1 0 1 2 3 4 time (s) velocity in z−direction (m/s) dz dz c Figure 5.40: Commanded and Actual Displacement in Z-Direction Without Error Com- pensation 130 0 1 2 3 4 5 6 7 8 9 10 −60 −40 −20 0 20 40 60 time (s) end deflection (mm) y−direction z−direction Figure 5.41: End De ection of Flexible Structure Position on Central Body's x-axis 131 0 1 2 3 4 5 6 7 8 9 10 −80 −60 −40 −20 0 20 40 60 80 100 120 time (s) end deflection (mm) y−direction z−direction Figure 5.42: End De ection of Flexible Structure Position on Central Body's x-axis Without Error Compensation 132 0 1 2 3 4 5 6 7 8 9 10 −150 −100 −50 0 50 100 150 200 250 time (s) end deflection (mm) y−direction z−direction Figure 5.43: End De ection of Flexible Structure Position on Central Body's x-axis Without Vibration Suppression 133 0 1 2 3 4 5 6 7 8 9 10 0 1000 2000 3000 4000 5000 6000 7000 time (s) commanded applied Figure 5.44: Commanded and Applied Forces in Z-Directions 134 Chapter 6 Concluding Remarks and Future Work In the preceding sections, a nonlinear model has been derived that illustrates the coupling that exists between the translational, angular, and modal degrees of freedom. A controller is then designed for this model for guidance and attitude control in a three dimensional reference plane and for vibration suppression in the exible structures. Using passivity and Lyapunov theory concepts, we have proven that when implementing this controller the trajectory tracking error is uniformly bounded in the presence of disturbances. In the analysis, we have accounted for disturbances due to on-o thruster implementation, sensor measurements, and parameter uncertainties. Further investigation is needed in the following areas. 1. Generalization of the Model (a) Additional Degrees of Freedom for Vibration The exible structures may also be susceptible to longitudinal and torsional vibration in addition to the lateral vibration. 135 Figure 6.1: Directions of Vibration in Three Dimensions (b) Flexible Central Body Another area of investigation is to design a controller for spacecraft with a central body modeled as a exible structure making the entire spacecraft ex- ible. In industry, one of the main concern in space missions is to minimize the mass, which, as mentioned earlier, can cause exibility. The Falcon 9 rocket at SpaceX was susceptible to longitudinal vibration during ascent, thus such a controller could be applied to ascent trajectory for a wide range of missions. Furthermore, such a controller could also be applied to entry, landing and descent of small spacecraft like Mars Phoenix on other planets. 2. Generalization of proof for three dimensional motion (a) Global Asymptotic Stability Erdong et al. proved global asymptotic stability during attitude control using passivity concepts of a exible spacecraft with disturbances to the angular degrees of freedom [15]. With actuators on the exible structures for damping vibration, it may be possible to extend the passivity concepts to prove global 136 asymptotic stability during tracking. Recall that attitude tracking can excite vibration in the exible structures. We will investigate if similar concepts could be applied to prove global asymptotic stability to the 3D model. (b) Damping Actuators During attitude tracking, any uncertainties in the model or external distur- bances can disturb the exible structures due to coupling, making it necessary to add actuators along the exible structures to compensate the disturbances. We assume that the number of actuators is equal to the number of modes we want to control in order to satisfy the matching conditions, and also because the commanded controller is dependent on the modal amplitudes thereby re- quiring that the input matrix be invertible. We want to reduce the number of actuators needed to control the modes that are excited, that is we want to nd a way to satisfy the matching conditions but minimize the number of actuators necessary for damping the vibration. 3. Realistic Simulation The spacecraft model and controller in the simulation will be updated to study a realistic spacecraft model and mission such as the Mars 98. The simulation will also include hardware inaccuracies, delays in the jet rings and termination, and a realistic model of external disturbances such as wind gusts. 137 References [1] Brij Agrawal, R.S. McClelland, and Gangbing Song. Attitude control of exible of spacecraft using pulse-width pulse frequency modulated thrusters. Space Technology, 17(1):15{34, 1997. [2] Hyochoong Bang. Manuver and Vibration Control of Flexible Structures by Lyapunov Stability Theory. PhD thesis, Texas A&M University, August 1992. [3] Hyochoong Bang, Cheol-Keun Ha, and Jin Hyoung Kim. Flexible spacecraft attitude manever by application of sliding mode control. Acta Astronautica, 57:841{850, 2005. [4] C. Benatzky, M. Kozek, and C. Bilik. Experimental control of a exible beam using a stack-bending actuator principle. In 20th Scientic Conference, 2006. [5] R.J. Benhabib, R.P. Iwens, and R.L. Jackson. Stability of large space structure con- trol systems using positivity concepts. Journal of Guidance and Control, 4(5):487{ 494, Sept.-Oct. 1981. [6] Dennis S. Bernstein. Matrix Mathematics. Princeton University Press, 2005. [7] Jovan D. Boskovic, Sai-Ming Li, and Raman Mehra. Robust tracking control design for spacecraft under control input saturation. Journal of Guidance, Control and Dynamics, 27(4):627{633, July-August 2004. 138 [8] Jonathan M. Cameron and Wayne J. Book. Modeling mechanisms with nonholo- nomic joints using the boltzmann-hamel equations. International Journal of Robotics Research, 16(1):47{59, February 1997. [9] Martin Corless. Guaranteed rates of exponential convergence for uncertain systems. Journal of Optimization Theory and Applications, 64:481{493, 1990. [10] Fabio Curti, Marcello Romano, and Riccardo Bevilacqua. Lyapunov-based thrusters' selection for spacecraft control: Analysis and experimentation. Journal of Guidance, Control and Dynamics, 33(4):1143{1160, July-August 2010. [11] Carlos Canudas de Wit. Theory of Robot Control. Springer-Verlag London Limited, 1996. [12] S. DiGennaro. Active vibration suppression in exible spacecraft attitude tracking. Journal of Guidance, Control and Dynamics, 21(3):400{408, May-June 1998. [13] S. DiGennaro. Output attitude control of exible spacecraft from quaternion mea- sures: a passivity approach. In Proceedings of the 37th IEEE Conference on Decision and Control, 1998. [14] S. DiGennaro. Passive attitude control of exible spacecraft from quaternion mea- surements. Journal of Optimization Theory and Applications, 116(1):41{60, January 2003. [15] Jin Erdong and Sun Zhaowei. Passivity-based control for a exible spacecraft in the presence of disturbances. International Journal of Non-Linear Mechanics, 45:348{ 356, 2010. 139 [16] D. M. Dawson Frank L. Lewis, Chaouki T. Abdallah. Control of Robot Manipulators. Prentice-Hall, Inc., 1993. [17] Shekhar Gosavi and A.G. Kelkar. Passivty-based robust control of piezo-actuated exible beam. In Proceedings of the American Control Conference, pages 2492{2497, Arlington, VA, June 2001. [18] Alexander Graham. Kronecker Products and Matrix Calculus. Ellis Horwood Lim- ited, 1981. [19] Qinglei Hu. Sliding mode maneuvering control and active vibration damping of three- axis stabilized exible spacecraft with actuator dynamics. Nonlinear Dynamics, 52:227{248, 2008. [20] Qinglei Hu and Guangfu Ma. Vibration suppression of exible spacecraft during attitude maneuvers. Journal of Guidance, Control and Dynamics, 28(2):311{380, March-April 2005. [21] Yong-An Huang, Zi-Chen Deng, and You-Iun Xiong. High-order model and slide mode control for rotating exible smart structure. Mechanism and Machine Theory, 43:1038{1054, 2008. [22] Ashok Iyer and Sahjendra N. Singh. Sliding mode of control of exible spacecraft under disturbance torque. International Journal of Systems Science, 21(9):1755{ 1771, 1990. 140 [23] Wu-Xing Jing, Shi-Jie, and Xu. Nonlinear attitude tracking control of a spacecraft with thrusters based on error quaternions. Chinese Journal of Aeronautics, 15:129{ 138, 2002. [24] USAF Joseph M. Fulton, Major. LQG/LTR Optimal Attitude Control of Small Flexible Spacecraft Using Free-Free Boundary Conditions. PhD thesis, University of Colorado, 2006. [25] S.M. Joshi, A.G. Kelkar, and J. T.-Y. Wen. Robust attitude stabilization of space- craft using nonlinear quaternion feedback. IEEE Transactions on Automatic Control, 40(10):1800{1803, October 1995. [26] S.M. Joshi, P.G. Maghami, and A.G. Kelkar. Design of dynamic dissipative compen- sators for exible space structures. IEEE Transactions on Aerospace and Electronic Systems, 31(5):1314{1324, October 1995. [27] John L. Junkins and James D. Turner. Optimal Spacecraft Rotational Maneuvers, volume 3 of Studies in Astronautics. Elsevier Science, 1986. [28] Atul Kelkar, Thomas E. Alberts, and Suresh Joshi. Dynamic dissipative compen- sators for multibody exible space structures. IEEE Transaction on Aerospace and Electronic Systems, 31(4):1325{1330, October 1995. [29] Atul Kelkar and Suresh Joshi. Control of nonlinear multibody exible space struc- tures. In Lecture Notes in Control and Information Sciences. Springer, 1996. 141 [30] Atul Kelkar, Suresh Joshi, and Thomas Alberts. Passivity-based control of nonlinear exible multibody systems. IEEE Transactions on Automatic Control, 40(5):910{ 914, May 1995. [31] Hassan K. Khalil. Nonlinear Systems. Prentice Hall, third edition, 2002. [32] Sylvia Kohn-Rich. Robust Fuzzy Logic Control of Mechanical Systems. PhD thesis, University of Southern California, May 2000. [33] Jack B. Kuipers. Quaternions and Rotation Sequences. Princeton University Press, 1999. [34] Jenchieh Lee. An Approach Experimentally Based Modeling and Simulation of Hu- man Motion. PhD thesis, University of Southern California, 2009. [35] Antonio Loria, Elena Panteley, and Henk Nijmeijer. A remark on passivity-based and discontinous control of uncertain nonlinear systems. Automatica, 37(9):1481{1487, 2001. [36] Leonard Meirovitch. Analytical Methods in Vibrations. The Macmillan Company, 1967. [37] Leonard Meirovitch. Methods of Analytical Dynamics. Dover Publications, Inc., 1970. [38] Ju. I. Neimark and N.A. Fufaev. Dynamics of Nonholonomic Systems. American Mathematical Society, 1972. 142 [39] Parviz E. Nikravesh. Computer-Aided Analysis of Mechanical Systems. Prentice Hall, 1988. [40] William J O'Connor. Control of exible mechanical systems: Wave-based techniques. In Proceedings of the 2007 American Control Conference, 2007. [41] Brad Paden, Brad Riedle, and Eduardo Bayo. Exponentially stable tracking con- trol for multi-joint exible link manipulators. In Proceedings of the 1990 American Control Confrence, pages 680{684, 1990. [42] Jianjun Shi. Control of Nonlinear Flexible Space Structures. PhD thesis, Iowa State University, 2005. [43] Malcolm D. Shuster. A survey of attitude representations. Journal of the Astronau- tical Sciences, 41(4):439{517, October-December 1993. [44] Sahjendra N. Singh and Ashok Iyer. Nonlinear decoupling sliding mode control and attitude control of spacecraft. IEEE Transaction on Aerospace and Electronic Systems, 25(5):621{633, 1989. [45] Hebertt Sira-Ramirez. The dierential algebraic approach in nonlinear dynamical feedback controlled landing maneuvers. IEEE Transactions on Automatic Control, 37:518{523, 1992. [46] Robert E. Skelton, Peter C. Hughes, and Hari B. Hablani. Order reduction for models of space structures using modal cost analysis. Journal of Guidance, 5(4):351{357, July-August 1982. 143 [47] Gangbing Song and Brij Agrawal. Vibration suppression of exible spacecraft during attitude control. Acta Astronautica, 49(2):73{83, 2001. [48] Yoon-Gyeoung Sung. Model Reduction and Robust Vibration Control of Flexible Space Structures. PhD thesis, The University of Alabama, 1997. [49] Sam W. Thurman. A Comprehensive Guidance and Control Theory for Spacecraft Using Pulse-Modulated Propulsion. PhD thesis, University of Southern California, 1995. [50] Sam W. Thurman and Henryk Flashner. New pulse-modulation technique for guid- ance and control of automated spacecraft. Journal of Guidance, Control and Dy- namics, 19:1007{1016, 1996. [51] Sam W. Thurman and Henryk Flashner. Robust digital autopilot design for space- craft equipped with pulse-operated thrusters. Journal of Guidance, Control and Dynamics, 19(5):1047{1055, September-October 1996. [52] James D. Turner and John L. Junkins. Optimal large-angle single-axis rotational maneuvers of exible spacecraft. Journal of Guidance and Control, 3(6):578{585, November-December 1980. [53] John Ting-Yung Wen and Kenneth Kreutz-Delgado. The attitude control problem. IEEE Transactions on Automatic Control, 36(10):1148{1161, October 1991. [54] Bong Wie and Peter M. Barba. Quaternion feedback for spacecraft large angle maneuvers. Journal of Guidance, 8(3):360{365, 1985. 144 [55] Bong Wie, H. Weiss, and A. Arapostathis. Quaternion feedback regulator for space- craft eigenaxis rotations. Journal of Guidance, 12:375{380, 1989. 145 Appendix A Quaternion Quaternions are an alternative representation of the attitude in the equations of motion of a spacecraft. They allow us to avoid singularities that can arise when using other representations like the Euler angles. The following overview of quaternions is based on Kuipers' Quaternions and Rotation Sequences [33] and Nikravesh's Computer-Aided Analysis of Mechanical Systems [39]. A.1 Quaternion Denition A quaternion q2R 4 can be expressed as q = 2 6 6 4 q 0 ~ q 3 7 7 5 (A.1) 146 where q 0 is the scalar part of the quaternion and ~ q = 2 6 6 6 6 6 6 4 q 1 q 2 q 3 3 7 7 7 7 7 7 5 2R 3 (A.2) is the vector part of the quaternion. A pure quaternion has a zero scalar part; therefore any vector ~ w2R 3 can be expressed as a pure quaternion as follows w = 2 6 6 4 0 ~ w 3 7 7 5 : (A.3) A.2 Quaternion Algebra Two quaternions are equal if and only if the have exactly the same components, that is if p = 2 6 6 6 6 6 6 6 6 6 6 4 p0 p1 p2 p3 3 7 7 7 7 7 7 7 7 7 7 5 and q = 2 6 6 6 6 6 6 6 6 6 6 4 q0 q1 q2 q3 3 7 7 7 7 7 7 7 7 7 7 5 (A.4) 147 then p =q if and only if p 0 =q 0 p 1 =q 1 p 2 =q 2 p 3 =q 3 : The sum of two quaternions p and q is dened by adding the components, as follows p +q = 2 6 6 6 6 6 6 6 6 6 6 4 p 0 +q 0 p 1 +q 1 p 2 +q 2 p 3 +q 3 3 7 7 7 7 7 7 7 7 7 7 5 : (A.5) The quaternion product, denoted by, of p and q is dened as pq = 2 6 6 4 p 0 ~ p T ~ p p 0 I 3 +J~ pK 3 7 7 5 2 6 6 4 q 0 ~ q 3 7 7 5 (A.6) = 2 6 6 4 p 0 q 0 ~ p T ~ q q 0 ~ p + (p 0 I 3 +J~ pK)~ q 3 7 7 5 (A.7) = 2 6 6 6 6 6 6 6 6 6 6 4 p 0 p 1 p 2 p 3 p 1 p 0 p 3 p 2 p 2 p 3 p 0 p 1 p 3 p 2 p 1 p 0 3 7 7 7 7 7 7 7 7 7 7 5 2 6 6 6 6 6 6 6 6 6 6 4 q 0 q 1 q 2 q 3 3 7 7 7 7 7 7 7 7 7 7 5 (A.8) 148 The complex conjugate of the quaternion is denoted q = 2 6 6 4 q 0 ~ q 3 7 7 5 : (A.9) The norm of a quaternion can then be dened as jqj = p q q: (A.10) The inverse of the quaternion is dened in terms of the norm and complex conjugate, that is q 1 = q jqj 2 (A.11) A.3 Unit Quaternion Recall from Chapter 2, that a unit quaternion of the form q = 2 6 6 4 q 0 ~ q 3 7 7 5 = 2 6 6 4 cos 2 ~ e sin 2 3 7 7 5 = 2 6 6 6 6 6 6 6 6 6 6 4 q 0 q 1 q 2 q 3 3 7 7 7 7 7 7 7 7 7 7 5 : (A.12) 149 is used for attitude representation. Note that ~ q; ~ e2 R 3 andj~ ej = 1. Additionally, the complex conjugate is equal to the inverse of the quaternion, that is q 1 =q = 2 6 6 4 q 0 ~ q 3 7 7 5 : (A.13) A.4 Quaternion Rotation Operator Let the unit quaternionq represent the rotation from the inertial coordinate frame to the body coordinate frame. The rotation of a vector ~ w2R 3 in the inertial frame to the body frame can be expressed as a product of quaternions, i.e. qwq = 2 6 6 4 1 T I b (q) 3 7 7 5 0~ w (A.14) where the transformation matrix T I b (q) = q 2 0 j~ qj 2 I 3 + 2~ q~ q T + 2q 0 J~ qK: (A.15) Henceforth, T I b will be denoted as T. The transformation matrix T can also be written as T = GL T (A.16) 150 where G = 2 6 6 6 6 6 6 4 q 1 q 0 q 3 q 2 q 2 q 3 q 0 q 1 q 3 q 2 q 1 q 0 3 7 7 7 7 7 7 5 = ~ q; q 0 I 3 +J~ qK (A.17) L = 2 6 6 6 6 6 6 4 q 1 q 0 q 3 q 2 q 2 q 3 q 0 q 1 q 3 q 2 q 1 q 0 3 7 7 7 7 7 7 5 = ~ q; q 0 I 3 J~ qK (A.18) These matrices are subject to the following properties i. The rows of G and L are orthogonal to q Gq = ~ q; q 0 I 3 +J~ qK q = 0 (A.19) and Lq = ~ q; q 0 I 3 J~ qK q = 0 (A.20) ii. The rows of G and L are orthogonal GG T = I 3 (A.21) and LL T = I 3 (A.22) 151 iii. Dierentiating Equations A.19 and A.20 yield G _ q = _ Gq (A.23) and L _ q = _ Lq (A.24) iv. The rows of _ G and _ L are orthogonal to _ q _ G _ q = _ ~ q; _ q 0 I 3 + r _ ~ q z _ q = 0 (A.25) and _ L _ q = _ ~ q; _ q 0 I 3 r _ ~ q z _ q = 0 (A.26) v. From Equations A.19 and A.20, G _ L T = _ GL T (A.27) vi. The derivative of the transformation matrix _ T = d dt GL T = _ GL T + G _ L T = 2 _ GL T = 2G _ L (A.28) 152 vii. The product G _ q can be expanded as G _ q = ~ q; q 0 I 3 +J~ qK _ q = _ ~ qq 0 ~ q _ q 0 +J~ qK _ ~ q (A.29) Transformation of G _ q into a skew-symmetric matrix yields J G _ qK G _ G T (A.30) Similarly, J L_ qK = L _ L T (A.31) viii. Additionally, the following identity can be derived using Equations A.23 and A.30 G _ G T = _ GG T (A.32) Similarly, using Equations Equations A.24 and A.31, we have L _ L T = _ LL T (A.33) 153 A.5 Identities with Arbitrary Vectors Let~ r be an arbitrary vector, G T ~ r = 2 6 6 4 ~ q T J~ qK +q 0 I 3 7 7 5 ~ r = 2 6 6 4 ~ q T ~ r J~ qK~ r +q 0 ~ r 3 7 7 5 = 2 6 6 4 ~ r T ~ q +J~ rK~ q +q 0 ~ r 3 7 7 5 = 2 6 6 4 0 ~ r T ~ r J~ rK 3 7 7 5 2 6 6 4 q 0 ~ q 3 7 7 5 = + rq: (A.34) Similarly, L T ~ r = rq (A.35) where r = 2 6 6 4 0 ~ r T ~ r J~ rK 3 7 7 5 : (A.36) A.6 Quaternion Time Derivative In Chapter 2, the derivative of the quaternion is dened as _ q = 1 2 q!; (A.37) 154 which can be written as ! = 2L _ q =2 _ Lq: (A.38) 155 Appendix B Kronecker Product and Matrix Derivative B.1 Kronecker Product For matrices A2R nm and B2R pq , the Kronecker product is dened as A B, a ij B 2R npmq : (B.1) If the matrix A is a function of the state vector p, i.e. A(p), then the matrix derivative is dened as @A(p) @p = 2 6 6 6 6 6 6 6 6 6 6 6 4 @A(p) @p 1 @A(p) @p 2 . . . @A(p) @p 1 3 7 7 7 7 7 7 7 7 7 7 7 5 2R snm : (B.2) 156 B.2 Matrix Derivative The derivative of matrix with respect to a vector is a tensor of order 3. With this denition of the matrix derivative, the product rule is @ @p [A(p)B(p)] = (I s A) @B @p + @A @p B (B.3) The following properties are a result of the Kronecker Product Property B.2.1. (I s _ p) _ p = ( _ p I n ) _ p = 2 6 6 6 6 6 6 6 6 6 6 4 _ p 1 _ p _ p 2 _ p . . . _ p s _ p 3 7 7 7 7 7 7 7 7 7 7 5 ,h _ p _ pi Property B.2.2. For any 2R n and 2R m ( I m ) = (I n ) (B.4) Property B.2.3. (A B) (C D) = (AC) (BD) for A2 R mn ; B2 R pq ; C2 R nk ; D2R qr Property B.2.4. For z2R n and A2R pq , the following properties hold (i). (I p z) A = A z (ii). A I q z T = A z T (iii). (z I p ) A =z A (iv). A z T I q =z T A 157 Property B.2.5. Let = 2 6 6 6 6 6 6 6 6 6 6 4 1 2 . . . m 3 7 7 7 7 7 7 7 7 7 7 5 where2R mn1 and i 2R n1 8i. Let A = 2 6 6 6 6 6 6 6 6 6 6 4 T 1 T 2 . . . T m 3 7 7 7 7 7 7 7 7 7 7 5 2 R mn . Then for any u2R n1 ; I m u T = Au Proof. I m u T = 2 6 6 6 6 6 6 6 6 6 6 4 u T u T . . . u T 3 7 7 7 7 7 7 7 7 7 7 5 2 6 6 6 6 6 6 6 6 6 6 4 1 2 . . . m 3 7 7 7 7 7 7 7 7 7 7 5 = 2 6 6 6 6 6 6 6 6 6 6 4 u T 1 u T 2 . . . u T m 3 7 7 7 7 7 7 7 7 7 7 5 = 2 6 6 6 6 6 6 6 6 6 6 4 T 1 u T 2 u . . . T m u 3 7 7 7 7 7 7 7 7 7 7 5 = Au Property B.2.6. For any 2R s ; I s T vec (I s ) = 158 where for the identity matrix I s = i 1 i 2 i s with i k being the column vectors of I s for k = 1; 2; :::; s; vec (I s ) = 2 6 6 6 6 6 6 6 6 6 6 4 i 1 i 2 . . . i s 3 7 7 7 7 7 7 7 7 7 7 5 : Proof. I s T vec (I s ) = 2 6 6 6 6 6 6 6 6 6 6 4 T T . . . T 3 7 7 7 7 7 7 7 7 7 7 5 2 6 6 6 6 6 6 6 6 6 6 4 i 1 i 2 . . . i s 3 7 7 7 7 7 7 7 7 7 7 5 = 2 6 6 6 6 6 6 6 6 6 6 4 T i 1 T i 2 . . . T i s 3 7 7 7 7 7 7 7 7 7 7 5 = 2 6 6 6 6 6 6 6 6 6 6 4 1 2 . . . s 3 7 7 7 7 7 7 7 7 7 7 5 = Property B.2.7. For any variable 2R n , @ @ = vec (I n ) (B.5) @ T @ = I n (B.6) Using the product rule Equation B.3, we have @v @p = I s T @ _ p @p + @ T @p _ p = @ T @p _ p = @ T @p v = I s ~ I @ _ p @p + @ ~ I @p _ p I s _ T @p @p @ _ T @p p = I s _ T @p @p = I s _ T vec (I s ) (B.7) 159 @v T @p = I s _ p T @ @p + @ _ p T @p = I s _ p T @ @p = I s v T T @ @p = I s _ p T @ ~ I T @p + @ _ p T @p ~ I T I s p T @ _ T @p @p T @p _ = _ (B.8) @v @ _ p = I s T @ _ p @ _ p + @ T @ _ p _ p = I s T @ _ p @ _ p = I s T vec (I s ) (B.9) @v T @ _ p = I s _ p T @ @ _ p + @ _ p T @ _ p = (B.10) The time derivative of the mass matrix can be written as _ M = s X i=1 @M @p i _ p i = _ p T I s @M @p = _ M T n because the mass matrix is symmetric, positive denite = @M @p T _ p T I s T = @M T @p T ( _ p I s ) = @M @p T ( _ p I s ) 160 Appendix C Modeling of Closed-Loop Systems in a Two Dimensional Space In this section of the appendix, we show a special case of control of a spacecraft in a planar reference frame using pulse-modulated thrusters. For a spacecraft moving in a planar reference frame, the spacecraft's rigid body motion has three degrees of freedom: two translational and one rotational and the exible structures are only susceptible to lateral vibration in they-direction. Additionally, we are able to generate the equations of motion in the inertial reference frame with the generalized coordinates using Lagrangian mechanics, because there are no singularities or constraints associated with the parame- ters. For control implementation, there are two thrusters with a xed gimbal angle for applying force and torque and a reaction wheel for small torque requirements, as can be seen in the diagram below. 161 Figure C.1: Translation and Rotation of Spacecraft in a Planar Reference Frame 162 The mode shapes for the exible structures used in deriving the equations of motion for the spacecraft in three dimensions are also used in the derivation of the equations of motion in two dimensions. C.1 Derivation of Equations of Motion using Lagrangian Mechanics To derive the equations of motion, we use the assumed modes methodology of discretizing the kinetic and potential energy of the exible structure, then use Lagrangian mechanics to determine the equations of motions. The kinetic energy of the rigid central body is T b = 1 2 m b _ x 2 + 1 2 m b _ y 2 + 1 2 J b _ 2 = _ p T 2 6 6 6 6 6 6 6 6 6 6 4 m b m b J b 0 3 7 7 7 7 7 7 7 7 7 7 5 _ p = _ p T M b _ p (C.1) where m b is the mass of the central body; J b is moment of inertia; the position of the central body's center of gravity (cg) relative to the inertial reference frame ~ r = 2 6 6 4 x y 3 7 7 5 ; (C.2) 163 is the angular displacement of the central body about its center of gravity (cg); _ p = 2 6 6 6 6 6 6 4 _ r _ _ 3 7 7 7 7 7 7 5 ; _ = 2 6 6 6 6 6 6 4 _ 1 . . . _ N 3 7 7 7 7 7 7 5 (C.3) and n forn = 1; 2;:::;N is the vector of modal amplitudes for then-th exible structure. The gravitational potential energy of the central body is V b =m b gy: (C.4) Figure C.2: Beam De ection Determining the kinetic energy of an appendage is a bit more involved. First we determine the position and velocity vector of dm n , thek-th element of then-the exible structure, then integrate the discretized square norm of the velocity along the exible 164 structure to determine the kinetic energy. The variables in Figures C.1 and C.2 are dened as follows. The position of dm n relative to the point of attachment between the exible structure and the central body ~ r (k) n = 2 6 6 4 x (k) n y (k) n 3 7 7 5 (C.5) Let R b be the radius of the central body and n be the constant angular position of the n-th exible structure relative to the central body's x-axis, then the position of dm n relative to the central body coordinate frame ~ r (k) n =b = 2 6 6 4 cos n sin n sin n cos n 3 7 7 5 ~ r k (C.6) where ~ r k = 2 6 6 4 R b +x (k) n y (k) n 3 7 7 5 (C.7) and R b = ~ R n for all n. The position of dm n relative to the inertial coordinate frame ~ r (k) n =I =~ r + 2 6 6 4 cos sin sin cos 3 7 7 5 2 6 6 4 cos n sin n sin n cos n 3 7 7 5 ~ r k =~ r + 2 6 6 4 cos n sin n sin n cos n 3 7 7 5 ~ r k =~ r + U n I ~ r k (C.8) 165 where n = + n and U n I is the transformation matrix from the coordinate frame xed to the n-th exible structure to the inertial coordinate frame. The position vector is dierentiated to obtain the velocity vector in the inertial reference frame of the n-th exible structure: ~ v (k) n =I = _ ~ r + _ U n I ~ r k + U n I _ ~ r k (C.9) = I 2 U~ r k U n I 2 6 6 6 6 6 6 4 _ ~ r _ _ ~ r k 3 7 7 7 7 7 7 5 (C.10) where _ U n I = _ 2 6 6 4 sin n cos n cos n sin n 3 7 7 5 = _ U: (C.11) 166 Assuming that the exible structure is only susceptible to vibration in the y-direction, the velocity reduces to _ ~ v (k) n =I = I 2 U~ r k (U n I ) 2 2 6 6 6 6 6 6 4 _ ~ r _ _ y (k) n 3 7 7 7 7 7 7 5 = I 2 U~ r k (U n I ) 2 T n 2 6 6 6 6 6 6 4 _ ~ r _ _ n 3 7 7 7 7 7 7 5 = I 2 U~ r k 0 (U n I ) 2 T n 0 2 6 6 6 6 6 6 6 6 6 6 6 6 6 6 4 _ ~ r _ . . . _ n . . . 3 7 7 7 7 7 7 7 7 7 7 7 7 7 7 5 (C.12) = N n _ p (C.13) where (U n I ) 2 is the second column of the matrix U n I ; _ y (k) n = T n _ n ; (C.14) n is the vector of mode shapes; and _ n is the vector of velocities of the modal amplitudes. The norm of the appendage's velocity is ~ v (k) n =I 2 = _ p T N n T N n _ p (C.15) 167 where N n T N n = 2 6 6 6 6 6 6 6 6 6 6 6 6 6 6 4 I 2 U~ r k 0 (U n I ) 2 T n 0 ~ r T k U T ~ r T k ~ r k . . . ~ r T k (U n I ) 2 T n . . . 0 0 n (U n I ) T 2 n (U n I ) T 2 U~ r k . . . n T n . . . 0 0 3 7 7 7 7 7 7 7 7 7 7 7 7 7 7 5 : (C.16) The kinetic energy of the nth appendage T n = 1 2 Z an ~ v (k) n =I 2 dm n = 1 2 _ p T Z an N n T N n dm n _ p = 1 2 _ p T M n _ p (C.17) The potential energy of the nth appendage V n = V gn +V sn = m n g y + R n + 1 2 L n sin n + 1 2 Z Ln 0 EI d 2 y n H dx 2 n H ! 2 dx n H = m n gy +m n g R n + 1 2 L n sin n + 1 2 ~ T n K n ~ n (C.18) where V gn is the gravitational potential energy; V sn is the elastic potential energy of the nth appendage; L n and m n is the length and mass of the n-th exible structure, respectively; the stiness matrix for the n-th exible structure K n =E Z Ln 0 I(x) n (x) T n (x)dx; (C.19) 168 E is the Young's modulus; and I(x) is the area moment of inertia. The total kinetic energy of the system T =T b + N X n=1 T n = _ p T M b _ p + N X n=1 _ p T M n _ p = _ p T M _ p (C.20) The total potential energy of the system V =V b + N X n=1 V gn +V sn =m b gy + N X n=1 m n g y + sin n R b + 1 2 L n + 1 2 T n K n n : (C.21) The Lagrangian of a system L =TV (C.22) is substituted into the Euler-Lagrange equation d dt @L @ _ p @L @p =Q (C.23) to generate the equations of motion. The Euler-Lagrange equation can be rewritten in terms of the kinetic and potential energy as d dt @T @ _ p @T @p + @V @p =Q (C.24) 169 We calculate the derivatives of the kinetic energy using the Kronecker product following the methodology presented in Control of Robot Manipulators [16]. Consider d dt @T @ _ p = M p + _ M _ p (C.25) and @T @p = 1 2 I s _ p T @M @p _ p: (C.26) where the subscript s is the length of p. The Coriolis/centripetal vector is denoted C (p; _ p) = _ M _ p 1 2 I s _ p T @M @p _ p (C.27) = _ M 1 2 U _ p (C.28) where U = I s _ p T @M @p : (C.29) The derivative of the system's mass matrix can be written as _ M = _ p 1 I s _ p s I s 2 6 6 6 6 6 6 6 4 @M @p 1 . . . @M @p s 3 7 7 7 7 7 7 7 5 = _ p T I s @M @p (C.30) = @M @p 1 @M @p s 2 6 6 6 6 6 6 4 _ p 1 I s . . . _ p s I s 3 7 7 7 7 7 7 5 = @M @p T ( _ p I s ): (C.31) 170 Also, note that (I s _ p) _ p = ( _ p I s ) _ p = 2 6 6 6 6 6 6 6 6 6 6 4 _ p 1 _ p _ p 2 _ p . . . _ p s _ p 3 7 7 7 7 7 7 7 7 7 7 5 (C.32) and for the Kronecker product (A B) T = A T B T : (C.33) Using Equations C.31-C.33, we can show _ M _ p = U T _ p: (C.34) Let's consider U T _ p = @M @p T (I s _ p) _ p = @M @p T ( _ p I s ) _ p = _ M _ p: (C.35) 171 This allows us to right the Coriolis/centripetal vector C (p; _ p) = _ M _ p 1 2 U _ p (C.36) = 1 2 _ M _ p + 1 2 U T _ p 1 2 U _ p (C.37) = C (p; _ p) _ p (C.38) where C (p; _ p) = 1 2 h _ M + U T U i (C.39) is the matrix of Coriolis/centripetal terms. Furthermore, the matrix N = 1 2 _ M C (p; _ p) = 1 2 U T U (C.40) is skew-symmetric. The partial derivative of the potential energy gives @V @p = Kp +G(p) (C.41) where the stiness matrix for the system K = 2 6 6 6 6 6 6 6 6 6 6 4 0 33 K 1 . . . K N 3 7 7 7 7 7 7 7 7 7 7 5 (C.42) 172 and the gravity vector G(p) = 2 6 6 6 6 6 6 6 6 6 6 6 4 0 mg g N P n=1 m n R b + 1 2 L n cos n ~ 0 3 7 7 7 7 7 7 7 7 7 7 7 5 : (C.43) The damping matrix D =T T (2 )T (C.44) where is the damping ratio; and T is the matrix of eigenvectors and is the diagonal matrix of eigenvalues of the generalized eigenvalue problem MT = KT: (C.45) The generalized equations of motion with respect to the inertial coordinate frame can be written as M (p) p + C (p; _ p) _ p + D _ p + Kp +G (p) = ~ Bu I a (C.46) where M (p) is the ss nonlinear positive denite mass matrix; C (p; _ r) is an ss matrix containing the Coriolis and coupling terms; D is the ss damping matrix determined by calculating the natural frequencies after linearizing Equation (C.46) using modal analysis; 173 K is thess positive semidenite stiness matrix determined by the mode shapes; ~ B = 2 6 6 4 I 3 7 7 5 is the input matrix; I is a 3 3 identity matrix; is the input matrix for the actuators on the exible structures; u I a = f I T a a ^ u f T is the actual input to the system in the inertial reference frame wheref I a is the force applied inertial reference frame; a is the applied torque; and~ u fa is the input applied to compensate for any vibration The matrix C is derived such that the following property holds. Property C.1.1. The matrix N (p; _ p) = 1 2 _ M (p) C (p; _ p) (C.47) is skew-symmetric where _ M (p) is the derivative of the mass matrix, and therefore has the property z T N (p; _ p)z = 0 (C.48) for any vector z2R n . C.2 Implementation of Commanded Control For the generalized form the equations of motion for the system, we quantify the un- certainty due to thrust implementation, which is a result of the inability to apply the commanded input continuously without delay. The commanded input is used to deter- mine which thrusters are ignited and the magnitude of the thrust from each thruster. 174 The commanded input also determines the force applied by actuators on the exible structures. The commanded input, which is calculated in the inertial coordinate frame, u I c = 2 6 6 6 6 6 6 4 f I c c u fc 3 7 7 7 7 7 7 5 (C.49) where f I c is the vector of commanded forces in the x and y-direction of the inertial coordinate frame; c is the commanded torque in the coordinate frame xed to the central body; and u fc is the commanded input for damping vibration in the exible structures. The commanded input is used to determineu I a , the input applied to the system. The rst step in determining the input applied to the spacecraft is to transform the commanded forcef I c into the central-body xed reference frame. The transformation matrix from the inertial reference frame to the central-body xed reference frame is U B I = 2 6 6 4 cos sin sin cos 3 7 7 5 : (C.50) 175 Note that the transformation matrix is a function , which is determined from sensor output thereby introducing uncertainty. The transformation matrix determined using the sensor measurements is denoted ^ U B I = 2 6 6 4 cos ^ sin ^ sin ^ cos ^ 3 7 7 5 = 2 6 6 4 cos ( +) sin ( +) sin ( +) cos ( +) 3 7 7 5 = 8 > > < > > : 2 6 6 4 1 0 0 1 3 7 7 5 + 2 6 6 4 0 0 3 7 7 5 9 > > = > > ; 2 6 6 4 cos sin sin cos 3 7 7 5 =: fI +g U B I (C.51) where the real angular displacement = ^ +; ^ is the angular displacement based on the the sensor output; is the sensor uncertainty. Note that terms denoted by ^ [ ] are dependent on sensor output. With respect the central body, the commanded force is ^ f B c =fI +g U B I f I c = U B I f I c + U B I f I c =f B c + U B I f I c (C.52) 176 The force ^ f B c to determine the thrust magnitude from the left and right thrusters. See Figure C.1 for thruster arrangement. The magnitude of the thrust required of each thruster is 2 6 6 6 6 6 6 4 ^ f Lc ^ f Rc rw 3 7 7 7 7 7 7 5 = 2 6 6 6 6 6 6 4 sin sin 0 cos cos 0 r l sin ( ) r r sin ( ) 1 3 7 7 7 7 7 7 5 2 6 6 4 ^ f B c c 3 7 7 5 =: U T B 2 6 6 4 ^ f B c c 3 7 7 5 (C.53) The thrusters can not be implemented exactly as commanded because of thrust imple- mentation error. Therefore the applied thrust is f B ta = 2 6 6 6 6 6 6 4 f La f Ra a 3 7 7 7 7 7 7 5 = 2 6 6 6 6 6 6 4 ^ f Lc ^ f Rc rw 3 7 7 7 7 7 7 5 + 2 6 6 6 6 6 6 4 f Lc f Rc rw 3 7 7 7 7 7 7 5 (C.54) where f t = f Lc f Rc rw is bounded uncertainty due to the modulation of the thrusters. Also, any other bounded uncertainties that are known to exist when applying thrust can be accounted for in this term. While the thrusters are applied with respect to spacecraft, for the purposes of analysis and simulation, the applied thrust must be transformed to the inertial reference frame which gives us f B a =f I c +U I B U B I f I c +U I B U B T f T t =:f I a +f a (C.55) 177 The arrangement of the thrusters also applies a torque which can be used for large angle maneuvers. The torque applied by the thrusters and reaction wheel ta =f La sin ( ) +f Ra sin ( ) + rw (C.56) The implementation of vibration suppression is similar to that of the exible structure in a three dimensional space as presented in Chapter 3. 178 Appendix D Trajectory Tracking in a Two Dimensional Space In the following stability proof, the control for trajectory tracking was designed for a model in a two-dimensional space. In this section we dene the controller and prove boundedness of the tracking error. Proposition D.0.1. Assume that we have the same number of actuators as the number of modes we want to control. If the commanded input u c =u 0 +u f +u = 2 6 6 6 6 6 6 4 f c c f c 3 7 7 7 7 7 7 5 (D.1) where f c represents the forces in the x and y direction; c is the torque applied to the body; f c vector of control inputs to the modal amplitudes. The nominal control term u 0 = M (p d ) + C (p d ; _ p d ) _ +G (p d ) (D.2) 179 where p d , _ p d , and p d denote the position, velocity, and acceleration of the commanded trajectory, respectively; _ = _ r d e; = r d _ e: (D.3) with position error e =pp d and velocity error _ e = _ p _ r d . The feedback control term u f =K d =K d _ r _ =K d ( _ e + e) ; (D.4) that is, = _ r _ = _ e e; (D.5) and it follows that _ = r = e _ e: (D.6) The error compensation term u =k n (;) (D.7) where k >kuk (D.8) for the bounded uncertaintykuk entering the system through the input, and n (;!) = 8 > > < > > : 0 if k!k< ! k!k if k!k (D.9) Additionally, the function n (;!) has the following properties 180 1. k!kn (;!) =kn (;!)k! 2. kn (;!)k 1 k!k Then the solutions to the system M (p) p + C (p; _ p) _ p + D _ p + Kp +G (p) = ^ B (u c +u) (D.10) are uniformly bounded. Proof. Consider the following system, M (r) r +C (r; _ r) +D _ r +Kr +G (r) = ^ B (u cm +u t ): (D.11) whereu cm =u 0m +u fm +u m is commanded input determined using sensor measurements and assumptions in modeling of the system and u t represents the thrust implementa- tion errors. First, we substitute the nominal control term into the system, but due to it's dependence on states measured by the sensors, we need to account for the sensor uncertainty. We denote terms that are dependent on sensor measurements as e m = e +e (D.12) _ e m = _ e + _ e (D.13) then _ m = _ r d e m = _ r d e e = _ e m = r d _ e m = r d _ e _ e = _ e: (D.14) 181 The measured nominal control term is u 0m = M (p d ) m + C (p d ; _ p d ) _ m +G (p d ) = M (p d ) + _ e + C (p d ; _ p d ) _ + e +G (p d ) = M (p d ) + C (p d ; _ p d ) _ +G (p d ) + M (p d ) _ e + C (p d ; _ p d ) e = [M (r) +M] + [C (r; _ r) +C] _ + [G (r) +G] +M (r d ) _ e +C (r d ; _ r d ) e = M (r) + C (r; _ r) _ +G (r) +u 0 (D.15) Substituting Equation D.15 into the equations of motion, and then using Equation D.5 and D.6 gives M (p) _ =C (p; _ p) D _ e Ke +u fm +u m +u t +u 0 : (D.16) The sensor uncertainty must also be accounted for when implementing the feedback control term. Let m = +, where = _ e + e, then the measured feedback term u fm =K d m =K d K d =K d +u f (D.17) Substituting Equation D.17 into the system gives M (r) _ = C (p; _ p) D _ e Ke K d +u m +u f +u 0 +u t = C (p; _ p) D _ e Ke K d +u m +u (D.18) 182 Let x T = T e T , and dene V (x) = 1 2 T M (p) + 1 2 e T (K + D)e (D.19) = 1 2 T e T 2 6 6 4 M (p) K + D 3 7 7 5 2 6 6 4 e 3 7 7 5 (D.20) = 1 2 x T Mx (D.21) to be a candidate for the Lyapunov function, where K = 2 6 6 4 0 K 3 7 7 5 (D.22) for K positive denite. Similarly, D is also positive denite. Note that the matrix M (p) is positive denite, thus V (x) 0. Furthermore 1 2 min (M)kxk 2 V (x) 1 2 max (M)kxk 2 (D.23) 183 where min and max denote the minimum and maximum eigenvalues of a matrix, re- spectively. Dierentiating (D.19) and then substituting (D.18) , gives _ V (x) = T M (r) _ + 1 2 T _ M (r) +e T (K + D) _ e = T [C (r; _ r) D _ e Ke K d +u m +u] + 1 2 T _ M (r) +e T (K + D) _ e = T K d ( _ e + e) T (D _ e + Ke) +e T (K + D) _ e + T (u m +u) = T K d _ e T D _ ee T Ke + T (u m +u) T K d e T Ke + T (u m +u) x T Kx + T (u a +u) min (K)kxk 2 + T (u a +u) (D.24) For the nominal system, i.e. the system with out uncertainty, M (p) _ =C (p; _ p) D _ e Ke K d ; (D.25) the derivative of the Lyapunov function Equation (D.19) is _ V = T K d _ e T D _ ee T Ke; (D.26) for all e; _ e2R n and _ V = 0 when , e, and _ e are all zero. Therefore, the nominal closed- loop system is globally asymptotically stable. The error compensation term u is derived 184 using Lyapunov Redesign. First we need to determine the upper bound on the uncertainty term E, T (u +u m ) (D.27) Again, due to sensor uncertainty the error compensation term can not be applied exactly; that is, the error compensation term is implemented as u m =k n (; m ) = 8 > > < > > : 0 ifk m k< m k m k ifk m k (D.28) and recall that m = +. Let = max . Recall that m = +, where is determined by the sensor measurement uncertainty. Let = maxkk be the maximum absolute uncertainty in due sensor measurement inaccuracies. If the thrusters are turned o, that is u m = 0, then E = T u = ( m +) T u k m kkuk +kkkuk + ,E off (D.29) 185 If the thrusters are on, then E, T (u +u m ) = ( m ) T (u n ;m ) = T m T u kn ;m k m k m k = T m T u T m T kn ;m k m k m k = T m T u kn ;m k k m k 2 k m k + kn ;m k k m k T m k m kkuk +kkkuk kn ;m kk m k + kn ;m k k m k kkk m k =k m kkuk +kkkuk + kn ;m kkk kn ;m kk m k k m k + + kn ;m k kn ;m kk m k = k m k (1kn ;m k) + (1 +kn ;m k) k m k kn ;m k + 2 = + 2 ,E on : (D.30) We see that the disturbances entering the system are greatest when the thrusters are on. Denote the upper bound on the disturbances as E =E on . If _ V (x)< min (K)kxk 2 +E < 0 (D.31) then the system is stable when kxk> E min (K) 1 2 =: (D.32) 186 If we dene the ultimate bound as r = max (M) min (M) ; (D.33) then x eventually converges to a ball B r with radius r for all t>t 0 ; therefore and e are also bounded. Now recall that = _ e + e; (D.34) this can be rewritten as _ e = e: (D.35) This equation represents a linear time-invariant feedback system with plantG(s) = 1 and feedback H(s) = 1 s . We see that both the plant and feedback are strictly positive real function, thus the output of the system _ e is bounded when is bounded. Similarly, equation (D.34) can be rewritten as e = 1 ( _ e) (D.36) and this too is a linear time-invariant feedback system with plant G(s) = 1 and feed- backH(s) =s and both plant and feedback are strictly positive real functions. Thus, for bounded input , the output e is bounded. Due to limits on the magnitude of force that can be applied by the thruster and actuators, the solutions are not globally uniformly bounded. Boundedness can only be 187 guaranteed when the commanded inputu c to the thrusters and actuators is less than the maximum force U max that can be applied. 188

Abstract (if available)

Abstract A methodology for the design of robust control systems for flexible spacecraft undergoing a large angle maneuver using pulse-modulated thrusters is proposed. The need for developing such a methodology stems from the need to reduce the weight of the modern spacecraft and to satisfy strict performance specifications resulting in vibrational modes within the control system bandwidth. Moreover, the discontinuous operation of pulse-modulation may excite large flexible motion that can lead to deterioration of performance, limit cycles, and even instability. Consequently, the current practice of designing the control system using a rigid body model and evaluating performance by simulation using a flexible model is inadequate for the modern spacecraft. In the proposed methodology structural flexibility is included in the control design model thus resulting in control systems with robust stability and improved performance. ❧ This dissertation has two main parts: formulation of flexible spacecraft dynamics performing large angle motion and control system design. The formulation of equations of motion is performed using quasi-coordinates that allow for parameterization of the attitude in terms of Euler-Rodrigues parameters that avoid singularities. Moreover, the formulation is in terms of variables measurable by the control system, such as angular velocities and is in a form that is convenient for establishing closed-loop robust stability. The control law structure proposed here results in guaranteed closed-loop performance in the presence of parameter uncertainties, measurement noise, and command implementation errors. A theorem that guarantees the bounds of trajectory following errors of a spacecraft performing a general translational and rotational motion in three-dimensional space is formulated and proven using Lyapunov stability theory. Finally, a model of a flexible spacecraft is used to demonstrate the performance of a control system developed using the proposed methodology are presented. A number of different maneuvers are studied by simulation showing satisfactory performance in the presence of uncertainties and implementation errors.

Linked assets

University of Southern California Dissertations and Theses

Conceptually similar

PDF

Control of two-wheel mobile platform with application to power wheelchairs

PDF

Active delay output feedback control for high performance flexible servo systems

PDF

Dynamic analysis and control of one-dimensional distributed parameter systems

PDF

Robust control of periodically time-varying systems

PDF

Nonlinear control of flexible rotating system with varying velocity

PDF

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

PDF

Design, modeling and analysis of piezoelectric forceps actuator

PDF

An approach to dynamic modeling of percussive mechanisms

PDF

Mathematical model and numerical analysis of a shape memory polymer composite cantilever beam for space applications

PDF

Modeling and vibration analysis of wheelchair-users

PDF

Modeling, analysis and experimental validation of flexible rotor systems with water-lubricated rubber bearings

PDF

Building cellular self-organizing system (CSO): a behavior regulation based approach

PDF

Flexible formation configuration for terrain following flight: formation keeping constraints

PDF

Controlled and uncontrolled motion in the circular, restricted three-body problem: dynamically-natural spacecraft formations

PDF

Periodic solutions of flexible systems under discontinuous control

PDF

Incorporation of mission scenarios in deep space spacecraft design trades

PDF

Large eddy simulations of laminar separation bubble flows

PDF

On the synthesis of dynamics and control for complex multibody systems

PDF

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

PDF

Using nonlinear feedback control to model human landing mechanics

Asset Metadata

Creator Reddy, Shalini A. (author)

Core Title Control of spacecraft with flexible structures using pulse-modulated thrusters

School Viterbi School of Engineering

Degree Doctor of Philosophy

Degree Program Mechanical Engineering

Publication Date 08/02/2012

Defense Date 06/19/2012

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

Tag Control,dynamics,OAI-PMH Harvest,spacecraft

Language English

Contributor Electronically uploaded by the author (provenance)

Advisor Flashner, Henryk (committee chair), Kukavica, Igor (committee member), Redekopp, Larry G. (committee member), Shiflett, Geoffrey R. (committee member), Yang, Bingen (Ben) (committee member)

Creator Email sreddy6@gmail.com

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

Unique identifier UC11290190

Identifier usctheses-c3-84663 (legacy record id)

Legacy Identifier etd-ReddyShali-1110.pdf

Dmrecord 84663

Document Type Dissertation

Rights Reddy, Shalini A.

Type texts

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

Access Conditions The author retains rights to his/her dissertation, thesis or other graduate work according to U.S. copyright law. Electronic access is being provided by the USC Libraries in agreement with the a...

Repository Name University of Southern California Digital Library

Repository Location USC Digital Library, University of Southern California, University Park Campus MC 2810, 3434 South Grand Avenue, 2nd Floor, Los Angeles, California 90089-2810, USA