<%BANNER%>

Robust Nonlinear Attitude Control With Disturbance Compensation


PAGE 1

ROBUST NONLINEAR ATTITUDE CONTROL WITH DISTURBANCE COMPENSATION By KEVIN J. WALCHKO A DISSERTATION PRESENTED TO THE GRADUATE SCHOOL OF THE UNIVERSITY OF FLORIDA IN PARTIAL FULFILLMENT OF THE REQUIREMENTS FOR THE DEGREE OF DOCTOR OF PHILOSOPHY UNIVERSITY OF FLORIDA 2003

PAGE 2

Copyright 2003 By Kevin J Walchko

PAGE 3

iii ACKNOWLEDGEMENTS I would like to thank Dr. Paul Mason for his guidance of my education and research. He has always made every effort to aid me in my work and life. Paul has held many roles in my life such as mentor, academic advisor, friend, and best man at my wedding. He has become more than an advisor over the years, and will forever be part of my family. I will always be in his debt. I would like to thank Dr. John Schueller who acted as my co-chair during this dissertation. He provided much advice and guidance in the bureaucracy of University of Florida. Schueller was always willing to help me with various paperwork and other such matters pertaining to fellowships, pay checks, and graduation. He is a good friend and I look forward to working with him in the future. I would also like to thank my family who sacrificed much at times so this could be possible. My wife Nina took great care of me and our new son Michael. She made sure I had plenty of time to do my work and never complained about anything. She is the source from which I draw all of my strength, ambition, and dedication. Finally, I would like to thank the horde, our pack of cats. They provided much amusement and stress relief through their wacky antics. I treasure each and every one of them, and will never forget those who have passed away.

PAGE 4

iv T ABLE OF CONTENTS page ACKNOWLEDGEMENTS ..............................................................................................iii LIST OF FIGURES ........................................................................................................viii LIST OF TABLES ............................................................................................................xi LIST OF SYMBOLS ......................................................................................................xiii ABSTRACT ...................................................................................................................... xv 1 INTRODUCTION ..........................................................................................................1 Formation Flying ..........................................................................................................1 General Problems Controlling Single Satellites ...........................................................2 External Disturbances .............................................................................................2 Parametric Uncertainty ...........................................................................................3 Literature Review ..........................................................................................................4 Fuel Slosh ...............................................................................................................4 Modelling fuel slosh dynamics..........................................................................4 Controlling fuel slosh.........................................................................................6 Do we need to worry about fuel slosh?..............................................................9 Missions involving fuel slosh disturbances.....................................................12 Thermally Induced Vibrations: Solar Snap ...........................................................12 Controlling thermal vibrations.........................................................................13 Missions which involved solar snap................................................................15 Parametric Uncertainty .........................................................................................16 Summary and Overview of the Thesis ........................................................................17 2 CONTROLLING THE UNKNOWN: FUEL SLOSH .................................................19 Introduction .................................................................................................................1 9 Modelling Fuel Slosh Dynamics .................................................................................20 Fuel Slosh Models ......................................................................................................23 Simulation Model ..................................................................................................23 Controller Model ...................................................................................................24

PAGE 5

v 3 DYNAMICS OF SOLAR ARRAYS ............................................................................26 Introduction .................................................................................................................2 6 Thermally Induced Vibrations ....................................................................................26 Heat Flux ...............................................................................................................27 Summing Forces ...................................................................................................28 Nonhomogeneous Boundary Conditions ..............................................................33 Separation of Variables .........................................................................................34 Solving the position function F(x)...................................................................35 Solving the time function G(t).........................................................................38 Summary of results..........................................................................................38 Thermally Induced Dynamics .....................................................................................39 Modal Equations ...................................................................................................39 Disturbance Torques .............................................................................................40 Thickness Temperature Gradient ..........................................................................42 Conclusion ..................................................................................................................44 4 SATELLITE ATTITUDE DYNAMICS ......................................................................46 Reference Frames ........................................................................................................46 Earth Centered Inertial ..........................................................................................46 Orbital Frame ........................................................................................................47 Body Fixed Frame .................................................................................................47 Fun With Vectors and Rotating Coordinate Systems ..................................................47 Spacecraft Attitude ......................................................................................................49 Quaternions ...........................................................................................................50 Rotations of Rigid Bodies in Space. .....................................................................51 Spacecraft Equations of Motion ..................................................................................51 Attitude Kinematics ..............................................................................................52 Spacecraft Dynamics ............................................................................................52 Simple Spacecraft .................................................................................................52 Angular Velocity ...................................................................................................56 Internal Disturbances ..................................................................................................56 Spacecraft with Reaction or Momentum Wheels .................................................56 Fuel Slosh .............................................................................................................60 Environmental Disturbances in Space ........................................................................61 Gravity Gradient ...................................................................................................62 Solar Snap .............................................................................................................62 Conclusion ..................................................................................................................62 5 ATTITUDE CONTROL OF SPACECRAFT ...............................................................64 Satellite Control Hardware .........................................................................................64 Control Computational Hardware .........................................................................64 Control Actuators ..................................................................................................65 Reaction and momentum wheels.....................................................................65

PAGE 6

vi Magnetic torque rods.......................................................................................66 Gas jets.............................................................................................................67 Proportional Derivative Control .................................................................................67 Sliding Mode Control Theory and Background ..........................................................69 What is Sliding Mode? ..........................................................................................69 Theory ...................................................................................................................69 Simple Example ....................................................................................................72 Other Sliding Mode Designs for Spacecraft ...............................................................73 Satellite Controller Design For This Work .................................................................75 Controller Derivation ............................................................................................75 Fuel Slosh Disturbance .........................................................................................77 Stability of Controller ...........................................................................................78 Saturation ..............................................................................................................80 Comparing The Two Types of Controller ...................................................................81 Conclusion ..................................................................................................................81 6 SIMULATION AND RESULTS ..................................................................................83 Simulation Dynamics and Parameters ........................................................................83 Satellite .................................................................................................................83 Satellite Dynamics ................................................................................................85 Fuel Tank ...............................................................................................................86 Reaction Wheels ...................................................................................................86 Solar Panels ...........................................................................................................87 Control System .....................................................................................................87 Finding an Optimal PD Gain ......................................................................................88 Fuel Slosh ...................................................................................................................9 3 IAE ........................................................................................................................94 ITAE ......................................................................................................................95 Settling Time .........................................................................................................95 CM ........................................................................................................................96 T orque and Momentum Limits .............................................................................96 Qualitative Analysis of Results .............................................................................98 Solar Snap ...................................................................................................................9 8 Performance ..........................................................................................................98 Qualitative Analysis of Results .............................................................................99 7 CONCLUSION ...........................................................................................................101 A ATTITUDE REPRESENTATIONS AND ROTATION MATRICES ......................104 Fixed Angle Rotations ...............................................................................................104 Euler Angles ..............................................................................................................105 Quaternions ...............................................................................................................107 Quaternion Algebra .............................................................................................108 Rotations of Rigid Bodies in Space. ...................................................................110

PAGE 7

vii Attitude Errors in Quaternions ..................................................................................111 Summary of Quaternions ..........................................................................................112 B MODAL EQUATIONS .............................................................................................113 Motion of Multi-Degree of Freedom System ...........................................................113 Dampening ................................................................................................................115 Simple Example ........................................................................................................116 C SIMULATION CODE ...............................................................................................117 Simulation Code ........................................................................................................117 Matlab .................................................................................................................117 C Programming Language ..................................................................................117 C++ Programming Language ..............................................................................118 cMathlib Header File ................................................................................................119 cRK4 Header File .....................................................................................................126 REFERENCES ...............................................................................................................127 BIOGRAPHICAL SKETCH ..........................................................................................134

PAGE 8

viii LIST OF FIGURES Figure pa ge 1-1 Finite element analysis of a fuel tank inside of a satellite. ......................................6 1-2 PPT for EO-1. ........................................................................................................ 10 1-3 Diagram of a PPT. .................................................................................................10 1-4 NEAR spacecraft orbiting the asteroid Eros. .........................................................11 1-5 Hubble Space Telescope diagram. .........................................................................13 2-1 Traditional fuel slosh pendulum model where M is the total fuel mass center, mp is the pendulum mass, ms is the stationary mass fraction, L is the pendulum arm length, and g is the local acceleration due to some external force. .......................22 2-2 Propagation of fuel slosh in an elliptical tank over time. These images demonstrate the difficulty in modelling the fluid movement with simple representations. ......22 2-3 Diagram of a single sloshing mode, with its mass attached to the fuel tank walls by springs and dampers aligned along the x, y, and z axes. ......................................23 2-4 Plot of the effect of multiple modes in the velocity of a simply supported beam. This is representative of the fuel slosh disturbance in flight data. ................................24 3-1 The Soviet satellite Sputnik which orbited the Earth in October 1957. ................27 3-2 Diagram of the penumbra and umbra regions which block a satellite from sunlight. ............................................................................................................................... .28 3-3 Another view showing the path of a satellite and how its orbit does effect the heat flux it is exposed to. ..............................................................................................28 3-4 Cross section of solar panel with coordinate systems superimposed. Note that the S’s in the figures represent ‘s and M’s are moments. ...........................................30 3-5 Free body diagram of a section of the beam of size dx where M(x,t) is a moment, V(x,t) is a shear force, f(x) is a force per unit length, and u(x,t) is the height of the mass element at position x at time t. .....................................................................32

PAGE 9

ix 3-6 Plot of [3-49]. Intersections of the two lines are solutions of the characteristic equation. ................................................................................................................36 3-7 A cross section of the beam in which a heat flux is applied to one side. Notice the beam starts off at a uniform temperature. .............................................................42 3-8 Plot of the temperature change of a beam. Notice that the temperature difference between the two sides remains constant eventually. .............................................42 3-9 Hubble Space Telescope which shows the solar panels bent due to a constant temperature gradient. ............................................................................................43 3-10 Another picture of the Hubble Space Telescope showing deflection of the solar panels during maintenance by the space shuttle. ..................................................44 4-1 The three frames of reference commonly used in spacecraft dynamics. ...............47 4-2 A reference frame in motion relative to an inertially fixed frame of reference. ....48 4-3 A rigid body satellite. ............................................................................................53 4-4 A rigid body spacecraft with a single reaction wheel. ...........................................57 4-5 A plot of common environmental disturbances a spacecraft is subjected too. ......60 5-1 Two types of reaction wheels produced by Ithaco. ...............................................66 5-2 Torque rods produced by Ithaco. ...........................................................................67 5-3 Spacecraft engine produced by Boeing. ................................................................68 5-4 Two possible solutions for the double integrator, where (a) is and (b) is These plots are typically referred to as phase portraits. They plot position and velocity in this example, but later (when we are doing controls with full state feed back) the portraits will represent error and error velocity. ...................................................70 5-5 Switching between the two controllers results in driving the system to zero in a stable manor. ...................................................................................................................71 5-6 Graphical representation of how effects the sliding surface’s orientation in the phase plane. The dotted lines represent the sliding surface (s = 0). ................................72 6-1 EO-1 satellite with solar panel deployed. ..............................................................83 6-2 Diagram of satellite’s internal components. ..........................................................84

PAGE 10

x 6-3 Composition of a typical simple solar panel. ........................................................86 6-4 Cost function surface plot for the optimal gain of the PD controller. ...................89 6-5 ITAE results from the optimal PD search. ............................................................90 6-6 IAE results from the optimal PD search. ...............................................................90 6-7 Settling time results from optimal PD search. .......................................................91 6-8 Control momentum or fuel used searching for the optimal PD controller. ...........92 6-9 A depiction of the step maneuver showing the orientation of the satellite at start and end. The angular distance and axis of rotation are also shown for better understanding of the terms. ...................................................................................94 6-10 Quaternion error and control effort of the PD controller rotating from -30 to 30. ............................................................................................................................... .96 6-11 Quaternion error and control of the SM controller rotating from -30 to 30. .......97 6-12 Quaternion error and control of the SMCNF controller rotating from -30 to 30. ............................................................................................................................... .97 6-13 Arcsecond error of the PD controller during solar snap. .....................................99 6-14 Arcsecond error of the SM controller during solar snap. ..................................100 6-15 Arcsecond error of the SMCNF controller during solar snap. ..........................100 A-1 Body reference frame attached to a rigid body. The x-axis points out the front of the vehicle. ................................................................................................................105 B-1 Simple multi-degree of freedom system composed of mass and spring elements. ..............................................................................................................................1 13 C-1 UML diagram of Mathlib class which contains C++ objects for vectors, matricies, and quaternion used in this work. .......................................................................119

PAGE 11

xi LIST OF TABLES T able page 1-1 Spacecraft with thermally induced disturbances. ...................................................15 3-1 Solar panel thermal material properties. .................................................................32 3-2 First four eigen modes of the beam ........................................................................35 5-1 Comparison of common active attitude control hardware. .....................................65 5-2 Summary of the switching law. ..............................................................................70 6-1 Technical Specifications for the EO-1 spacecraft ..................................................84 6-2 EO-1 fuel slosh parameters. ....................................................................................86 6-3 EO-1 wheel specifications. .....................................................................................86 6-4 EO-1 solar panel parameters. ..................................................................................87 6-5 Initial and final euler angles where roll, pitch, and yaw are equal. ........................93 6-6 IAE results with values having units of radians. ....................................................94 6-7 ITAE results with values having units of radian-seconds2. ...................................95 6-8 Settling time results with values having units of seconds. .....................................96 6-9 CM results with values having units of Nms. .........................................................96 6-10 Maximum wheel momentum achieved where units are in Nms. .........................98 6-11 Solar snap results of the three controllers. ............................................................99 A-1 Properties of a rotation matrix. ............................................................................105 A-2 Comparison of the two major types of rotations, sequential and first and third axes rotation ................................................................................................................106

PAGE 12

xii A-3 Quaternion Algebra Summary .............................................................................109

PAGE 13

xiii LIST OF SYMBOLS Thermally Induced Vibrations Symbols Absorptivity Coefficient of thermal expansion. Natural frequency of mode n Density of a material. Thermal diffusivity. A Cross sectional area of the solar panel. c Heat capacity. EY oung’s Modulus F(x)Position function from the separation of variables method. F n Position function for mode n G(t)Time function from the separation of variables method. h Thickness of the solar panel. k Thermal conductivity. M T (t)Thermal moment. T otal heat flux on surface. Heat flux from sun. Earth emitted heat flux. Earth reflected heat flux. u(x,t)Displacement of solar panel. V(x,t)Shear force. Satellite Dynamics and Fuel Slosh Symbols Angular velocity of satellite system relative to inertial frame. Angular velocity of satellite system relative to orbital frame. Angular velocity of orbital frame relative to inertial frame. Tn qoqsqeqaBI /BO /OI /

PAGE 14

xiv C Rotation matrix. c First moment of inertia. c fn Damper coefficient for sloshing mode n F External forces and torques. h Angular momentum. J Second moment of inertia. k fn Spring coefficient for sloshing mode n L Length of pendulum. M Mass matrix of satellite system. mpPendulum mass. msStationary mass fraction. p Linear momentum. q Quaternion The imaginary part of the quaternion, 1-3 R c Distance of satellite from center of Earth. T Kinetic energy of satellite system. V Composite linear and angular velocity of satellite system. Controls Symbols W eighting term on the sliding surface. K d Derivative controller gain. K p Proportional controller gain. PDProportional Derivative controller s Sliding surface. SMSliding Mode controller. SMCNFSliding Mode controller with a colored noise filter T d Gain ratio between the proportional and derivative gains. Satellite Dynamics and Fuel Slosh Symbols q ˆ

PAGE 15

xv Abstract of Dissertation Presented to the Graduate School of the University of Florida in Partial Fulfillment of the Requirements for the Degree Doctor of Philosophy ROBUST NONLINEAR ATTITUDE CONTROLWITH DISTURBANCE COMPENSATION By Kevin J Walchko May 2003 Chair: Paul Mason Cochair: John Schueller Major Department: Mechanical and Aerospace Engineering Attitude control of small spacecraft is a particularly important component for many missions in the space program: Hubble Space Telescope for observing the cosmos, GPS satellites for navigation, SeaWiFS for studying phytoplankton concentrations in the ocean, etc. Typically designers use proportional derivative control because it is simple to understand and implement. However this method lacks robustness in the presence of disturbances and uncertainties. Thus to improve the fidelity of this simulation, two disturbances were included, fuel slosh and solar snap. Fuel slosh is the unwanted movement of fuel inside of a fuel tank. The fuel slosh model used for the satellite represents each sloshing mode as a mass-spring-damper. The mass represents the wave of fuel that propagates across the tank, the damper represents the baffling that hinders the movement, and the spring represents the force imparted to the spacecraft when the wave impacts the tank wall. This formulation makes the incorporation

PAGE 16

xvi of multiple modes of interest simple, which is an advance over the typical one sloshing mode, pendulum model. Thermally induced vibrations, or solar snap, occur as a satellite transitions from the day-to-night or night-to-day side of a planet. During this transition, there is a sudden change in the amount of heat flux to the solar panels and vibrations occur. Few authors have looked at the effects of solar snap. The disturbance dynamics were based on the work by Earl Thorten. The simulated effects compared favorably with real flight data taken from satellites that have encountered solar snap. A robust sliding mode controller was developed and compared to a more traditional proportional derivative controller. The controllers were evaluated in the presence of fuel slosh and solar snap. The optimized baseline proportional derivative controller used in this work showed little effort was needed to obtain better performance using sliding mode. In addition, a colored noise filter was developed to compensate for the fuel sloshing disturbance and incorporated into the sliding mode controller for greater performance increase at the expense of requiring a little more control effort.

PAGE 17

1 CHAPTER 1 INTRODUCTION This chapter provides the motivation and a literature review for this work. Additional literature will be presented in subsequent chapters as appropriate. Formation Flying As scientists aim for higher goals of discovery, satellites prove to be an increasingly important tool. Bigger, more powerful satellites can search the cosmos more easily for distant planets, map the Earth, or even try to discover the origins of the universe itself. However, when a satellite reaches a certain size it effectively can no longer be launched into space. The only answer to bigger, more powerful satellites is either construction in space or formation flying. Construction is space is difficult and costly, but perhaps one day when the international space station (or its successor) is finished this will be a viable solution. Formation flying is the only near term solution. Thus, the idea of not one big satellite but rather many little satellites, each acting as an individual element of a much larger sensor array, is actively being pursued by National Aeronautical and Space Administration (NASA) and other agencies. These missions greatly improve the sensing capabilities available to NASA, while keeping the associated cost and risk down. One of the main problems with formations and individual spacecraft in orbit is the amount of fuel consumed. Small differences in altitude between satellites will result in different velocities and orbits. This will require some satellites to expend more fuel in

PAGE 18

2 order to continually maintain formation. In addition, if they are in different inclinations, they are in different orbital planes (planes that cut through the center of the Earth) and thus will naturally move farther apart or closer together depending on where they are in the orbit. Therefore fuel is required for formation maintenance, which inevitably leaves fuel tanks partially filled. Partially filled tanks can lead to fuel slosh which can degrade the pointing accuracy and therefore affect the science mission. Fuel slosh compensation is not an easy task and unfortunately is often neglected. Satellites do not typically change their orbits unsupervised, but with formations there will have to be some level of autonomy where the satellites must make corrections to their orbits. Satellites will be required to change orbit to close gaps in formations or to avoid collisions from orbital debris. Gaps in the formation will occur as satellites run out of fuel or have failures. The control solution must be robust to account for uncertainty, produce stable motion, be easy to understand, and easy to implement. There has been much research on satellite formation control [1-8]; however these works assume satellites with robust individual control systems. General Problems Controlling Single Satellites To make a formation more stable, we must improve individual satellite control. This section will provide a brief overview of some of the problems with satellite control. Later greater detail will be presented on two specific problems: fuel slosh and solar snap. External Disturbances NASA engineers must contend with many types of disturbances when controlling a spacecraft. Low altitude satellites must fight gravitational and aerodynamic torque which can deteriorate their attitude accuracy and over time degrade their orbit. These

PAGE 19

3 gravitational torques can be present on large satellites where gravity has more of an effect on one end of a satellite than the other. In deep space or geosynchronous orbits, solar based disturbance becomes the dominate environmental disturbance since there is no friction or dampening. Also in these orbits, the gravity effects become negligible compared to effects produced by the sun. One of the major solar disturbances is solar snap which is thermally induced vibrations. V ibrations due to thermal changes have been known for a long time, but here on Earth there are few situations where these vibrations are of any interest. Another disturbance due to the sun is called radiation torque. Photons emitted from the sun slam into the satellite (a force) which produces a torque about the center of mass of the craft. This radiation torque grows and wanes with solar flare activity, but none the less pushes and prods satellites until they slip behind the dark side of a planet or moon. Parametric Uncertainty Control systems are typically designed using mathematical models of the system to be controlled. Thus there is a direct correlation between the performance of a standard classical controller and the accuracy of the mathematical model of the plant. Unfortunately in satellite control, the inertial matrix is only known within 5% of its in orbit value. This model inaccuracy must be addressed in order to achieve the higher level of performance expected from the missions. Thus engineers are currently examining new control schemes which help to combat this problem and improve performance in the presence of modelling error. Another problem that effects the inertia matrix is fuel slosh and fuel expenditure. Fuel can compose 40% of the total mass of a satellite prior to launch [9]. During orbit insertion,

PAGE 20

4 some of this fuel will be burned to perform various v maneuvers to achieve a desired orbit. However, some of the fuel will be kept so that later orbital corrections can be made during the life span of the spacecraft. This remaining fuel poses a problem as it moves around inside of the satellite and thus changes the inertial matrix on orbit. Literature Review This section will provide a review of spacecraft that have had problems with fuel slosh and solar snap. The controllers developed to compensate for these disturbances will also be discussed. An analytical description of these disturbances will be provided in a later chapter. Fuel Slosh Modelling fuel slosh dynamics Accurate analysis of the fuel requires the use of the Naiver-Stokes equations, which are extremely computationally intensive. Since most satellites have a limited amount of processor power assigned to control, this is not practical for controller implementation. Therefore the most used models are the simplest ones such as the pendulum, mass-springdampener, and moving fuel tank walls [11]. For small motions, the objective is to reproduce the sloshing mode resonance frequencies. Typically a mass is attached to the wall of the fuel tank by a spring. The movement of this mass accounts for the free surface wave motion and the influence on the spacecraft’s dynamics is represented by the force the spring exerts on the wall of the fuel tank. Dampeners are commonly introduced to represent the viscous dampening of the fuel. There is typically one set of mass-springdampeners for each sloshing mode that is being represented. Swirl effects of the fuel can also be included as reaction wheels with torsional springs attached. These models however

PAGE 21

5 are only accurate in representing the sloshing modes when small linear or angular motions are performed. For accuracy during larger motions, better models must be utilized. However these models are not ideal for control purposes. The typical pendulum model used for fuel slosh is derived from work done with rockets. Typically in this scenario, the rocket’s main engine is firing and produces a local acceleration that pools the fuel at one end of the tank. The amount of fuel that remains pinned to this end of the tank is the stationary fuel mass fraction. The rest of the fuel is able to oscillate about the center line of the fuel tank in a pendulum-like motion. These oscillations are typically small since the fuel has only a small momentum compared to the forces produced by the local gravity. Unfortunately this model is not as useful in the situation of satellite attitude control, since there is no external force produced by a thruster that is constantly firing throughout the life of the vehicle. Thus a model is needed in which the fuel has more freedom of movement, can exhibit a larger range of movement, and produce accurate results. The work developed by Agrawal [10] for use on the INTELSAT IV, dual spin spacecraft, is a boundary layer model of the fuel for a spinning spacecraft. The author reached the conclusion that the analytical model gives the most accurate prediction of the fuel dynamics. The finite element model consisted of 81 nodes which represented a spherical fuel tank, shown in Figure 1-1. This fuel tank had a radius of 0.42 m, was 1.31 m from center of mass of the craft, and spun at 30 rpm. The solution for the fluid motion is obtained by solving three equations: an inviscous fluid problem, a boundary layer problem, and a viscous correction problem. The inviscid and viscous correction solutions are obtained by using finite element methods, while the boundary layer problem is solved analytically. The numerical results showed two slosh modes which could be accurately

PAGE 22

6 modeled by the pendulum. However there were also first azimuth and elevation modes, lower modes to be inertial modes and higher modes to be a combination of inertial and slosh modes. These other modes would not be accurately modeled by a pendulum. Controlling fuel slosh Fuel slosh is an example of under-actuated control. The objective is to control both the rigid body dynamics of the spacecraft and the fluid dynamics of the fuel only having access to effectors that act directly on the spacecraft. The fuel’s dynamics are controlled through the coupling of the two. One of the main problems with controlling this type of disturbance is that one can not measure the position, orientation, etc. of fuel in a satellite. One can only measure the effects the fuel slosh on the total system. The states and parameters of the fuel must be estimated, and a controller must try to account for the fuel’s influence. Therefore it is difficult to control what can not be measured directly. Thus several types of passive methods are used to control slosh such as baffles [12], slosh Figure 1-1. Finite element analysis of a fuel tank inside of a satellite.

PAGE 23

7 absorbers [13] to dissipate energy through sloshing, and large fuel tanks are broken up into smaller ones. However these steps add complexity and weight to a satellite’s design. The added complexity increases the construction time, the cost, and the chances of mechanical failure of the satellite. The added weight also increases the cost of launching the satellite since it costs approximately $10,000 per pound to put things into space. Fuel slosh directly affects the performance of the attitude control system, in both pointing and tracking maneuvers. One of the simplest ways to handle fuel slosh is to identify the sloshing frequencies and eliminate the influence on the feedback with a notch filter. Unfortunately this denies the satellite the ability to operate in certain frequency ranges since the filter attenuates the slosh frequencies in the feedback. Hence, normal feedback without sloshing effects will also be attenuated. Also as fuel is expended and the sloshing frequency changes, the filter would need to be updated to notch out the new sloshing frequencies again. Cho and McClaroch [14] developed a controller which utilized a common 2D pendulum model. Their equations only dealt with motion in a plane, but “a significant extension of what has been done previously, since it simultaneously controls both the rigid vehicle motion and the fuel slosh dynamics.” One of the major problems they set out to solve, which others were unsuccessful in solving, was the transverse motion of the craft. Their method utilized a Lyapunov function approach, and was able to control the system. They utilized additional control effectors to counter this movement (i.e., gas jet thrusters) while also controlling the pitch and suppressing the slosh dynamics. Kuang and Leung [15] use a simple model of a liquid-filled spacecraft subject to external disturbances, and control it with a state feedback controller. The spacecraft’s H

PAGE 24

8 elliptical fuel tank is assumed to be completely filled with an ideal liquid in uniform vortex motion. Because of these assumptions, the slug of liquid can be represented by a finite number of variables. These variables can then be solved using Helmholtz equations. For the simulations, the liquid slug is located at the center of mass of the vehicle and surrounded by a viscous layer. An attitude controller was then developed using the theory of Hardy-infinity optimal control theory. The controller was capable of attenuating the vortices in the fuel tank. The pointing accuracy remained constant regardless of small changes in the viscosity. The authors plan to next look at partially filled fuel tanks, which is a more difficult problem. Thus far we have discussed more complex models to model liquids and said that in order to provide the level of control that a designer wants, one must use a complex model. An alternate example is Grundelius and Bernhardsson [16], who developed a minimumtime optimal controller for moving open rectangular fluid containers down an assembly line. Proper control of the containers is important. If the accelerations applied to the containers are too high, then the contents of the containers may spill out of the open top. This is obviously a loss to the company in material and results in a more expensive process. Also, if the top of the container is to be eventually sealed, the contents may spill on to an area where glue may be applied, thus making for a bad seal. This could result in customers disliking the product because of “cheap” design. When the speed of the system began to exceed two meters per second, fluid began to slosh out of the containers and the controller was not able to achieve the desired result. The authors incorporated a minimumenergy cost function into the design. Using only a simple linear model for the position of the surface of the liquid and the container, a controller that used a quadratic penalty on the H

PAGE 25

9 slosh and a quadratic penalty on the terminal state was able to produce the desired results without a complex hydro-dynamical model. Do we need to worry about fuel slosh? Fuel slosh is associated with liquid fuel which is used by thrusters. These types of thrusters are efficient and can provide a significant amount of control authority. Slosh forces can be eliminated via the use of alternate actuators. Let us take a look at some of the various types of control actuators available to satellites. W isniewski and Blanke [17] developed a three axis magnotorquer controller for satellites subject to a gravity gradient while in a polar orbit about the Earth. Magnotorquers are an attractive form of control for small inexpensive satellites since they use the Earth’s magnetic field to control the orientation of the spacecraft. However magnotorquers produce very small torques, and their principle of operation is nonlinear and difficult to use. They are difficult due to the fact that control torques can only be generated perpendicular to the Earth’s magnetic field. Also the size of the spacecraft might prohibit the use of magnotorquers since they produce low torque. Chen et al. [18] developed an optimal controller which utilizes both thrusters and magnotorquers to unload extra momentum that has accumulated and saturated the reaction wheels of a spacecraft. Their method separates the required torques for the thrusters and magnotorquers using a simple optimal algorithm. Simulations conducted by Chen showed a 20% fuel savings using this method. Liquid [19] apogee motors are common in geosynchronous spacecraft design. This results in liquid fuel constituting almost half of the vehicles’ mass. In space, liquid fuel motion greatly influences the attitude stability and control. The level and nature of these

PAGE 26

10 influences are dependent on the geometry and location of the fuel tanks, the ratio of spacecraft mass to fuel mass, the fill level of the fuel tanks, and the physical properties of the fuel itself. They also depend on the flexible structure of the spacecraft and the behavior of the estimation and control schemes that are present. Unfortunately, accurate prediction of the fuel dynamics is a difficult problem due to the complexity of the hydrodynamical equations of motion. LEASAT, which was launched in 1984, experienced attitude instability during the pre-apogee injection phase which immediately followed the activation of despin control. The instability was the result of the interaction between the lateral sloshing modes and the attitude control. Pulsed plasma thrusters (PPT) are becoming a more popular control solution for small spacecraft [20]. This fact is due to the small electrical consumption and the use of a solid teflon fuel. The fuel is attractive because it is stable, compact, and does not provide complications due to sloshing liquids. EO-1 provided a successful demonstration of an attitude control system using these PPT’s shown in Figure 1-2 and Figure 1-3. PPT’s are capable of very high pulse width modulation frequencies (PPT’s are only capable of onoff or bang-bang type of control). Although PPT’s have some very favorable properties, they still have the liability of using a non-renewable fuel and are only capable of operating Figure 1-2. PPT for EO-1.Figure 1-3. Diagram of a PPT.

PAGE 27

11 on small spacecraft. Another technology that shows promise for satellites ranging from 100 grams to 500 kg is MEMS engines [21]. MEMS (micro-electrical-mechanical systems) is an exciting new technology field which currently is beginning to make its way out of the R&D labs and into commercial products. MEMS are microscopic mechanical-electrical systems that are developed using a combination of micromachining and standard integrated circuit design. Some of the common MEMS components on the market are accelerometers, gyros, and GPS receivers. Here MEMS rocket engines are developed for space applications. The MEMS engines are the size of a human hair and have a thrust-to-weight ratio hundreds of times greater than the space shuttle. Although there are many good alternatives to liquid thrusters for small spacecraft, liquid fueled thrusters still provide the best source of high thrust capability for both small Figure 1-4. NEAR spacecraft orbiting the asteroid Eros.

PAGE 28

12 and large spacecraft. Therefore fuel slosh will always be present in these types of systems and thus affect high precession pointing and tracking accuracy. Missions involving fuel slosh disturbances We now understand why fuel slosh exists. Next let us look at missions where slosh has been observed. The Near Earth Asteroid Rendezvous (NEAR) spacecraft was sent to orbit and observe the Eros asteroid on December 23, 1998. The spacecraft, shown in Figure 1-4, had to fly by the asteroid following an unsuccessful firing of its main engine a few days earlier. A subsequent successful firing of the engine put NEAR on course to rendezvous with Eros to begin its planned year long orbital mission starting in mid-February 2000. The reason the main engine firing did not work is the on board control system terminated the burn when sensors registered higher than expected lateral accelerations. These accelerations were later determined to be caused by the spacecraft’s fuel. The NEAR’s mission had to be delayed for over a year because the control system was incapable of handling the coupled dynamics. Situations like NEAR do not just happen around asteroids, but with the advent of the International Space Station similar problems could arise. As resupply ships travel to and from the station, similar problems could occur and a ship could crash into the station. This scenario is similar to the accident that occurred with MIR and one of its resupply vessels. Although that incident was not caused by fuel slosh, it still does represent a potential problem. Thermally Induced Vibrations: Solar Snap Thermally induced vibration or solar snap occurs when dynamics are excited by a rapid temperature change. This type of vibration has been known for some time but rarely

PAGE 29

13 appears in applications on Earth. However in the harsh, irradiated vacuum of space, without the presence of an atmosphere to dampen the vibrations, this disturbance can adversely affect the performance of high precision instruments and equipment. This section, like the proceeding one, will present missions that where adversely affected by solar snap and control systems that were designed to compensate for it. Controlling thermal vibrations The Hubble Space Telescope (HST), shown in Figure 1-5, was launched on 24 April 1990 with the hopes of making numerous new scientific discoveries. The HST was a 13 ton, free-flying spacecraft with a pointing precision requirement of 0.007 arc-sec over a 24 hour period (which is the most stringent requirement ever imposed). During its initial checkout period, a pointing “jitter” was discover. The jitter was produced by a thermally induced vibrations imitating from the solar panel array. Later in December 1993 the solar Figure 1-5. Hubble Space Telescope diagram.

PAGE 30

14 arrays were replaced by the Space Shuttle Endeavour in an attempt to reduce the jitter. The solar panels exhibited an end-to-end bending oscillation when the spacecraft transitioned between sunlight and shadow. They also had a sideways oscillation on the day side of the Earth. At its worst, the 20 ft. solar arrays would deflect as much as 3 ft. This would have the effect of rendering any long time exposures of over 25 min. pointless. The original controller was a digital PID controller with an FIR filter in the rate path to attenuate highfrequency, main-body bending modes. This design was incapable of handling the induced disturbances from the solar arrays, which caused a peak pointing jitter of 0.1 arc-sec. ([22], p. 583) W ie and Liu [25] redesigned the attitude control system using an controller. This controller was designed to attenuate the two major modes from the solar panels which operated at 0.12 and 0.66 Hz. Through a trial-and-error process, they were able to tune the weights of the controller and properly achieve the mission requirements. Singhose et al. [23], while working at M.I.T., designed a method to minimize vibrations in flexible spacecraft which utilized command profiles. The commands are designed by the following: 1.Formulate the equations of motion for a system subject to a sequence of impulses. 2.The equations are solved for the sequence of impulses that results in the smallest amount of residual vibrations in the system. 3.The sequence is convolved with the desired input to generate the command sequence which is issued to the real system. The command sequence is a pulse width modulation (i.e., on-off sequence) of commands for the spacecraft’s thrusters. Later at Georgia Tech., Ooten and Singhose [24] extended the work to sliding-mode control. The drawback for this work is, although flexible spaceH

PAGE 31

15 craft are mentioned, all simulations involve only a mass-spring-damper system. Also, solving the equations of motion (step 2) for a nonlinear, flexible satellite would be very difficult. Missions which involved solar snap Originally spacecraft were very simple, but as technology progressed into the 1960s more appendages were being put on. These structures offered new challenges in controlling flexible structures and vibrations. Below is a table of spacecraft that had their performance compromised by thermal vibrations. These satellites do not include military and Soviet spacecraft. T able 1-1: Spacecraft with thermally induced disturbances a Spacecraft Date Description Ulysses19901991 Spin-stabilized solar probe experiences unexpected nutation due to thermally induced vibrations of beryllium-copper axial boom Upper Atmosphere Research Satellite 1991 “Thermal snap” at orbital night/day transition attributed to rapid heating of large single solar array Hubble Space Telescope 1990Jitter phenomena attributed to thermally induced vibrations of FRUSA-type solar arrays LANDSAT-4/5 1980s Thermal elastic shock experienced during orbital night/day transition due to disturbance of large single solar array. Communications Technology Satellite 1978Three axis stabilized satellite experienced “thermal elastic shock” during orbital night/day transition V oyager 1977 Low frequency oscillations of PRA booms during LEO operations attributed to “thermal flutter” Apollo 151971Thermally induced vibrations of twin Bi-STEM booms filmed during 64th lunar orbit by astronauts. Explorer 45 1971 Unexpected nutation of spin stabilized satellite attributed to thermal bending of four radial booms. NRL 161, 163, and 1641969Thermally induced vibrations of STEM booms used for gravity-gradient stabilization leads to large satellite motions. OGO-IV 1967 NASA Orbiting Geophysical Observatory Satellite STEM booms experience strong thermally induced vibrations during orbital night/ day transition. OV-101966US Air Force gravity-gradient satellite using STEM booms equipped with tip masses experiences thermal flutter and complete inversion in orientation.

PAGE 32

16 “In spite of the practical importance of thermally induced vibrations to the successful operation of spacecraft, recurrent vibration problems from the 1960s to the present day suggest that the phenomena is either underemphasized or not well understood in spacecraft design” (Thorten [26], p. 344) Parametric Uncertainty A very simple way to model the effects of fuel and thermal vibrations is to represent them as external disturbances (i.e., deterministic noise). Depending on the mass fuel ratio, fuel slosh affects the mass distribution of the vehicle which affects the inertia matrix. Since control systems are model based, this will directly affect the level of control. If a control system could be designed that provided generic control, or model independent control, then the effects of these disturbances could be negated. This idea of generic control was the subject of work done by Walchko and Mason [27] utilizing fuzzy logic. The idea was that fuzzy logic is rule based and not model based control. Therefore fuzzy logic would be capable of providing better control in the presence of uncertainty. Although several different configurations of fuzzy logic were developed (PD, PID, and sliding mode) none produced the desired results when there was a significant amount of uncertainty. In fact it was shown that the two best were sliding mode Explorer XX 1964 NASA spin-stabilized satellite experiences spin decay due to interaction between STEM boom thermal bending and solar radiation pressure. 1964-83D and 196322A 1964/ 1963 Dynamic thermal bending of STEM booms on APL gravity-gradient satellites. Aloutte 1 1962 Canadian spin-stabilized satellite experiences spin-rate decay due to thermal structural response of STEM booms similar to Explorer XX. a.Information contained in this table comes from Thorten [26] on page 345. T able 1-1: Spacecraft with thermally induced disturbances a Spacecraft Date Description Continued

PAGE 33

17 and fuzzy sliding mode. However the non-fuzzy sliding mode had the added advantage that stability could be shown for the controller, whereas fuzzy logic lacks any general method to prove stability Petroff et al. [28] attempted to address this problem of fuzzy logic stability. They looked at simple linear single input single output systems and were able to develop a stability analysis. However, this method did not naturally extend to multi-input multioutput systems and definitely not to nonlinear systems. Thus although fuzzy logic appeared to be attractive as a model independent control architecture, it was extremely weak with respect to stability analysis. This fatal flaw ultimately disqualified it for our applications. Summary and Overview of the Thesis What is needed are controllers that are easy to adapt from one satellite to another but still retain their robustness and efficiency. Problems due to simple control systems being implemented on satellites are numerous. This usually results in a redesign of the control system which is costly and a waste of time and resources. Now couple these problems with the difficulties of formation flying and things get very risky. This work will primarily focus on fuel slosh and solar snap. Both of these phenomena are interesting in that they have rarely been observed. Rather, engineers have examined the data from satellites that have experienced unexplained problems and deduced these were the causes. Neither of these problems is generally looked at when a control system is developed, but rather they are dealt with only when they appear. This attitude proves to be a problem when applied to formation control. The possibility of one satellite going

PAGE 34

18 unstable due to these phenomena then introduces the possibility of that satellite affecting the entire formation, thus endangering the mission. Chapter 2 will provide a more in depth discussion of fuel slosh. Various models will be described that are used to model fuel slosh for simulation purposes and models that are used in controllers. Chapter 3 will cover the interactions of solar panels with satellite dynamics. Thermally induced vibrations are described by the solution to their partial differential equations and modal equations. Chapter 4 will cover satellite attitude dynamics and control. This chapter will introduce the reader to the fundamental dynamics of satellites. Chapter 5 will introduce the basics of sliding mode control and other attitude control algorithms and control hardware. The specific controller for this work will also be derived and discussed. Chapter 6 will provide an overview of the simulations. All parameters, equations, and assumptions will be presented here. Also all results will be discussed in this chapter. Chapter 8 will provide conclusion to the work done.

PAGE 35

19 CHAPTER 2 CONTROLLING THE UNKNOWN: FUEL SLOSH The introduction covered much of the problems associated with fuel slosh and highlighted missions where this phenomena occurred. This chapter will expand on the introduction by focusing more on the fuel slosh disturbance. Two models were developed, one for use in the simulation to provide the disturbance torques on the satellite and another model to help the controller account for the fuel slosh disturbances. Introduction The sloshing behavior of liquids and its modeling has been studied in many different areas: fuel slosh in aerospace applications, movement of fluid filled containers in industrial/manufacturing applications, earth quakes, vehicle and ship dynamics. Unfortunately the choice of a model is not a simple one. Most analytic models that try to accurately represent the dynamics are three dimensional partial differential equations. These models are heavily dependant on boundary conditions and computationally expensive to solve. Thus they are typically not used in controller design. However, the nonlinear dynamics that appear from these complex models are important when designers wish to move large amounts of fluid rapidly. The aerospace industry has looked at controlling the effects of fuel movement in aircraft fuel tanks for years. Most of the ideas (i.e., baffled fuel tanks, breaking one big tank up into many smaller ones, etc.) can be traced back to standard aircraft design.

PAGE 36

20 The interaction of the fuel and the spacecraft is more complex in three-axis stabilized systems than spin-stabilized systems [29]. This is due to the fact that this system has constraints that reduce the number of degrees of freedom, simplify the model structure, and is numerically easier to compute. However, spacecraft will typically operate in several modes of operation over their life span and thus it is important to be able to control fuel slosh in any mode of operation. Modelling Fuel Slosh Dynamics Models for fuel slosh have not changed much since Abramson started in 1961 [9]. Models over the years have been: (single and multi) mass-spring-damper, pendulum, liquid slug, and CFD/FEA models. Typically, the controller uses a model to account for the fuel present inside of the spacecraft. Unfortunately, the models listed are too simple or based on a static representation of the fuel and do not attempt to modify the parameters of the model during run time. The coordinate system for the fuel is typically the body coordinate system [29]. This has the result of making the spacecraft’s motions appear as a body-force field in the Naiver-Stokes equations. The boundary conditions at the wall also become stationary with only the free surface remaining transient while the viscous forces provide a dampening of the liquid. The fuel motion is influenced by four classes of forces: gravity, inertia, viscous, and capillary forces. The early work with slosh dynamics and spacecraft centered around liquid filled rockets and their stability [29]. This topic is still of interest and the assumptions [9] for this problem typically are: • The fuel has small displacements, velocities, and slopes of the liquid-free surface.

PAGE 37

21 • The fuel tank is rigid. • The fuel is assumed to be a nonviscous liquid. • The fuel is also modelled as an incompressible, homogenous fluid. Propellant slosh is the induced motion of liquid due to acceleration of the container. This motion along with its reactive forces can deteriorate the pointing performance. This is especially true for systems with large fuel/weight ratio. Due to the nonlinear dynamics of sloshing, it is not always possible to compensate for its affects using a standard attitude controller. T ypically a pendulum system, shown in Figure 2-1, is used to model a single sloshing mode. Since slosh dynamics typically display several modes, several pendulums should be included in the model of the system. Thus the more dominate modes are modelled, the more realistic the simulation. A problem with this type of model is determining the values of the pendulum lengths, masses, springs, and dampeners. One way to do this is to look at actual flight data where you can determine the modal frequencies. Assuming the amount of mass in the slosh mode, the value of the spring constant can be solved. The damper value is picked so that the system’s oscillations die out after an appropriate amount of time. The main problem with liquid slosh dynamics is estimating the hydrodynamic pressure distribution, forces, and moments. One reason that this is so difficult is the dynamic boundary conditions at the free surface varies with time in a manor not known a priori, as depicted in Figure 2-2. Hydrodynamic pressure in many rigid tanks is composed of two parts. The first part is the fluid that moves with the tank in unison. The second component is the sloshing at the free surface. This component is typically modelled as a

PAGE 38

22 mass-spring-damper or pendulum. Suggested values for the pendulum shown in Figure 21 are given by EL-Sayad et al. [30] as [2-1] where L is the length of the pendulum, R is the radius of the tank, and mp is the mass at the end of the pendulum. L R 1.84 ---------1.84 h R --coth = mpMtotalR 2.2 h ---------1.84 h R --tanh = Figure 2-1. Traditional fuel slosh pendulum model where M is the total fuel mass center, mp is the pendulum mass, ms is the stationary mass fraction, L is the pendulum arm length, and g is the local acceleration due to some external force. Figure 2-2. Propagation of fuel slosh in an elliptical tank over time. These images demonstrate the difficulty in modelling the fluid movement with simple representations. g L mp M ms

PAGE 39

23Fuel Slosh Models T wo different models for the fuel will be utilized in this work. The first model will be included in the satellite dynamics and provide the proper disturbance torques on the system. The second model will be computationally simpler, but still contain the important aspects of the fuel’s disturbance on the system. Simulation Model The model utilized is based on the mass-spring-damper models by Sidi [9] and Hughes [31] shown in Figure 2-3. Each sloshing mode will be composed of a mass, a spring, and a damper. The mass represents the mass of fuel in a specific mode. The spring represents the force exerted on the tank wall by a sloshing mode. The damper represents the baffling in the fuel tank that dissipates the sloshing movement. The actual dynamics of this disturbance will be derived later when the dynamics of the satellite are covered. This is due to the fact that both the dynamics of the spacecraft and fuel are tightly coupled and it is easier to present them all at once. Figure 2-3. Diagram of a single sloshing mode, with its mass attached to the fuel tank walls by springs and dampers aligned along the x, y, and z axes.

PAGE 40

24Controller Model The model used in the controller needs to be capable of representing several sloshing modes, but can not be computationally expensive. Figure 2-4 provides a frequency domain visualization of the proposed model. The model represents the sum of the slosh modes as a band pass filter. This will allow us to incorporate the effects of the multiple sloshing modes without having to specifically identify their individual characteristics. A simple band pass filter is the combination of a low pass and high pass filter in series. The standard transfer function for this filter is as follows: [2-2] [2-3] as as 1 + -------------1 bs 1 + -------------as abs2ab + () s 1 ++ ---------------------------------------------= a 1 freqhighpass---------------------------= b 1 freqlowpass-------------------------= Figure 2-4. Plot of the effect of multiple modes in the velocity of a simply supported beam. This is representative of the fuel slosh disturbance in flight data.

PAGE 41

25where the transfer function on the left is a high-pass filter and the one on the right is a lowpass filter. The variables a and b are the inverse of the high and low pass cutoff frequencies. Putting [2-2] into state space form yields [2-4] [2-5] The interesting thing to note about 2-5 is the output given by the model is velocity and not position. This is due to the derivative term (i.e., s in the numerator) in 2-3. This is a reasonable characteristic, since the fuel slosh effects show up in the rate terms and not the position terms. x x ab + () ab ---------------– 1 ab ----– 10 x x 1 0 u + = y 1 b -0 x x =

PAGE 42

26CHAPTER 3 DYNAMICS OF SOLAR ARRAYS As discussed in the introduction, thermally induced vibrations are a source of problems for satellite attitude control. This chapter will provide the theoretical development of thermally induce vibrations in long, thin beams. These results will finally be used to generate the solar snap disturbance in the simulations. Introduction Research on thermally induced vibrations began back in the mid-1950’s, before Sputnik orbited the Earth in October 1957. However it was not until the 1960’s that applications for this new idea appeared. Initially spacecraft appendages (i.e., solar panels, communication arrays, booms, etc.) where thought to be simple systems like those used by Sputnik (Figure 3-1), but slowly it became apparent that they were much more complicated. Strange and unexplained spacecraft behavior was now being attributed to this new nonlinear behavior being observed. Thermally Induced Vibrations First, this section will derive and solve the partial differential equations associated with solar snap. Then differential equations will be obtained from the PDEs which are utilized in the simulations. The derivation that follows will parallel work done by Thorten [26] and Rietz [67]. However, there was great difficulty in actually following their work. Fortunately, all key equations were independently derived. The results in this work were also similar to

PAGE 43

27derivations performed on cantilever beams excited by sources other than thermal vibrations by Greenburg [32], Kreyszig [33], Meirovitch [34], and Zill and Cullen [35]. Heat Flux Thermal vibrations are excited by a sudden change in heat flux to the object. In this work, heat flux changes greatly during the day-to-night and night-to-day transitions. Here the satellite passes into and out of the Earth’s shadow. The Earth’s shadow is actually composed of two regions shown in Figure 3-2. The penumbra is a region of partial shadow and the umbra is a region of full shadow. When a satellite is in low Earth orbit (LEO), the penumbra is small and the satellite essentially transitions from the day side to the night side (i.e., umbra). This sudden change from dayto-night or night-to-day results in larger thermal vibrations. This is due to the sudden change in heat flux, much like turning on or off a light switch. The amount of heat flux ( q ) seen by a satellite is composed of solar flux from the sun ( ), Earth emitted radiation flux1 (), and Earth’s reflected heat flux () from the sun. qsqeqaFigure 3-1. The Soviet satellite Sputnik which orbited the Earth in October 1957.

PAGE 44

28This amount changes with the satellite’s position in orbit and its orientation relative to the solar flux vector (see Figure 3-3). [3-1] [3-2] where is the surface absorptivity for solar radiation and is the angle between the surface normal and the solar flux vector. Summing Forces The solar panel is modelled as an Euler-Bernoulli cantilever beam subjected to a thermal moment. The beam, a cross section is shown in Figure 3-4, will have the following assumptions: • Plane cross sections before bending will remain planes after bending. 1.The radiation referred to here is black body radiation.qqsqeqa++ = qs1350 as () cos = as Figure 3-2. Diagram of the penumbra and umbra regions which block a satellite from sunlight. Figure 3-3. Another view showing the path of a satellite and how its orbit does effect the heat flux it is exposed to.

PAGE 45

29• Lateral displacements in v and w are due to bending only. Deformations from shearing and changes in transverse beam dimensions are negligibly small. • The beam is a Euler-Bernoulli cantilever beam. • The system is uncoupled, meaning that changes in beam orientation do not effect the amount of heat flux entering the system. • The temperature gradient is only one dimensional (i.e., through the thickness of the beam). The moments on the system are defined as [3-3] [3-4] Linear strain displacement is as follows: [3-5] Hooke’s uniaxial law is [3-6] Then substituting [3-5] into [3-6] yields: [3-7] Now substituting [3-7] into [3-3] produces [3-8] where MzxyA dA= MTE TyA dA= xx u x d d u0y x2 2 v – z x2 2 w – == xE xE T – = xE x d d u0y x2 2 v – z x2 2 w – E T – = MzE x d d u0y x2 2 v – z x2 2 w – E T – yA dAEI x2 2 v – EIyzx2 2 w – MT– ==

PAGE 46

30 and [3-9] Next, summing the forces in the y direction in Figure 3-5 results in [3-10] Note that from here on out, the terms in parenthesis will be dropped unless it is felt that they are necessary for clarity. [3-11] [3-12] Summing the moments about the center of mass of the beam yields: [3-13] Izy2A dA= Iyzy2A dA= Fyfx () – Vxt () – Vxt () x Vxt () dxmu xt () +++ 0 == A t2 2 u dx x V dx + fdx = A t2 2 u x V + f = MM x M dx 1 2 -Vdx 1 2 -V x V dx + dxM – +++ 0 == Figure 3-4. Cross section of solar panel with coordinate systems superimposed. Note that the S’s in the figures represent ‘s and M’s are moments.

PAGE 47

31[3-14] [3-15] Substituting [3-8] into [3-5] and ignoring the middle term (since we are only looking at motion in the u direction and not the w direction). [3-16] Now taking [3-12] and substituting [3-16] gives us: [3-17] [3-18] This equation can be further simplified by assuming that no uniform load ( f(x) ) on the solar panel exists and realizing the thermal moment is not a function of x thus its second derivative is zero. [3-19] Incorporating these, [3-18] becomes [3-20] where a subscript of x or t refers to differentiation with respect to that variable. The boundary conditions of the system are as follows: [3-21] [3-22] x M V + 0 = x M – V = x EI x2 u2MT+ V = A t2 2d d u x x EI x2 2 u MT+ + f = A t2 2d d u EI x4 4 u x2 2 MT++ f = f 0 = MTxx0 = 4EI A -----= utt4uxxxx+ 0 = uxt () 0 = ux0 t () 0 =

PAGE 48

32[3-23] [3-24] T able 3-1: Solar panel thermal material properties. V ariable Description Units k thermal conductivity q heat flux density c heat capacity thermal diffusivity coefficient of thermal expansion MzL () EIuxx– MT– 0 == VLt () EIuxxxMTx+ 0 == W mK -------W m2K ---------kg m3-----J kgK --------m2s -----1 K --Figure 3-5. Free body diagram of a section of the beam of size dx where M(x,t) is a moment, V(x,t) is a shear force, f(x) is a force per unit length, and u(x,t) is the height of the mass element at position x at time t. M x M dx + V x V dx + f(x) V M u x y

PAGE 49

33Nonhomogeneous Boundary Conditions Unfortunately the boundary conditions for the system are nonhomogeneous and in order to solve the equations this must be corrected. Getting rid of the nonhomogeneous boundary conditions will be accomplished by changing the dependent variable. [3-25] Substituting [3-25] into [3-23] yields: [3-26] The second derivative of u must be equal to zero to make the system homogenous, thus must negate the right hand side of the equation. [3-27] Now we must solve for thus integrating [3-27] twice results in an equation with two constants that can be solved by looking at the initial conditions. [3-28] [3-29] [3-30] Thus [3-25] can be rewritten as [3-31] Our new beam equation and boundary conditions using [3-31] are shown below. Notice how this is now a simpler homogeneous equation to solve. uxt () vxt () x () + = uxxvxxxx+ MTEI ------– == xxMTEI ------– = x () MT2 EI -------– x2C1xC2++ = 0 () 0 C2== x0 () 0 C1== x () MT2 EI -------x2– = uxt () vxt () MTx22 EI -----------– =

PAGE 50

34[3-32] [3-33] [3-34] [3-35] [3-36] Separation of Variables Now [3-32] with its homogeneous boundary conditions can be solved via separation of variables. We will assume that the function can be broken up into two separate functions. Each of these two functions are dependent on only one variable, either position ( x ) or time ( t ). [3-37] Substituting [3-37] into [3-32] gives us an equation which can be separated. [3-38] [3-39] [3-40] Since each side of the equation is dependent on a single variable, it must be equal to a constant. We will choose this constant to be which is squared to simplify the derivation. Each of the two sides can now be solved, which gives us one differential equation in x and one differential equation in t [3-41] vtt4vxxxx+ 0 = v 0 t () 0 = vx0 t () 0 = vxxLt () 0 = vxxxLt () 0 = vxt () vxt () Fx () Gt () = FxxxxG 4FGtt+ 0 = GttG -----– Fxxxx4F -----------2== 4EI A -----= 2Gtt2G + 0 =

PAGE 51

35[3-42] Solving the position function F(x) We will begin by solving the F(x) equation first. Assuming that the solution to this differential equation is an exponential, the following solution is found. [3-43] [3-44] [3-45] [3-46] [3-47] Combining the last two equations into a matrix give us: [3-48] The determinate of this matrix gives us the characteristic equation from which we can solve for the values of that result in zero. The characteristic equation is shown below and a plot (shown in Figure 3-6) which can be used to graphically solve the equation. [3-49] It is now possible to solve the equations in the matrix as a ratio Thus it is possible to solve for 3 of the constants in terms of the fourth. T able 3-2: First four eigen modes of the beam n 1 2 3 4 1.8751044.6940917.85475710.995541 Fxxxx24F – 0 = Fx () B1nx () cosh B2nx () sinh B3nx () cos B4nx () sin +++ = F 0 () B1B3+ 0 == Fx0 () B2B4+ 0 == FxxL () B1nL () cosh nL () cos + () B2nL () sinh nL () sin + () + 0 == FxxxL () B1nL () sinh nL () sin – () B2nL () cosh nL () cos + () + 0 == nL () cosh nL () cos + nL () sinh nL () sin + nL () sinh nL () sin – nL () cosh nL () cos + B1B20 = nL nL () nL () cos cosh1 + 0 = nL B2B1

PAGE 52

36[3-50] [3-51] This represents the nth spatial solution or mode shape. In order to calculate the complete response of the beam, a summation over all modes must be done. [3-52] The coefficient is found by the orthogonality property of the modes. This coefficient is what defines the amplitude of the vibrations while only defines the shape of the vibrating modes. Now because it is assumed that the eigenvectors form an orthonormal basis, the dot product of any one eigenvector with another is zero. For example, if we take the above equation and dot both sides with the first eigenvector only the first eigenvector Fnx () nx () cosh nx () cos nnx () sinh nx () sin – () + – = nnL ()2L2---------------EI ----= n 4 A n 2EI -------------= nB2B1----nL () sin nL () sinh – nL () cos nL () cosh + ------------------------------------------------------== Fx () anFnx ()n 1 = = anFnx () Figure 3-6. Plot of [3-49]. Intersections of the two lines are solutions of the characteristic equation.

PAGE 53

37and its coefficient will remain. We are now free to assume that the solution to the differential equation is a form of the first eigenvector and solve for the unknown coefficient. [3-53] where the thermal moment () comes from Johnson and Thorten which is similar to Reitz [66] shown below respectively. [3-54] [3-55] The in Thorten’s equation is the temperature difference through the thickness of the solar panel. An equation for assuming the solar panel can be modelled as a thin plate, will be derived later in the chapter. The temperature difference equation is dependent on the physical shape of the solar panel, while the thermal moment equation is more general. A quick numerical example, to calculate the coefficients for the first mode, set to 1.875104 and L to a position on the beam ( w ). [3-56] Then solve for the coefficient and finally for F(x) Although the derivations show that the summation symbol goes from zero to infinity, in actuality you would only sum over a finite number of modes. anMt0 () x 22 EI ----------------------Fnx () x d0 LFn 2x () x d0 L-------------------------------------------------= Mtt () MTt () E wh22 ----------------T = MTt () 48 EIqo 4 --------------------496 ----en22 t h2---------------–n4-----------------n 1 = – = T T nL 11.875104 ()2w2----------------------------EI ----= 11L L -------1.875104 w ---------------------== 11.875104 () sin1.875104 () sinh – 1.875104 () cos1.875104 () cosh + --------------------------------------------------------------------------------= an

PAGE 54

38Solving the time function G(t) The solution of the time part of the separation of variables method is a little complicated. The change of the dependent variable is what complicates finding the result. Our solution for the time variable has the form shown below. [3-57] The solution of [3-57] will be in the form shown below, where the term is the particular solution which can be determined if we have an equation for the thermal moment. The general form of the solution is shown below. [3-58] Summary of results We have taken a long, winding road in order to solve the equations of motion for the beam. The final solution of the beam, which describes the motion in one direction, will need to be numerically solved is shown below. [3-59] [3-60] [3-61] [3-62] [3-63] Gttt () 2Gt () + MTtt– = GpGt () A nt () sin B nt () cos Gp++ = uxt () vxt () x () + Fx () Gt () Mtt () x22 EI -----------------– == MTt () E wh22 ----------------T = Fx () annx () cosh nx () cos – nnx () sinh nx () sin – () + []= nnL ()2L2---------------EI ----= nnL L -------= nnL () sin nL () sinh – nL () cos nL () cosh + ------------------------------------------------------= Gt () A nt () sin B nt () cos Gp++ =

PAGE 55

39 [3-64] Again, summation occurs only over the number of modes of interest. Thermally Induced Dynamics Modal Equations The preceding derivation results in the path of the beam at some point x over time. However for the simulation modal equations are developed which are ordinary differential equations (ODE). This section follows the ideas in McConnell [36] and Thompson [37] and will cover the formulation of the modal equations. Further explanation of the theory pertaining to modal equations can be found in Appendix B or Tongue [71]. Each mode of the system will be a second order equation as shown below. [3-65] The subscript n denotes that the equation is for the nth mode. The q ’s are the generalized coordinates of the mode and Q is the generalized force applied to that mode. Note that in this formulation, the generalized force for the second mode will only effect the second mode and have no effect on the first mode or any other mode. This idea is due to orthogonality and as in the previous derivation, the modes are assumed to be orthogonal and do not interact with each other. The relationship between our physical parameter of distance x and the generalized coordinate q is as follows: [3-66] anMt0 () x 22 EI ----------------------Fnx () x d0 LFn 2x () x d0 L-------------------------------------------------= q nCnq nn 2qn++ Qn= xFnqn n=

PAGE 56

40where is our position function1 for mode n from the previous derivation, note that again it’s dependence on x is not shown for simplicity. [3-67] where P(x) is some external distributed load per unit length. The mass, spring, and dampener matrices are shown below: [3-68] [3-69] [3-70] The dampening matrix is assumed to be the same magnitude as the mass and spring coefficient matrices. Thus the dampening matrix is just a scaled version of the two where and are scalar weighting coefficients that are tuned until the desired amount of dampening is present in the system. The natural frequencies () are defined by [3-51] and are shown again below: [3-71] Disturbance Torques In order to include the solar snap disturbance torques in a simulation, not all of the preceding steps need to be performed. This section will follow along the work of Johnson and Thorten [68]. Instead of solving [3-25], a differential equation is developed which 1.Thompson refers to as a normal mode and uses the variable see page 438 in [37]Fn FnnQp1 Mp-----Px () Fnx d0 L= Mn AFn 2x d0 L= KnEI x2 2 Fn 2x d0 L= Cnn 2+ = nnnL ()2L2---------------EI ----=

PAGE 57

41draws on the material already presented. This solution is easier to understand and integrate into a computer simulation. The equations Johnson and Thorten give for the coupled satellite solar panel system using a generalized form of Lagranian equations is as follows: [3-72] An expression for the disturbance torque can be obtained by moving terms associated with movement of the solar panel to the right side of the equation.The disturbance torque () is a composite of two torques. The first torque () is generated from the quasi-static terms, or the nonhomogeneous boundary conditions. The second torque () is derived from the vibrational aspects of the beam. [3-73] This can be seen by substituting the second derivative of [3-25] into the disturbance torque equation obtained from [3-72]. [3-74] [3-75] [3-76] [3-77] Isat ARox + () u xt () x d0 L+ 0 = gssgqsgvgssgqsgv+ = gss ARox + () u xt () x d0 L– = u xt () t2 2d d uxt () t2 2d d vxt () xt () + () == gqs ARox + () t2 2d d xt () x d0 L– A2h T4 I ------------------RoL33 ----------L44 ----+ T – == gv ARox + () t2 2d d vxt () x d0 L– ARox + () nx () x d0 Lqnn 1 = N==

PAGE 58

42where is the shape function (which use to be but is changed to match Johnson and Thorten’s work) and is the generalized coordinate for mode n and is the temperature difference through the thickness ( h ) of the solar panel. Thickness Temperature Gradient The equation for the temperature in a thin plate of thickness h can be found in Thorten’s book [26] on pg. 86 or Incropera and DeWitt [69]. A diagram of the solar panel cross section is shown in Figure 3-7, the equation for the temperature change through the thickness is given below. [3-78] The temperature difference through the thickness of the beam can be found by subtracting the temperature at the heated surface from the temperature at the insulated surface at time t After subtracting the two equations and some algebra, the temperature gradient at time t steady state at time and its acceleration becomes: nFnqn T Tyt () qoh k -------t h2----1 2 -y h -1 2 -+ 21 6 -– 2 2----1 – ()nn2------------n y h -1 2 -+ en22 t – h2-------------------cosn 1 = – + = Figure 3-7. A cross section of the beam in which a heat flux is applied to one side. Notice the beam starts off at a uniform temperature. Figure 3-8. Plot of the temperature change of a beam. Notice that the temperature difference between the two sides remains constant eventually. y

PAGE 59

43[3-79] [3-80] [3-81] Notice that Figure 3-8 and [3-80] show the temperature gradient of the solar panel will eventually reach a steady state. The two sides of the beam never simultaneously achieve the same temperature while there is a heat flux constantly applied. This means there is always a moment on the solar panels, and thus always a deflection of the solar panels. Examples of this are shown in Figure 3-9 and Figure 3-10 where the solar panels of Hubble are deflected due to the constant temperature gradient. Hubble’s solar panels are more complex than what is represented in this work. Hubble’s solar panels incorporate bistems, solar blankets, and spreader bars which are very thin, light weigh structures that are highly susceptible to solar snap. Hubble’s solar panels have been replaced several times Tt () qoh k -------1 2 -2 2----2 – n2----en22 t – h2------------------n 135 … ,,, = + = T () qoh 2 k -------= Tt () d2 T dt2-----------= Figure 3-9. Hubble Space Telescope which shows the solar panels bent due to a constant temperature gradient.

PAGE 60

44over the years, and the current ones are smaller and stiffer than the original design. This is an attempt to reduce the influence of solar snap on Hubble. Even though there is a constant temperature gradient through the panels, this does not mean there is constant vibrations. The vibrations are driven by a change in the heat flux, and thus a change in the temperature gradient. These changes primarily occur during the day-to-night and night-to-day transitions of the satellite around the Earth. Conclusion Thermally induced vibrations or solar snap is a disturbance which occurs when a spacecraft transitions into and out of the Earth’s shadow. A temperature gradient occurs in the long thin solar panel which in turn creates a thermal moment. This moment causes the panel to vibrate until the transients die out. During this time frame, the scientific mission has to be temporarily halted and the spacecraft put into a stable orientation. Figure 3-10. Another picture of the Hubble Space Telescope showing deflection of the solar panels during maintenance by the space shuttle.

PAGE 61

45The theory and math behind this phenomenon is not a mystery and the aim of this chapter was to provide the reader with a complete derivation of this problem for a simple solar panel. Here the assumption was made that the temperature gradient only occurs though the thickness of the beam which reduced the problem to one dimension.

PAGE 62

46CHAPTER 4 SATELLITE ATTITUDE DYNAMICS Satellite attitude control can be difficult due to the various types of disturbances and limited control capabilities of satellites. In order to develop a controller, a good understanding of the dynamics of the system is needed for proper simulations to be conducted. This chapter will give an overview of the dynamics of rigid bodies in rotating reference frames. Reference Frames Several different coordinate systems will be used to develop the equations of motion for a satellite. This section will cover the reference frames used in this work and are shown in Figure 4-1. Earth Centered Inertial The first reference frame is the Earth centered inertial (ECI) frame. This frame is a non-rotating frame relative to the fixed stars1. Another way to say this is, the ECI frame is an inertially fixed frame where Newton’s laws are valid. This frame has its origin located at the center of the Earth with its z-axis along the Earth’s mean axis of rotation. The y-axis and x-axis are pointing to some convenient set of stars. Note that the Earth is allowed to rotate in this fixed frame. 1.Actually the stars move, and thus our ECI frame moves with them. This movement is small enough to be ignored in most cases.

PAGE 63

47Orbital Frame The orbital reference frame or local vertical, local horizontal (LVLH) frame is the reference frame from which all of our spacecraft’s angle will be defined. The orbital frame’s origin is located at the center of mass of the spacecraft. The z-axis points towards the center of the Earth. The y-axis is perpendicular to the orbital plane. The x-axis points in the general direction of travel. When the orbit is circular, then the velocity vector and the x-axis are co-linear. Body Fixed Frame The body reference frame is located with its origin located at the center of mass of the spacecraft. The x, y, and z axis are located along the x, y, and z principle axis of the spacecraft. Fun With Vectors and Rotating Coordinate Systems This section will present some simple vector identities and properties which maybe useful to the reader. Although all of the important steps are shown in the following Figure 4-1. The three frames of reference commonly used in spacecraft dynamics.

PAGE 64

48derivations for spacecraft dynamics, at times some algebraic simplifications are not shown. Conventions used: • All vectors and matrices are bold face. • 1 and 0 represent the identity matrix and the zero matrix respectively. The dynamics for a spacecraft are written with respect to several frames of reference. T wo important reference frames will be the ECI inertially fixed frame and a body fixed frame which is allowed to rotate in space. Differentiating a vector relative to a rotating reference frame is given by: [4-1] abc + () abac + = ab ba – = ab () c abc () = ab ba = abc () bac () cab () – = ab aTb = t AIt d d AB AB + = Figure 4-2. A reference frame in motion relative to an inertially fixed frame of reference. p

PAGE 65

49where is the relative angular rotation rate between the inertially fixed frame I and the rotating body frame B Now for the point p located in the rotating system at position r its velocity is: [4-2] [4-3] where is the velocity of the origin of the rotating system relative to the fixed system and is the velocity of the point in the rotating system. A similar derivation can be shown for the acceleration of the point which results in: [4-4] When dealing with rigid bodies, with the rotating frame representing a body fixed frame, the point p does not move within this frame. This is due to the fact that p is fixed in the body. Thus and the proceeding equations reduce to: [4-5] where the velocity and acceleration are functions of the rotation rate of the rotating frame and its velocity and acceleration. Spacecraft Attitude Attitude refers to the orientation a spacecraft occupies in 3D space. Although there are several different ways to represent this (i.e., euler angles, Gibbs vector, fixed angles, Euler symmetric parameters, etc.), typically space applications utilize quaternions. This section does not attempt to provide the extensive understanding needed to employ quaternions but vIt d d rIt d d RIt d d r 'I+ R t d d r 'B r 'B ++ R vB r 'B ++ ==== rRr + = R vBaIaBR 2 vB r () + r +++ = vBaB0 == vIR r + = aIR r () r ++ =

PAGE 66

50rather a simple introduction. Further information can be found in Wertz [37] or Crane and Duffy [42], or Appendix A. Quaternions Quaternions were invented by William Rowan Hamilton in 1843. Prior to his discovery, it was believed impossible that any algebra could violate the laws of commutativity for multiplication. His work introduced the idea of hyper-complex numbers. Here real numbers can be thought of as hyper-complex numbers with a rank of 1, ordinary complex numbers with a rank of 2, and quaternions with a rank of 4. Hamilton’s crucial rule that made this possible: [4-6] Hamilton supposedly developed this rule while on his way to a party. When he realized what the solution was, he took out his pocket knife and carved the answer into a wooden bridge. This rule would forever change mathematics as was known at the time. Now mathematicians could look at algebra where communitivity did not work. This is where Gibbs and others developed algebra of vector spaces, and quickly eclipsed Hamilton’s work until recently. Quaternions, also known as Euler symmetric parameters, are more mathematically efficient ways to compute rotations of rigid and non-rigid body systems than traditional methods involving standard rotational matrices or Euler angles. Quaternions have the advantage of few trigonometric functions needed to compute attitude. Also, there exists a product rule for successive rotations that greatly simplifies the math, thus reducing processor computation time. Quaternions also hold the advantage of being able to interpolate between two quaternions (through a technique called spherical linear i2j2k2ijk 1 – ====

PAGE 67

51interpolation or slerp) without the danger of singularities, maintaining a constant velocity, and minimum distance travelled between points. The major disadvantage of quaternions is the lack of intuitive physical meaning. Most people would understand where a point was in cartesian space if they were given [1 2 3]. However, few would comprehend where a point was if given the quaternion [1 2 3 4]. Rotations of Rigid Bodies in Space. Quaternions are able to represent a rotation of a rigid body in space. To perform a rotation () of a rigid body about an arbitrary moving/fixed axis ( e ) in space, the quaternion representation of this operation is [4-7] [4-8] [4-9] Notice that only one sine and one cosine function call is needed to calculate a quaternion, where an euler rotation matrix would require three sine and three cosine function calls, one each for roll, pitch, and yaw. Since trigonometric function calls are computationally expensive, this is a great savings. Spacecraft Equations of Motion In this section, the kinematic and dynamic equations of motion for a spacecraft are presented. q q ˆ q4 T= Norm q () 1 = q ˆ e 2 -sin = q4 2 -cos = e e1e2e3 T= Norm e () 1 =

PAGE 68

52Attitude Kinematics A spacecraft’s orientation in space is represented by a quaternion. The quaternion kinematic equations of motion are given by Wertz [38] as, [4-10] where and [4-11] The terms with theby them signify a skew matrix. Spacecraft Dynamics First the equations of motion will be developed for a satellite and then expanded to include other dynamics. Thus a lengthy derivation is needed to derive something as simple as Euler’s equation. This, however, is necessary so that when reaction wheels and fuel slosh are included, understanding how these dynamics are incorporated will be easier. Simple Spacecraft To model a satellite, we starting off with a couple of definitions for the total mass of the spacecraft, first moment of inertia, and the second moment of inertia respectively. [4-12] where is a point mass in the spacecraft, N is the total number of points. These definitions are in reference to Figure 4-3 (Adapted from Hughes [31] p. 43). The second moment of inertia ( J ) is typically just referred to as the moment of inertia. The moment of inertia has two important properties: symmetry and positive definite. q 1 2 -q 1 2 -q () == – T– 0 = q () q4I33 q ˆ + q ˆT– = mmn n 1 = N= c mnrn n 1 = N= J mnrn Trn1rnrn T– ()n 1 = N= mn

PAGE 69

53The momentum and velocity of each point in the spacecraft is [4-13] where is the external forces acting on point p and is the force exerted by point m on point n The linear momentum ( p ) of the total spacecraft is [4-14] This equation could be further simplified by realizing that a spacecraft1 is a rigid body, and thus Thus [4-15] 1.Actually since spacecraft are composed of light weight, flexible materials, some authors model them as flexible structures.mnv nfnfmn m 1 = N+ = vnvor n+ = R vo= f n f mn p pn n 1 = Nmnvn n 1 = Nmnvor n+ ()n 1 = Nm voc + ==== c 0 = p ffn n 1 = N== Figure 4-3. A rigid body satellite.

PAGE 70

54where [4-16] which comes from Newton’s third law, The angular momentum of the spacecraft is: [4-17] Now differentiating this equation with respect to time, results in the following equation. [4-18] Realizing that results in the following equation. [4-19] where is the total external torque about the point O. Now if the equations are defined about the center of mass of a rigid body spacecraft, then , and the cross product with momentum in [4-19] disappears. The equations collapse to the more familiar ones for linear and angular momentum. [4-20] Now using the previous equations involving rotating reference frames, we can write the linear and angular momentum equations for a spacecraft in a rotating body fixed frame. [4-21] [4-22] fmn m 1 = Nn 1 = N0 = f mnfnm+ 0 = hornpnn 1 = N= h or npnrnp n + ()n 1 = Nvnvo– () mnvnrnfn + []n 1 = N== vnmnvn 0 = h ovopgo+ – = goc 0 = vovn= p m vo= h ogo= p m voc c ++ = p – p f + =

PAGE 71

55[4-23] As before, these equations can be simplified using certain assumptions. However, deriving the equations without the assumptions makes including internal dynamics (i.e., fuel slosh and reaction/momentum wheels) easier. These equations put into a matrix form are as follows: [4-24] [4-25] where M is the mass matrix of the system. The kinetic energy T of the system is given by: [4-26] [4-27] The dynamics can be obtained from the kinetic energy by using the following: [4-28] [4-29] These equations are called the quasi-Lagrangian equations and are shown in Hughes [31] and Meirovitch [39]. Meirovitch derives the equations using quasi-coordinates, which is somewhat difficult to understand and very long. However, this method has been used by some authors such as Miller et al. [40] to model various flexible appendages on a satellite and moving submasses inside of a satellite. h o – ho vopgo+ – = pho T= V v T= F fgo T= MV = M m 1c – c J = T 1 2 -VTM V = T 1 2 -m vo Tvo Tcvo1 2 -TJ + + = v T p = T ho= t d d vo T vo T + f = t d d T T vovo T + + go=

PAGE 72

56Angular Velocity The spacecraft’s orbit in this work is assumed to be circular. Since we are using an L VLH reference frame for our spacecraft, the angular velocity of the system must account for the orbital rate in addition to the maneuvers that the satellite is conducting. [4-30] [4-31] where B refers to the body frame, O refers to the orbital frame, I refers to the fixed inertial frame, refers to angular velocity of B relative to O, C is a rotation matrix, is the gravitational parameter of the Earth, and is the distance the spacecraft is from the center of the Earth. Internal Disturbances V arious factors complicate the dynamics of a satellite’s dynamics. Moving masses such as reaction wheels and sloshing fuel create additional linear and angular momentum that needs to be accounted for. This section will show how the above equations for a simple spacecraft need to be modified in order to account for these effects. Spacecraft with Reaction or Momentum Wheels Starting off with reaction wheels, we will introduce a method to account for any number of moving wheels in any configuration. However, since this work only deals with reaction wheels, we will not concern ourselves with momentum wheels which change their orientation and not their rate of spin. BI / BO / OI /+ = OI /Cn – = n 0 Rc 3----0T= BO / Rc

PAGE 73

57Now we will add the contribution of the reaction or momentum wheels to the spacecraft dynamics. Looking at Figure 4-4 (adapted from Hughes [31] p. 66), we can define the following terms: [4-32] where terms with the subscript sat refer to the satellite and terms with the subscript w refer to the wheels. The linear momentum is given by the following equations: [4-33] [4-34] [4-35] where b is the vector from the point O to the wheel. The angular equations of motion are given by the following equations: [4-36] mmsatmw+ = ccsatmwb + = JJsatIwa mwbTb1bbT– () ++ = p satmsatv csat + = p wmwv mw b + = p psatpw+ m v c + == hsatJsat csatv + = Figure 4-4. A rigid body spacecraft with a single reaction wheel.

PAGE 74

58[4-37] [4-38] [4-39] where a is a vector that describes the direction of the wheel axis when dealing with a single wheel. When there are multiple wheels, a becomes a matrix where the columns are the individual direction vectors for each wheel. For example, a system with n wheels would look like: [4-40] where [4-41] The wheel speed is a vector where each element is the speed of a specific wheel. Thus for a system with n wheel, the vector would be: [4-42] The equations of motion can be derived from kinetic energy for the system by using the quasi-Lagrangian. The kinetic energy for the spacecraft with wheel dynamics is: [4-43] Putting these equations into a matrix form results in the following equations. [4-44] hwIs w= hhsatha+ = hcvJ ahw++ = a aw 1{} aw 2{}… awn{} = awie1e2e3 T= Norm awi() 1 = w ww 1w 2… wn T= T 1 2 -m vTv 1 2 -TJ 1 2 -Iw w T wvTc Iw w TaT + – ++ = phoha T= V vo w T= F fgga T=

PAGE 75

59 [4-45] The equations of motion, like before, are as follows: [4-46] [4-47] [4-48] Although these equations look the same as the previous equations without wheel dynamics, remember the equations for linear and angular momentum are different. The linear and angular momentum of the spacecraft and wheels are now coupled in these equations. Equations for spacecraft with wheel dynamics are often simpler than what are shown here. However, if the wheels are assumed to be located at the center of mass (which is physically impossible) then the equations are greatly simplified, [4-49] since and Further more, since and substituting this into [4-49] results in: [4-50] [4-51] [4-52] This final equation is the common equation found in Wie [22] and Sidi [9]. MV = M m 1c – 0 c J Iwa 0TIwaTIwaT= p – p f + = ho – ho vopgo+ – = ha ga= ho – ho go+ = p m vo= vop vom vo 0 == hJ hw+ = c 0 = J hw + J hw+ () go+ – = J J hw + hw – go+ – = J J hw + ga– go+ – =

PAGE 76

60Fuel Slosh Fuel slosh, as stated in a previous chapter, is the unwanted movement of fuel inside a spacecraft. Each sloshing mode is modelled as a mass, with the forces exerted on the spacecraft body represented by a spring, and the effects of baffling in the fuel tank Figure 4-5. A plot of common environmental disturbances a spacecraft is subjected too.

PAGE 77

61represented by a dampener. The dynamics of this are incorporated into our simple satellite model by the following equations which will represent one sloshing mode. [4-53] [4-54] [4-55] where n is the position of the sloshing mode, is the mass of the fuel, and b is the vector from the center of mass of the spacecraft to the center of the fuel tank. These equations reflect how the linear, angular, and fuel slosh momentums change with the addition of fuel slosh’s mass-spring-damper. The mass matrix now becomes: [4-56] As with the wheel dynamics, the differential equations of motion for the linear and angular momentum remain the same. The new differential equation for the fuel slosh is: [4-57] where the and are the dampening coefficient and spring constant for the sloshing mode. Environmental Disturbances in Space There are many environmental disturbances that plague a spacecraft in orbit. These disturbances constantly effect the performance of the spacecraft which necessitates the use of a control system that counter acts them. The disturbance torques produced by some of p m vc mfn + – = hc J mfbn ++ = pfmfvb n + – () = mf M m 1c – mfc J mfb mfmfb – mf= p f m –f vrd – () cfn – kfn – = cf kf

PAGE 78

62these are shown in Figure 4-5 (adapted from Hughes [31] p. 271). A further explanation of the important disturbances follows. Gravity Gradient Euler and d’Alemert in 1749 first pointed out how the Earth’s gravitational field was not uniform. Later in 1780, Lagrange used these ideas to explain why the moon always had the same side facing the Earth. Since the gravitational field over a rigid body is not uniform, the center of mass is not the center of gravity. Thus there are torques about the center of mass which depend on the orientation of the spacecraft. The gravity gradient torque disturbance is found in [22], [31], and [38] as [4-58] where is the orbital rate of the spacecraft, is the third column in the orbital frame to body frame rotation matrix. When a satellite is in low Earth orbit, gravity gradient torques become the predominate environmental disturbance. Solar Snap In a previous chapter, modelling solar snap was discussed. The disturbance torque due to solar snap is a combination of the quasi-static displacement of the solar panel and its vibration. [4-59] Conclusion This chapter showed the basic concepts of spacecraft and disturbances dynamics. These dynamics produce a system where the satellite body, wheels, and fuel are all coupled. The disturbances presented in this chapter were internal (fuel slosh) and external ggg3 e 2o3ˆ Io3ˆ = eo3ˆ gssgqsgv+ =

PAGE 79

63(gravity gradient). A complete system can now be produced which has all of the dynamics of interest: satellite, wheel, fuel slosh, and thermally induced vibrations.

PAGE 80

64CHAPTER 5 ATTITUDE CONTROL OF SPACECRAFT This chapter will introduce the reader to attitude control of a satellite and specifically to a variable structure method called sliding mode control. First however we will cover some of the typical hardware used in a satellite to control its attitude. Next the traditional proportional derivative controller and sliding mode controllers will be introduced. Later in this chapter we will compare and contrast the two controllers. Satellite Control Hardware This section will provide a brief overview of hardware used by a satellite to control attitude. Control Computational Hardware Satellite’s are typically equipped with only the bare minimum in computational ability in an effort to reduce the cost of the vehicle. Unfortunately the complexity and performance of the control system is constrained only to the simpler algorithms. In fact some spacecraft, such as the lunar explorer, had no on board control systems. Instead it was remotely controlled from Earth, which was a very brittle solution at best. This was hailed as a great step forward in reducing the cost of space flight, but what happens if you lose the connect between ground control and the satellite? Any control algorithm used must be computationally inexpensive so it does not overburden the processor, and is typically designed to operate at 1 Hz. One of the simplest, easiest, and most used controller available is the venerable proportional

PAGE 81

65derivative (PD) controller. This is one reason that the PD has flourished in spacecraft control. Control Actuators In order to overcome environmental disturbances, satellites are typically equipped with various types of actuators that can provide an attitude control effort. Some of these are active control systems which utilized electrical power or consumable resources to effect the spacecraft’s attitude. While others are passive, and rely on conservation of energy to effect the satellite’s motion. Passive methods of control are generally relegated to dampening out disturbances during flight. These methods do not produce the same levels of torque as their active counterparts, but are typically employed along with active controls in an effort to reduce the power consumption of the spacecraft. Some of the most common types of active control actuators will now be discussed. A quick comparison of their control torques is given below. Reaction and momentum wheels Reaction and momentum wheels are cylindrical masses that rotate about their center of mass driven by a motor. They are produced in a variety of types and sizes, some of which are shown in Figure 5-1. These devices are capable of changing or stabilizing the orientation of a spacecraft with respect to its axis of rotation through conservation of angular momentum. These are the best for delivering smooth continuos control efforts. T able 5-1: Comparison of common active attitude control hardware. Control Hardware Control Force Range (Nm) Magnetic Torque Rods10E-6 (geostationary) 2.5E-3 (400 km) Reaction Wheels 0.05-2 Thrusters0.1-30

PAGE 82

66Wheels that have a nominal spin rate of zero1 are called reaction wheels. Wheels that have a constant momentum are called momentum wheels. Moment wheels are used to create a bias in the system which resists small parasitic environmental torques. Sometimes these wheels are also mounted in a gimbaled system that has one or two degrees of freedom. This allows the momentum wheels (with their constant momentum) to change the direction of the their momentum vector. Magnetic torque rods Magnetic torque rods are another method of attitude control and are shown below in Figure 5-2. They are based on the idea of generating torque by magnetic fields. Here the rods produce a magnetic field which interacts with the magnetic field of the Earth. Obviously this method of attitude control produces small control forces compared to 1.This statement is not always true. Sometimes wheels are biased by a small amount to fight friction in the wheels and bearings. These zero-crossing problems, dead zone around zero rpm, can effect the accuracy of a satellite when it is gathering scientific data.Figure 5-1. Two types of reaction wheels produced by Ithaco.

PAGE 83

67reaction wheels. Also the maximum amount of control effort produced by this method is dependant on the altitude of the spacecraft since the Earth’s magnetic field is stronger the closer the spacecraft is to the center of the Earth. In addition to small torque capabilities, electromagnetic shielding to protect satellite components must be considered. This has the additional disadvantage of adding weight and increasing the complexity of the design process. Gas jets Gas jets or control jets (shown in Figure 5-3) are the most powerful hardware available for attitude control, but these are not the most desirable to use for two reasons. First they consume fuel which is a limited, non-renewable resource, and second they do not have variable control but rather “on-off”. Thus bang-bang control schemes are created to provide pulse-width modulation control efforts to mimic a variable type of device. Proportional Derivative Control First let’s start off with one of the simplest to design and computationally efficient to implement, the proportional derivative (PD) controller. Wie and Barba [41] developed several computationally efficient control schemes for large angle maneuvers. Many of Figure 5-2. Torque rods produced by Ithaco.

PAGE 84

68these stabilizing schemes utilize quaternion and angular velocity feedback. The use of quaternion representation allows for more realistic, large angle maneuver control schemes. These schemes are formulated based on Lyapunov analysis, which produces a range of positive stable gains for that control law. Thus in order to meet desired performance, engineers must iterate through a significant number of gain combinations to obtain the desired response in a clean simulation. However, even when a satisfactory response is finally obtained, there are no guarantees how the satellite will behave in the presence of disturbances, noise, or uncertainties. Bong Wie [22] defines three different PD controllers utilizing quaternions in his book where the origin of the quaternion system is (0,0,0,1). [5-1] [5-2] [5-3] [5-4] where is the control effort and all k ’s and c ’s are positive. uKq ˆ – C – = K k 1 = C diagc1c2c3,, () = K k q4 3----1 = C diagc1c2c3,, () = K kq4() 1 sgn = C diagc1c2c3,, () = u Figure 5-3. Spacecraft engine produced by Boeing.

PAGE 85

69Sliding Mode Control Theory and Background This section will present an overview of sliding mode to the reader. Further in depth explanation can be obtained from Slotine [42] or DeCarlo et al. [43]. What is Sliding Mode? Nonlinear model based control systems offer a level of dynamic capabilities which linear techniques are incapable of providing when dealing with parameter uncertainties and unmodeled dynamics. Sliding mode, which has been studied in the Soviet Union for many years, is categorized as a variable structure control. This control structure has excellent stability, robustness, and disturbance rejection characteristics. Sliding mode is not a new control technique either, many researchers have utilized its robust properties to control a variety of systems such as missiles [44], mechanical systems [45], robotic manipulators [46][47][48],and submarines [49][50][51][52]. Theory First let’s look at a simple example of controlling a system. Given the system below, [5-5] [5-6] a solution to this differential equation can be found and is plotted in Figure 5-4. Notice that neither of the two solutions drive the system to zero at the origin, but rather remains oscillating with a stable behavior. The trick is to switch between the two solutions depending on which quadrant the system is in. [5-7] y ut () = ut () kyt () – = ut () k1yt () if yy 0 < – k –2yt () otherwise =

PAGE 86

70By switching, we can drive the system to the origin in a stable manor, as shown in Figure 5-5. Slotine took this idea and developed a systematic method for designing sliding mode controllers. Looking at Figure 5-5, we can develop the following table: This table can be summarized by the switching law. [5-8] where [5-9] T able 5-2: Summary of the switching law. Quadrant 1++ 2 + 3-+ 4 y yut () k2– k1k1k2– s () sgn 1 syy () 0 > – 1 syy () 0 < = syy () y y + =(a) (b)Figure 5-4. Two possible solutions for the double integrator, where (a) is and (b) is These plots are typically referred to as phase portraits. They plot position and velocity in this example, but later (when we are doing controls with full state feed back) the portraits will represent error and error velocity. y y y y kk1= kk2=

PAGE 87

71is called your sliding surface and the control effort is given by: [5-10] where K is typically a scalar or diagonal matrix of positive definite gains. In order to design a controller, we will now replace our states () with our errors ( ). This will now (as previously stated) drive our system errors to zero in a stable manor. Second, we will modify the control effort equation to include a feed forward aspect as follows. [5-11] The term is called the equivalent control effort. It is model based, thus if we have a perfect model of the system, we can calculate the required control effort to produce the desired results. However, that is never the case. There are always disturbances and unmodelled dynamics that influence our system. Thus the job of the second term in [5-11] is to reduce the error and error velocity produced by the equivalent control in a stable manor. Also note that the is a weighting factor between the error and error velocity. This term can be adjusted depending on which is more important, and its effects are depicted ut () K – s () sgn = yy x x ut () u ˆ K – s () sgn = u ˆ Figure 5-5. Switching between the two controllers results in driving the system to zero in a stable manor. y y

PAGE 88

72below in Figure 5-6. Typically this term is constant, but there is no reason it could not be dynamic. Simple Example For a simple example, we will develop a sliding mode controller for a mass-springdampener system. The dynamic equations for this system is [5-12] Our sliding surface and equivalent control are [5-13] [5-14] [5-15] [5-16] Now looking at [5-16] we notice that if there is no error in the modelling of the system or disturbances, then this control effort is sufficient to control the system. Also if we had a model of any disturbances in the system, we could include them in our equivalent control mx cx kx ++ u = sx x + = s x x + 0 == m1 –u ˆ cx – kx – () x d– x + 0 = u ˆ mx d x – cx kx ++ = Figure 5-6. Graphical representation of how effects the sliding surface’s orientation in the phase plane. The dotted lines represent the sliding surface (s = 0). y y s = 0, 1 = y y s = 0, 1 45 =

PAGE 89

73term. However in this example, there are no disturbances included in the system or we do not have a model for any disturbances seen by the system during run-time. Unfortunately there will always be error in modelling, noise, or disturbances. Thus we must include another term to account for these problems. [5-17] Other Sliding Mode Designs for Spacecraft Many authors have suggested many types of control architectures for satellite attitude control [53-59]. However, the one of the most interesting and simpler to implement types of control is sliding mode. The use of sliding mode control is not new to satellite attitude control. Lo and Chen [60] designed a sliding mode controller scheme which avoids the inverse of the inertia matrix, and smooths the control effort of the controller. Such a strategy provides for a more efficient use of fuel. Their simulations involved small uncertainties in the spacecraft inertia matrix and a sinusoidal disturbance. Simulations showed favorable results. The control effort was [5-18] where the for all intensive purposes. Boskovic et al. [61] and [62] developed a sliding mode controller, but specifically designed the controller such that it did not saturate the control effort. They were able to produce accurate results, and even tested the controller on a system with a larger inertia matrix than what the controller was designed with. Unfortunately they drew an incorrect conclusion that the controller is independent of the inertia matrix it is presented with. Sliding mode will be able to control a system with a larger inertia matrix, but may not control a system with a smaller inertia matrix. This is because the control effort is model uu ˆ Ks () sgn – = u J Jv d Jq – ks () sgn – + = v d d

PAGE 90

74based, and thus it is possible to make a satellite go unstable by issuing too much control effort. The simulations utilized a square wave disturbance and small uncertainties in the inertia matrix. The controller developed was as follows: [5-19] Crassidis et al. [63] developed an optimal variable-structure (sliding mode) controller for large-angle maneuvers on satellites. A cost function was developed, which when minimized, resulted in the optimal control effort. The the minimization lead to a two point boundary value problem. The cost function and the resulting optimal control effort were [5-20] [5-21] Simulations of the controller on the MAP satellite showed favorable results. Coleman and Godbole [64] conducted a performance trade study between fuzzy logic, PID, and sliding mode control. The controllers were tuned for and tested on three similar linear plants. In all three cases, the sliding mode controller out-performed the fuzzy logic controller. This is a typical result, where sliding mode tends to provide a superior model based performance compared to fuzzy logic, assuming the model is accurate enough. However all of these authors utilized simple satellite dynamics with a simple external disturbance in their simulations. Fuel slosh, which is coupled with the dynamics of the satellite, is a more complex problem. It is difficult to accurately model the real disturbance and requires a control system that can account for this uncertainty. u J Jq – diag s1() sgn s2() sgn s3() sgn – k1k2k3= () 1 2 -q Tq T + () t dts = u J J 1 2 -q 4() q () qd() d qd() q () – [] sgn + dks () sgn – + =

PAGE 91

75Satellite Controller Design For This Work This section will cover the sliding mode controller developed for this research. Controller Derivation This formulation, which is a full-state feedback technique, utilizes the existing dynamic model and compensates for uncertainty while formulating a control effort that tracks the desired trajectories. The equations of motion for a satellite’s attitude (in quaternion space) are given by [22] [5-22] [5-23] where [5-24] Also f in [5-22] represents the modelling error and potential disturbances in the system. This term will provide a feed forward aspect to the controller and should yield better results. The q is a vector composed of the three imaginary elements of a quaternion. However, in this formulation, we can take the state variable to be or small euler angle orientation, and equate this to q [5-25] Thus the sliding surface ( s ) is given by: J – J fu ++ = q 1 2 -q 1 2 -q4 + = 0 3– 230 1– 2– 10 = q q1q2q3= xyz2 qxqyqz =

PAGE 92

76[5-26] where the error terms are defined as: [5-27] In order to calculate the equivalent control effort (), we need to take the derivative of the sliding surface and set it equal to zero. [5-28] Now multiplying both sides by the inertia matrix (J) will allow us to substitute in the equations of motion ([5-22] and [5-23]). [5-29] [5-30] The hat over particular variables denotes that they are estimates, since accurate values for the inertia matrix and disturbances may not be known. Since the sliding surface does not move, its time derivative is zero, and thus the equation becomes [5-31] Now solving for the equivalent control effort () from [5-31] yields: [5-32] The quaternion error rate in the above equation is defined as: [5-33] s t d d + n 1 –x x x + q + === x xxdesired– = u ˆ s q + = Js J ˆ J ˆ q + J ˆ J ˆ q J ˆ desired– + == Js J ˆ f ˆ u ˆ J ˆ q J ˆ desired– +++ = Js J ˆ f ˆ u ˆ J ˆ q J ˆ desired– +++ 0 == u ˆ u ˆ – J ˆ f ˆ – J ˆ q J ˆ desire d + – = q 1 2 - q 1 2 -q 4 + =

PAGE 93

77[5-34] The estimated disturbance in the above equation is equal to the estimated solar snap and fuel slosh disturbance. [5-35] where is the disturbance due to solar snap and is the disturbance due to fuel slosh. Finally, the sliding mode control effort is given by: [5-36] The K term in [5-36] will be discussed when stability of the controller is covered. This switching term provides the appropriate control action to quickly drive the trajectories back onto the sliding surface. Since real actuators do not have an infinite bandwidth, the sign function tends to cause excessive chatter in the control effort. The sign function is typically replaced by a saturation function. The saturation function not only smooths the control response, it also reduces the amount of fuel used by the controller. Fuel Slosh Disturbance The effects of the fuel slosh are incorporated into the controller by the term in the controller. A discrete state space system is used to track the effects of the fuel slosh: [5-37] [5-38] where the state vector is: 0 3 – 2 3 0 1 – 2 – 1 0 = f ˆ f ˆssf ˆfs+ = f ss f fsuu ˆ K sgn s () – = f ˆfsx fsAxfsBu + = f ˆfsCxfs=

PAGE 94

78[5-39] where x is the position of the fuel slosh and is its velocity. The input to the system u is the previous control effort used by the satellite. Also the state space matrices are discrete time representations to reduce the computational load of the controller. Stability of Controller The following Lyapunov candidate function is proposed. This function is clearly positive definite. [5-40] T aking its time derivative and substituting in previous equations of motion and control effort gives: [5-41] where the true dynamics of the system are inserted into the equations are [5-42] The modelled dynamics, which define the controller output are given by: and [5-43] After some algebra, the dynamics cancel each other out, and the only terms left are the difference between the true disturbance forces and the estimated/modelled ones and the last terms on the control effort equation. Thus [5-41] becomes: [5-44] where xfsx x = x V 1 2 -sTJs = V sTJs = Js – J uf Jq +++ = uu ˆ K s () sgn – = u ˆ J ˆ J ˆ q – = V sTF K s () sgn – () =

PAGE 95

79[5-45] Note that any differences between the true inertia matrix and the real inertia matrix, plus any other modelling error or disturbance, is contained in f Where f is the true disturbances, modelling errors, unmodelled higher order terms, etc. For a stable system, the lyapunov derivative needs to be negative definite. All variables are known except for the scalar gain value K Two of the three remaining terms are positive definite, K and F Thus we need to do further simplifications and algebra to solve for the gain. [5-46] [5-47] where [5-48] This operation is a cross between an absolute value operation (if it were a scalar) and a norm of a vector (if the sign function was not present). Thus this operation will be denoted by the angled brackets and not confused with either absolute value or norm of a vector. A quick example is shown below, and note that the result will always be scalar and positive. [5-49] Thus the scalar gain K is [5-50] Fff ˆ – = sTFsTK s () sgn – 0 < sTF K s – 0

PAGE 96

80In order to guarantee the above inequality, we will add a little offset () that will always be positive. [5-51] Now substituting this gain back into [5-47] results in: [5-52] This equation has two results depending on whether s is positive or negative. If s is negative then the two ‘s sum together and the entire equation is negative definite. If s is positive, then the two cancel each other out and the equation is still negative definite. or [5-53] Saturation The sign function used in the sliding mode controller leads to an excessive amount of chattering in the control effort. A common modification to the sliding mode controller is to change the sign function to a saturation function with the following properties. [5-54] Thus our sliding mode control equation changes to [5-55] where the is a scalar value that can be used to adjust when the saturation function will saturate. All other sliding mode equations remain the same. K sTF s ----------+ = sTF sTF s ----------+ s – 0 < sTF 2 sTF s – 0 < s – 0 < sata () 1 a < – 1 < a elsea () sgn = uu ˆ Ksat s -– =

PAGE 97

81Comparing The Two Types of Controller Slotine’s sliding mode formulation and a PD controller have some very similar attributes when the sliding mode is in the sliding region. Simplifying the system some, assume that the equivalent control effort () is zero, since this is really a feed forward term. Thus the control effort is now: [5-56] The sliding surface s is very similar to a PD control effort. [5-57] where is a rate error and is a position error. When the errors are on the sliding surface (meaning when the errors are small, they are in the saturation bounds), the saturation function returns a scaled version of the values passed to it (as opposed to a 1 or -1 when outside of the saturation bounds). Thus the control effort becomes scaled errors and scaled rate error multiplied by a gain K [5-58] [5-59] Then the PD parameters can be written as follows: [5-60] Conclusion This chapter introduced PD and sliding mode theory for satellite attitude control. An in depth discussion of sliding mode design, stability, and a comparison between sliding mode u ˆ u Ksat s -– = s t d d + n 1 –x x x + == x x u K s -– = u K --x x + – = kp K ------= Ti1 -=

PAGE 98

82and PD was provided. This laid the foundation for the similarities and differences between the two controller which will soon manifest themselves when the results are presented.

PAGE 99

83CHAPTER 6 SIMULATION AND RESULTS This chapter will discuss the simulation, the gains chosen for the PD controller, and finally discuss the results. The first section will describe the specifics of the satellite and the simulation. All of the variable from proceeding chapters will be filled in and explained. Next the simulation, equations of motion, and the controllers that were developed are covered. Finally the results of the simulation are discussed. Simulation Dynamics and Parameters The simulation used to evaluate the performance of the controllers was written in C++. The first simulations were written in Matlab because of the ease of working with vectors and matrices. The problem with Matlab, however, was the amount of time each simulation took. Once the simulation were rewritten for C++, they were magnitudes faster. The reader is referred to Appendix C, where some of the mathematical code developed is presented. Satellite The satellite shown in Figure 6-1 and Figure 6-2 is Earth Observing 1 (EO-1), and will be the satellite used in all simulations. EO-1 was primarily constructed by Swales Figure 6-1. EO-1 satellite with solar panel deployed.

PAGE 100

84Aerospace and is part of the New Millennium Program which is dedicated to validating new technologies. The spacecraft was launched aboard a Boeing Delta II rocket from V andenberg Air Force Base on November 21, 2000. Some of the new technology being tested is the Advanced Land Imager, hyperspectral imaging spectrometer called Hyperion, X band phased array, pulsed plasma thrusters, and formation flying with Landsat 7 spacecraft which is already in orbit. T able 6-1: Technical Specifications for the EO-1 spacecraftDry Mass568 kg V olume of Bus 1.5m x 1.5m x 2m Inertia Tensor[443 179 429] kg m^2 ACS Zero Momentum, 3 axis stabilized Command and Data Handling Bus Architecture Mongoose V, Rad Hard at 12 Mhz, RISC Architecture Solar Arrays 3 panel / Si w/GaAs / Articulating Bus StructureHexagonal with aluminum honeycomb Propulsion 1 fuel tank with 4 thrusters Propellent Capacity23 kg Mission Design Life 1.5 years Figure 6-2. Diagram of satellite’s internal components.

PAGE 101

85Satellite Dynamics The true equations of motion for the satellite’s linear and angular momentum are the now familiar ones: [6-1] [6-2] The differential equation for the linear momentum of sloshing mode i is: [6-3] where for this work there are two sloshing modes. The wheel dynamics are: [6-4] The mass matrix of the spacecraft, which contains the linear and angular momentums for the various components, is: [6-5] [6-6] [6-7] p p – = h hvpg + – – = p i mfi vrfi – () cfni kfi– – ni= ha ga= MV = phop1p2ha T= V v n1 n2 s T= M m 1c – 1 mf 11 mf 20 c J mf 1b mf 2b Isa 1 mf 11 mf 1– b 1 mf 1– 00 1 mf 21 mf 2– b 0 1 mf 2– 0 0 Isa 00 Is1 =

PAGE 102

86Fuel Tank The location of the fuel tank in the body frame is and the remaining parameters are shown below. Reaction Wheels The reaction wheels used on EO-1 are customized versions Ithaco TW-4A12. Their specifications used in this work are shown below. Each of the three reaction wheels is located on one of the three prime axes, half a meter from the center of mass. T able 6-2: EO-1 fuel slosh parameters. Mode 1 Mode 2 Fuel Mass () 15 kg8 kg Natural Freq.() 0.1 Hz 0.3 Hz Dampening Ratio () 0.50.5 Neutral Location from C.M. (b) 1 m 1 mT able 6-3: EO-1 wheel specifications.Mass2.55 kg Speed 5100 rpm DimensionsDia. 20.5 cm Height 6.4 cm Limits T orque 0.025 Nm Momentum 4 Nms Location from C.M..5 m b 001T= mf Figure 6-3. Composition of a typical simple solar panel.

PAGE 103

87Solar Panels The parameters used for the solar panel in this work are shown below in Table 6-4. Since a solar panel is actually a composite of various materials and components as shown in Figure 6-3. Thus the parameters are primarily based on the largest component, the aluminum honey comb core. Johnson and Thorten [68] provides a formula for calculating the material constants taking into account that the solar panel is a composite structure. However the formula appeared in a journal article with no explanation where it came from, nor of the constants used to weight the different parameters. These parameters have been determined based on what little information that is publicly available from the manufacturer of EO-1’s solar panel. Note that the dampening ratio used was a very small number since structural dampening is assumed. Control System The control system issues control corrections to the spacecraft at 1 Hz with the colored noise filter running at a faster 4 Hz. The controllers operate with a 5% error in the inertia matrix. T able 6-4: EO-1 solar panel parameters.MaterialAluminum 6061 honey comb Dimensions (w x L x h) 1.428 m x 3.8 m x 0.03 m Density () 300 kg/m^3 Natural Freq. () 0.48, 3.0, 8.4 Hz Structural Dampening Ratio () 0.1 Structural Stiffness (EI) 2E4 N m^2 Specific Heat (c)920 J/kg K Thermal Conductivity (k) 1.2 W/m K Coefficient of Thermal Expansion () 22.7E-6 1/K Location from C.M. () 1 m TRo

PAGE 104

88PD: SM: Finding an Optimal PD Gain In order to compare the PD controller to the sliding mode controller, an attempt to find the best gains was conducted. A simulation was run which looped through various values of the PD controller. [6-8] where is the control effort, is the proportional gain, is the derivative gain, is the attitude error in the system, is the rotational rate error, and The dynamics used for the satellite where basically Euler’s equation, which are commonly used to represent a simple spacecraft. [6-9] where is the angular momentum of the satellite, is the angular momentum of the wheels, and are the external disturbance torques. These dynamics were used because they are simpler and allowed the search to proceed at a greater rate. Also no mathematical model ever exactly matches the real system. Thus it was seen as undesirable to train the satellite on the same model that the experimenting would be conducted on. The reasoning, again, is that is real life whatever model the controller was trained on before the spacecraft was put into space would have some error, big or small. The satellite was commanded to perform a step maneuver which rotated the satellite through a rotation from euler angles [-30 -30 -30] to [30 30 30]. Each time the values for and where changed and various performance criteria were calculated. The Kp1.5 = Td0.1 = 0.1 = 0.3 = 0.1 = u kpe kde + kpe Tde + () == ukpkde e Tdkdkp---= hsat hsat – hwNw– Next+ – = hsathwNextKpTd

PAGE 105

89performance measures used were the Integration of Absolute Error (IAE), Integration of T ime times Absolute Error (ITAE), and settling time () of the system. [6-10] [6-11] The error used in the calculations came from the quaternion error () and the absolute value bars indicate the euclidian norm of the quaternion error vector. Then a cost function was calculated using these individual measures to determine the optimal controller gains. [6-12] tsIAEqerrorqerrord0 t= ITAEqerrortt d0 t= qerror CostCostitaeCostiaeCostts = Figure 6-4. Cost function surface plot for the optimal gain of the PD controller.

PAGE 106

90Figure 6-5. ITAE results from the optimal PD search. Figure 6-6. IAE results from the optimal PD search.

PAGE 107

91This cost function was used because these were some of the performance measurement used to evaluate the different controllers. Thus optimizing the PD controller’s gains for the performance measurements would help it perform better. Each of the individual measures were also normalized by dividing them by the maximum value found after the search was completed. For example, the ITAE would have been calculated as follows: [6-13] [6-14] Thus the values for the cost functions that interest us most are the highest values (i.e., 1.0). The smaller error in the system and the shortest settling time will all have a value of one. Also the cost function is now unit less since there is a division by the maximum value (i.e., the units cancel). itae 1 ITAE ------------= Costitaeitae maxitae () ------------------------= Figure 6-7. Settling time results from optimal PD search.

PAGE 108

92The resulting surface plot of the cost function is shown in Figure 6-4. Here the optimal gains were: Kp = 1.5 and Td = 0.1. Looking at the individual components of the cost function, the IAE and ITAE in Figure 6-6 and Figure 6-5 produce similar surfaces. Although the ITAE produced a narrower ridge of gains, both criteria prefer lower gain values that run along a ridge approximately located at Td = 0.1. Next looking at the settling time results in Figure 6-7, a similar ridge to the error surfaces was formed. The ridge again showed a preference for lower gains like the error surfaces. Here, however, the ridge resides slightly lower than the error gains; running approximately along Td = 0.06. The cost function for finding an optimal gain was composed of the error terms and the settling time. Originally, however, a term for the amount of power / fuel / control momentum was intended on being included. The results of this performance criteria is Figure 6-8. Control momentum or fuel used searching for the optimal PD controller.

PAGE 109

93shown in Figure 6-8. This term was not used since it favors higher Td values, but when looking at pervious surfaces, it becomes apparent that performance would suffer greatly. Fuel Slosh This section and the one following it will talk about the performance of three different controllers. The controller used in this work were the proportional derivative (PD), sliding mode (SM), and sliding mode with a colored noise filter to compensate for fuel slosh (SMCNF). A step simulation was used to demonstrate the controller’s ability to react to the two main disturbances: fuel slosh and gravity gradient torque. The simulation had all euler angle set equal to each other, but the initial and final conditions where differentiated by sign. Table 6-5 shows the three sets of initial and final angles used along with the angular distance between the two orientations. The results of the simulations were compared by looking at four performance criteria. They were ITAE, IAE, Ts (settling time), and CM (Control Momentum in Nms). The CM was determined by integrating the euclidian norm of the control effort ( u(t) ) used by the controller from time 0 to the end of the simulation (). [6-15] T able 6-5: Initial and final euler angles where roll, pitch, and yaw are equal. Start (deg) Stop (deg) Angular Distance (deg) -101034.46 -30 30 98.99 -6060165.64 tfinalCMut () t d0 tfinal=

PAGE 110

94IAE The IAE is the sum of error in the system. The amount of error in a system is proportional to the amount of control effort used to get reduce the error. In a dynamic system, such as the satellite, there are other factors that contribute to the amount of error and control effort used, such as momentum and external disturbances. The results for the simulations were: The sim column in Table 6-6 refers to the euler angles. Thus the 10 row means the satellite rotated from euler angles [-10 -10 -10] to [10 10 10]. Notice that the numbers decrease from left to right. This means that the SM controller produced less error than the T able 6-6: IAE results with values having units of radians. Sim PD SM SMCNF 10230198197 30 887 861 858 60229221692142 Figure 6-9. A depiction of the step maneuver showing the orientation of the satellite at start and end. The angular distance and axis of rotation are also shown for better understanding of the terms.

PAGE 111

95PD controller, and the SMCNF controller produced less error than the SM and PD controllers. The only exception to this was the results for the PD controller’s IAE at 60. ITAE The ITAE is a measurement similar to the IAE, but ITAE keeps track of when an error occurs. Errors that occur early in a simulation are weighted very small while errors that occur later in the simulation are weighted larger. This performance criteria has the effect of favoring systems with short settling times, but does not penalize controller that have large overshoots very early in the simulation. The numbers in Table 6-7 again decrease from left to right, meaning that the SMCNF controller performed better than the previous two. Settling Time The settling time is the amount of time the system takes to reach a steady state value. The spacecraft reached steady state when its quaternion error was below 100 arcsec1. T able 6-7: ITAE results with values having units of radian-seconds2. Sim PD SM SMCNF 10424927562745 30 18036 16476 16377 60845047192468917 1.An arc second (arcsec) is 1/3600 of a degree or one degree is equal to 3600 arcsec. This is a very fine unit of measurement used for high precision orientation.T able 6-8: Settling time results with values having units of seconds. Sim PD SM SMCNF 10509278278 30 576 370 362 60918730721

PAGE 112

96The results shown in Table 6-8 show that the SM produces faster setting times than the PD, but the SMCNF does not always produce faster setting times than the SM. CM The control momentum (CM) is a measurement of the amount of energy used by each controller. The results for the three simulations are shown in Table 6-9. The SM and SMCNF controllers utilized more fuel than the PD controller, while the addition of the CNF to the SM produced a very small decrease in the amount of fuel consumed. To r que and Momentum Limits The final point of interest for this section was the maximum control efforts. The momentum wheels are limited to 0.025 Nm of torque by the simulation. The simulation T able 6-9: CM results with values having units of Nms. Sim PD SM SMCNF 104.835.325.35 30 10.62 10.64 10.86 6019.3718.0117.55 Figure 6-10. Quaternion error and control effort of the PD controller rotating from -30 to 30.

PAGE 113

97does nothing to limit the momentum that the wheels can produce, but the maximum wheel momentum produced by the controllers was recorded. T able 6-10: Maximum wheel momentum achieved where units are in Nms. Sim PD SM SMCNF 101.591.661.67 30 3.44 3.43 3.45 605.154.724.70 Figure 6-11. Quaternion error and control of the SM controller rotating from -30 to 30. Figure 6-12. Quaternion error and control of the SMCNF controller rotating from -30 to 30.

PAGE 114

98The wheels on EO-1 have a momentum limit of 4 Nms. Thus during the final maneuver of 60, the wheels would have saturated and the control system. The thruster are used to unload momentum and the wheels are slowed down to some nominal value. This maneuver would require a much more complex control system involving a thruster with bangbang control which is outside of the scope of this work. Qualitative Analysis of Results Plotting the quaternion error for the 30 degree maneuver showed a smooth decrease of the error for each controller, as shown in Figure 6-10, Figure 6-11, and Figure 6-12. Thus there was little overshoot and waste of energy in the system. Although there was a definite difference between the PD and sliding mode controllers, the two sliding mode controllers showed little difference in their error plots. This was expected since the addition of the colored noise filter only slightly increases the performance of the sliding mode controller. Solar Snap The final simulation done involved looking at the effects of solar snap on the controllers. Performance Here the satellite tried to maintain a nadir1 orientation as it orbits the Earth. The solar panel was then subjected to a heat flux which induced vibrations. 1.A nadir orientation means the satellite is pointing towards the center of the Earth.T able 6-11: Solar snap results of the three controllers. PD SM SMCNF IAE (rad)3.71.00.9 ITAE (rad-sec^2) 14.9 2.3 2.1 CM (Nms)0.520.570.50 ts (sec) 139 46 45

PAGE 115

99The results in Table 6-11show that the SMCNF produced less error and a quicker settling time than the other two controllers, but at the expense of using more control effort. Qualitative Analysis of Results Again looking at the error produced by the three different controllers reveals the sliding mode controllers have better performance than the PD controller, as shown in Figure 6-13, Figure 6-14, and Figure 6-15. The error curves for the two sliding mode controllers looks identical. Figure 6-13. Arcsecond error of the PD controller during solar snap.

PAGE 116

100Figure 6-14. Arcsecond error of the SM controller during solar snap. Figure 6-15. Arcsecond error of the SMCNF controller during solar snap.

PAGE 117

101CHAPTER 7 CONCLUSION This work revolved around comparing two different types of attitude controllers: proportional derivative (PD) and sliding mode (SM). In order to make the comparison as good as possible, the PD controller was optimized to produce the best effort possible in this simulation. The SM was tuned by trial and error to achieve an acceptable results. The tuning took only a short amount of time, and easily produced better results than the PD controller. Thus the SM parameters used in this work may not be globally optimal. The simulation written to test these controllers was written in C++. The benefits of programming in a language such as C++, as opposed to using prototyping software such as Matlab, was speed. The solar snap disturbances required very small time steps to avoid numerical instability, and under Matlab this translated into very long simulation times. Thus code written in C++ allowed equations to be written in clear mathematical expressions, but produced very fast execution times. A simulation in Matlab would take approximately 10 minutes, while the same simulation written in C++ would take 20 seconds. Since simulation times were now short, much experimentation could easily be conducted in a short amount of time. Fuel slosh disturbances are a difficult problem to handle. They are impossible to measure or estimate without making assumptions about the form of the sloshing modes (e.g. pendulum, mass-spring-damper, liquid slug, etc.). These assumptions greatly simplify the problem to something that can be realized in an embedded system. Otherwise

PAGE 118

102very complex models which require massive amounts of computing power would need to be employed which is impossible in an embedded system like a satellite. The SM controller was shown to provide better control than the venerable PD controller. The SM controller was designed to handle inaccuracies in the system model. Here the inaccuracies are movement of the fuel which change the inertia matrix during maneuvers. In order to improve the performance of the SM controller, a colored noise filter (CNF) was added to account for the disturbances. The results shown in the previous chapter illustrated how the CNF improved the SM controller’s performance over the PD. The SM and SMCNF controllers produced good results when encountering solar snap activity. Typically satellites stop their scientific mission during times of possible solar snap disturbance because the errors generally exceed the allowable range. Here the SM and SMCNF produce much quicker settling times and less than half the error of the PD controller. Thus the science mission could be run with less down time due to solar snap and could possibly even be run during solar snap disturbances if the SM or SMCNF errors are with in an acceptable range. This work clearly showed the superiority of SM or SMCNF over a PD controller. The question now becomes, which of the two types of sliding mode controllers was better. That really depends on the mission for which they would be applied. The SM is less computationally expensive than the SMCNF, and for satellites where fuel sloshing disturbances are small, computational power is very limited, and high precision is unnecessary, the SM controller would most likely be the better choice. This would be especially true where engineers are unfamiliar with sliding mode, and this is their first opportunity to work with it. However, in the situation where the opposite is true, then the

PAGE 119

103benefits of the CNF at the expense of increased computational cost and increased development time become attractive. Thus, in this situation the SMCNF would be a better controller to implement.

PAGE 120

104APPENDIX A ATTITUDE REPRESENTATIONS AND ROTATION MATRICES When dealing with a subject such as inertial navigation, controlling spacecraft, or computer vision, it is common to have to deal with multiple reference frames and have to relate the position of one point in one frame to its position in another reference frame. Thus it is important to know how to rotate points from one frame to another efficiently and understand what limitations may be present in various methods. This Appendix will attempt to familiarize the reader with quaternions, rotations, attitude in space, and euler angles. Comparisons will be drawn between euler angles and quaternions to show their strengths and weaknesses. However it is assumed that the reader is already somewhat familiar with euler angles. Fixed Angle Rotations T ypically when trying to relate a common point between different reference frames, a rotation matrix will be constructed. This matrix will allow easy transformation between frames. [A-1] Here a point (P) is rotated from where it is originally defined, in reference frame a, to another reference frame b. Not that the rotation matrix R does not include any translation, but only rotation. Note also that R has the following properties: PbRa bPa=

PAGE 121

105When a rigid body is rotated about a fixed inertial reference frame, it is referred to as a fixed angle rotation. These types of rotations are common in robotics where one is trying to determine the joint angles of a serial robot relative to some point in space. What is important to understand, is rotations of this type are made about the inertial axis which is inconvenient. A better method would be rotations about the body axes. These type of rotations are known as Euler Angles. Euler Angles Euler angles are typically though of in terms of roll, pitch, and yaw angle about body axes. These terms are shown graphically below for a rigid body in space. T able A-1: Properties of a rotation matrix. Rotation Matrix R R1 –RT= Ra bxayazaxbybzb== R 1 = columniR () rowiR () 1 == Ra cRb cRa b= Figure A-1. Body reference frame attached to a rigid body. The x-axis points out the front of the vehicle.

PAGE 122

106[A-2] There are two different groups of euler rotations out of the possible twelve. The two groups are distinguished by their singularity locations. Typically aerospace and navigational applications use the 1-2-3 rotation where this singularity only effects fighter aircraft who dive or climb at such steep angles. Some space applications, such as orbits around the Earth, use the 3-1-3 rotation sequence. In order to perform these operations on the rigid body, a rotation matrix (R) is used to find the new orientation of the spacecraft given the old orientation. Given below are the attitude matrices for rotations about the x, y, and z-axes. [A-3] [A-4] [A-5] T able A-2: Comparison of the two major types of rotations, sequential and first and third axes rotation Singularities in Pitch Possible Rotation Sequence T ype I+/-901-2-3, 1-3-2, 2-1-3, 2-3-1, 3-1-2, 3-2-1 T ype II 0/180 1-2-1, 1-3-1, 2-1-2, 2-3-2, 3-1-3, 3-2-3 roll pitch yaw = Rx () 100 0 () cos () sin 0 () sin – () cos = Ry () () cos0 () sin – 010 () sin0 () cos = Rz () () cos () sin0 () sin – () cos0 001 =

PAGE 123

107Now, subsequent rotations about these primary axes can be accomplished by multiplying the matrices together. Thus, successive rotations about the z-axis, x-axis, and y-axis are given by: [A-6] Quaternions Quaternions where invented by William Rowan Hamilton in 1843. Prior to his discovery, it was believed impossible that any algebra could violate the laws of commutativity for multiplication. His work introduced the idea of hyper-complex numbers. Here real numbers can be thought of as hyper-complex numbers with a rank of 1, ordinary complex numbers with a rank of 2, and quaternions with a rank of 4. Hamilton’s crucial rule that made this possible: [A-7] Hamilton supposedly developed this rule while on his way to a party. When he realized what the solution was, he took out his pocket knife and carved the answer into a wooden bridge. This rule would forever change mathematics as was known at the time. Now mathematicians could look at algebra where communitivity did not work. This is where Gibbs and others developed algebra of vector spaces, and quickly eclipsed Hamilton’s work until recently. Quaternions, also known as Euler symmetric parameters, are more mathematically efficient ways to compute rotations of rigid and non-rigid body systems than traditional methods involving standard rotational matrices or Euler angles. Quaternions have the advantage of few trigonometric functions needed to compute attitude. Also, there exists a product rule for successive rotations that greatly simplifies the math, thus reducing Rzyx ,, () Rz () Ry () Rx () = i2j2k2ijk 1 – ====

PAGE 124

108processor computation time. Quaternions also hold the advantage of being able to interpolate between two quaternions (through a technique called spherical linear interpolation or slerp) without the danger of singularities, maintaining a constant velocity, and minimum distance travelled between points. The major disadvantage of quaternions is the lack of intuitive physical meaning. Most people would understand where a point was if they were given [1 2 3]. However, few would comprehend where a point was if given the quaternion [1 2 3 4]. This section does not attempt to provide the extensive understanding needed to employ quaternions but rather a simple introduction. Further information can be found in W ertz [37] or Crane and Duffy [42]. Quaternion Algebra The quaternion is composed of a scalar and a vector part. The scalar is a redundant element that prevents singularities from occurring since the four elements are all dependent upon each other.1[A-8] where the v stands for vector and the r is a scalar real. Given below is a table that will summarize the important mathematical operations of quaternions. 1.The order of the quaternion elements is not standardized. I have chosen the element order that NASA uses which is the imaginary part first then the real. Some authors put the real first then the imaginary which is similar to complex numbers. Thus before you use any equations that utilize quaternions, make sure to understand how the author is arranging the elements.q qxqyqzqr Tq1q2q3q4 Tvr ixjykzr +++ ====

PAGE 125

109Multiplication is an important operation for quaternions, thus we will elaborate on it a little. [A-9] where [A-10] [A-11] [A-12] [A-13] or rewritten in matrix form. T able A-3: Quaternion Algebra Summary Operation Formula Add / Subtract Scalar Multiplication Norm Quaternion Multiplication aa.Notice that the order of the quaternions in the matrix multiplication have been reversed. Conjugate Inverse qq qxqx' qyqy' qzqz' qrqr' T= q qx qy qz qr T= qqx 2qy 2qz 2qr 2+++ = q q q4q3q –2q1q –3q4q1q2q2q –1q4q3q –1q –2q –3q4q = q q –xq –yq –zqr T= q1 –q q ------= p qrr4ir1jr2kr3+++ == r4p4q4p1q1– p2q2– p3q3– = r1p4q1p1q4p2q3p3q2– ++ = r2p4q2p1q3– p2q4q3q1++ = r3p4q3p1q2p2q1– p3q4++ =

PAGE 126

110[A-14] The last part is the same representation present in the table of quaternion operations above. Notice how only the signs change when we reverse the order of the multiplication in the matrix-vector form of the quaternion multiplication. This is because of the cross product relationship with the vector part: Rotations of Rigid Bodies in Space. Quaternions are able to represent a unique rotation in space. To perform a rotation () of a rigid body about an arbitrary moving/fixed axis (e) in space, the quaternion representation of this operation is: [A-15] [A-16] Notice that only one sine and one cosine function call is needed to calculate a quaternion, where euler would require three sine and three cosine function calls, one each for roll, pitch, and yaw. Since trigonometric function calls are computationally expensive, this is a great savings. Rotation of a point in space with a quaternion is done as follows: [A-17] [A-18] p q r1r2r3r4p4p3– p2p1p3p4p –1p2p –2p1p4p3p –1p –2p –3p4q1q2q3q4q4q3q –2q1q –3q4q1q2q2q –1q4q3q –1q –2q –3q4p1p2p3p4=== qp pq – = q13 –e 2 -sin = q4 2 -cos = e e1e2e3 T= PbqPaq = Paxyz 0T=

PAGE 127

111Since the quaternion has four elements, the point P must have a zero appended on to it in place of the real part of the quaternion. Thus the point P values constitute the imaginary (vector) part of the quaternion. This rotation can also be done similar to euler and fixed angle rotations by creating a rotation matrix from the quaternions. The attitude matrix (A) or rotation matrix (R) for this is [A-19] Successive rotations can be accomplished by multiplying the attitude matrices together. Hopefully it can be seen that quaternions are better than other Euler operations to determine attitude. Quaternions lack singularities due to one redundant element. They also lack the computationally intensive trigonometric functions, and contain a simplified way to determine successive rotations about an arbitrary axis. Attitude Errors in Quaternions Since quaternions can represent attitudes and not just rotations, it is important to be able to calculate the difference between where you want to point and where you are pointing. For example, say you had to control where a satellite was to point and you need an error term to feed into your control system. [A-20] Aq () R q1 2q2 2q3 2– – q4 2+ 2 q1q2q3q4+ () 2 q1q3q2– q4() 2 q1q2q3q4– () q1 2– q2 2q3 2– q4 2++ 2 q2q3q1q4+ () 2 q1q3q2q4+ () 2 q2q3q1q4– () q –1 2q2 2– q3 2q4 2++ == qEqS 1 –qTqT 4qT 3q –T 2qT 1q –T 3qT 4qT 1qT 2qT 2q –T 1qT 4qT 3q –T 1q –T 2q –T 3qT 4q –S 1q –S 2q –S 3pS 4==

PAGE 128

112where is the quaternion error, is the attitude of the satellite, and is the target attitude. Note that the quaternions are in reversed order in the matrix and the spacecraft quaternion has the imaginary part negated due to the inversion. Summary of Quaternions 1.Quaternions have singularities at angles of rotation of 180 degrees, but this is generally not a problem. Typically if an application encounters a rotation greater than 180 degrees, it is shorter to go in the other direction. 2.The quaternion is a compact method of representing a unique rotation with only four elements with out the redundant information present in euler rotation’s nine elements. 3.Quaternions are difficult to visualize without extracting out the axis of rotation and rotation angle. Although the quaternion can be thought of as a sphere in 4D space, this does the majority of people little good, since they can not perceive a 4D surface. 4.The use of quaternions allows you to compute the shortest distance between two different attitudes. If lerp is used, the resulting velocity between the two attitudes will not be constant but more bell shaped. If slerp is used, the resulting velocity will be constant. This process is difficult using euler angles, which contain many singularities that must be avoided and the solution is not guaranteed to be the shortest distance. Also the solutions that result from euler angles are not uniformly distributed along the surface of a sphere. Thus it again is not always possible to transition from one orientation to another using the shortest distance, because these empty areas must be skirted around to arrive at the desired location. qEqSqT

PAGE 129

113APPENDIX B MODAL EQUATIONS Motion of Multi-Degree of Freedom System Looking a the simple multi-degree of freedom (MDOF) system in Figure B-1, we can write the equation of motion as [B-1] Assuming [B-1] will have the simple harmonic solution [B-2] where is the displacement and is the amplitude of the nth mass. Thus substituting [B-2] into [B-1] gives [B-3] Mx Kx + f = f 0 = xnt () Anej t= xnAn Ma – Ka + 0 = a A1A2… AN= 2= Figure B-1. Simple multi-degree of freedom system composed of mass and spring elements.

PAGE 130

114Now dividing through by the mass and simplifying the equation a little we have the characteristic equation. [B-4] This is assuming that The roots of this equation are eigenvalues which provide the resonance frequency of the system. The roots must be real and positive if M and K are positive definite. Since the eigenvalues provide the resonance frequencies, the eigenvectors () provide the mode shapes of the system. Thus the nth mode shape of the system is [B-5] The mode shapes are arbitrary to a constant value, thus the modal vector is typically normalized. Two common methods to scale the mode shape are by using the mass matrix or by setting the maximum value of the modal vector to one. [B-6] The mode shapes also follow the following property of orthogonality. [B-7] where is the Kronecker delta function which is defined as [B-8] I K M ---– 0 = a 0 nA1 nA2 n… AN n= TM 1 = max 1 = TM nmMMn== TK nmKKn== nm0 nm 1 n m = =

PAGE 131

115The mode shapes form a basis which can be used to solve for x. To accomplish this, we will now define a modal matrix () whose columns are made up of modal vectors. The solution to the system is now assumed to take the form [B-9] where q(t) is a generalized coordinate or the modal displacement. Substituting [B-9] into [B-1] and pre-multiplying by the transpose of the of the modal matrix resulting in the equations of motion. [B-10] This can be simplified to: [B-11] This can further be manipulated by dividing the equation through by the mass matrix. This results in an equation with a matrix of natural frequencies and the modal force is now defined and inserted in the equation. [B-12] [B-13] The initial position and velocity of the system are defined by: [B-14] Dampening Throughout these derivations, we have assumed that there was no dampening in the system. The dampening coefficient matrix is assumed to be a linear combination of the mass and spring constant coefficient matrix. [B-15] xt () qt () = TM q TK q + f = Mnq Knq + f = q n 2q + Fn= FnMn 1 – f = qn0 () n TMx 0 () = q n0 () n TMx 0 () = Mnq Cnq Knq ++ f =

PAGE 132

116[B-16] or represented another way: [B-17] [B-18] Simple Example Thus for a system with two modes, the state space representation would be: [B-19] Again notice that the coefficient matrices are diagonal due to orthogonality of the modes in the system. Also note that the input forces for the nth mode only effects the nth mode. In order to derive this system, only the diagonal mass matrix, natural frequencies, and the normal modes need to be calculated. Cn Mn Kn+ = q Cnq n 2q ++ Fn= Cnn 2+ = q 1q 2C10 0 C2q 1q 21 20 0 2 2q1q2++ F1F2=

PAGE 133

117APPENDIX C SIMULATION CODE This Appendix will describe the satellite simulation, which was eventually written in the C++ programming language. Simulation Code Matlab The simulation originally was created in Matlab to utilize the quick development environment. The main problem with Matlab, however, was the slow speed at which the simulations ran. Thus a decision was made to port the code over to C to gain the advantage of speed. C Programming Language All of the mathematical benefits of Matlab had to be created for use in C. The book Numerical Recipes in C [73] was a great aid in developing some of the vector and matrix operations. The functions written to perform the various operations were designed to be as efficient and fast as possible. Thus most mathematical function required two inputs and an output, so that processor time was not wasted making expensive calls to malloc to create a return value. The C programming language produced fast simulations, however the dynamics equation code was difficult to read and follow due to the functional architecture of C. For example of how C’s functional aspects obscure the mathematics, see the example below.1 #include “mathlib.h” 2

PAGE 134

118 3 int main(void){ 4 vector_t h_dot, w, h, tmp, u; // create vectors 5 6 initVector(h_dot,3); // initialize the vectors to size 3 7 ... // initialize remaining vectors 8 9 // Euler’s equation: h_dot = w x h + u 10cross(w,h,tmp); // calc cross product and store in tmp vector 11 add(tmp,u,h_dot); // add the two vectors together and store in h_dot 12 13}Although Euler’s equation is calculated properly using the two functions cross() and add(), where the first two arguments are the inputs and the third is the output, this is not as clear as Matlab code and therefor can lead to errors. Thus the code was moved from C to C++ which is an object oriented language where concepts such as vectors and matrices can exist. C++ Programming Language C++ provided cleaner mathematical code without sacrificing performance. This was due to C++’s object oriented nature and the ability to overload operators such as: +, -, *, and /. Basically the C code previously written was modified for C++ by following some of the examples in [72]. An example of the previous C code rewritten in C++ is shown below.1 #include “mathlib.h” 2 3 int main(void){ 4c V ector h_dot(3), w(3), h(3), u(3); // create vectors of size 3 5 6 // Euler’s equation: h_dot = w x h + u 7 h_dot = cross(w,h) + u; 8 9} 10Notice that the temporary vector tmp is no longer required. This is due to the fact that the math library developed here automatically creates and reuses temporary vectors and

PAGE 135

119matrices for mathematical operations. This means that the new operator1does not need to be called to return temporary values. This is accomplished by creating a cache of vectors, matrices, and quaternions that are need using an STL (standard template library) data structure called a vector2. Further explanation of STL can be found in [75], which covers programing with the wide variety of data structures (i.e., linked lists, hash tables, trees, etc.) available in STL. cMathlib Header File This header file provides the foundation of the programming work. Here there are definitions for vectors (cVector), matrices (cMatrix), quaternions (cQuaternion), and error handling (cMLError). 1.The new operator is equivalent to malloc in C. For example, double *a = new double in C++ is basically the same as double *a = (double*)malloc(sizeof(double)) in C. 2.Note the STL vector is not a mathematical vector, but rather a dynamic array data structure.Figure C-1. UML diagram of Mathlib class which contains C++ objects for vectors, matricies, and quaternion used in this work.

PAGE 136

120The cVector, cMatrix, and cQuaternion are mathematical objects and are easy to comprehend. The cMLError or class Math Library Error is an object that is created when an error occurs. The C++ language has a method of error handling called catch and throw. Basically when an error occurs in a function, an error (in this case cMLError) is thrown by the offending piece of code. The code that called the function then catches the error and handles the problem appropriately. Further information on catch and throw can be found in Deitel and Deitel [74].1 #ifndef CMATHLIB_H 2 #define CMATHLIB_H 3 4 //--C++ --5 #include 6 #include 7 #include 8 9 //--C --10#include 11 #include 12#include 13 14//--defines --15#define ML_CHECK 16#define ERROR_STRING_SIZE 255 17#define ml_data double 18 19#define R2D (180.0/M_PI) 20#define D2R (M_PI/180.0) 21 22class cBaseMath; 23class cMatrix; 24class cQuaternion; 25class cVector; 26template < class T > class cArray; 27 28class cMLError { 29 public: 30 cMLError(int,char*); 31 ~cMLError(void); 32 friend std::ostream &operator<<(std::ostream&,cMLError&); 33 void operator+=(cMLError&); 34 void operator+=(char*); 35 36 enum {NONFATAL,FATAL}; 37 char errMsg[255];

PAGE 137

121 38 int type; 39}; 40 41 42/*! 43 \warning NEVER pass a temp matrix, vector, or 44 quaternion to or from a function. A 45 temp should ALWAYS be passed to a user 46 created matrix, vector, or quaternion. 47*/ 48class cBaseMath { 49 public: 50 cBaseMath(void); 51 ~cBaseMath(void); 52 static cQuaternion* getTmp(void); 53 static cVector* getTmp(int); 54 template static cArray* getTmpArray(int); 55 static cMatrix* getTmp(int,int); 56 void releaseTmp(void) const; 57 inline void setKey(int a){key = a;} 58 inline int getKey(){return key;} 59 void setName(char*); 60 char* getName(void){return name;} 61 inline bool valid(void){return p == NULL ? false : true;} 62 inline int getType(void){return type;} 63 ml_data get(const int i) const; // {return p[i];} 64 ml_data& set(const int i); //{return p[i];} 65 66 //protected: 67 enum {ML_VECTOR,ML_MATRIX,ML_QUATERNION,ML_ARRAY}; 68 static std::vector cache; 69 int key; 70 int type; 71 bool available; 72 ml_data *p;// pointer to data 73 char *name; 74 static int keyNum; 75}; 76 77 78template < class T > 79cArray* cBaseMath::getTmpArray(int size){ 80 cBaseMath *b = NULL; 81 cArray *a = NULL; 82 unsigned int index = 0; 83 84 for(index = 0; indextype == cBaseMath::ML_ARRAY && b->available ){ 87 a = (cArray*) b; 88 if (a->getSize() >= size){ 89a->available = false; 90//printf(" using array[%d](%d)\n",a->getKey(),size); 91break; // get out of for loop

PAGE 138

122 92 } 93 else{ 94a = NULL; 95 } 96 } 97 } 98 99 if( a == NULL ){ 100 a = new cArray(size); 101 a->setKey(keyNum++); 102 a->available = false; 103 //printf(" + array[%d](%d)\n",a->getKey(),size); 104 cache.push_back(a); 105 } 106 107 memset(a->array,0,sizeof(T)*size); 108 109 return a; 110} 111 11 2 /////////////////////////////////////////////// 113 11 4 template < class T > 11 5 class cArray : public cBaseMath { 11 6 public: 11 7 cArray(int); 11 8 ~cArray(void); 11 9 inline int& operator[](int a){return array[a];} 120 inline int getSize(void){return size;} 121 //protected: 122 T *array; 123 int size; 124}; 125 126template < class T > 127cArray::cArray(int s){ 128 type = ML_ARRAY; 129 key = 0; 130 size = s; 131 132 array = NULL; 133 array = new T[size]; 134 135 if(array == NULL){ 136 throw cMLError(cMLError::FATAL,"could not create array"); 137 } 138} 139 140template < class T > 141cArray::~cArray(void){ 142 delete array; 143 size = 0; 144} 145

PAGE 139

123 146 147 148/////////////////////////////////////////////// 149 150class cVector : public cBaseMath { 151 friend cVector& operator+(const cVector&, const cVector&); 152 friend cVector& operator-(cVector&, cVector&); 153 friend cVector& operator-(cVector&); 154 friend cVector& operator*(cVector&, ml_data); 155 friend cVector& operator*(ml_data, cVector&); 156 friend cVector& operator/(cVector&, ml_data); 157 friend cVector& operator/=(cVector&,ml_data); 158 friend cVector& operator*=(cVector&,ml_data); 159 friend cVector& operator*=(ml_data,cVector&); 160 friend cVector& operator+=(cVector&, cVector&); 161 friend cVector& operator-=(cVector&, cVector&); 162 friend std::ostream& operator<<(std::ostream&,cVector&); 163 friend ml_data dot(cVector&,cVector&); 164 friend cVector& cross(cVector&,cVector&); 165 friend cMatrix& outer(cVector&,cVector&); 166 167 public: 168 cVector(int,char* =NULL,ml_data* =NULL); 169 cVector(void){size=0;name=NULL;p=NULL;} 170 cVector(cVector&); 171 ~cVector(); 172 cVector& operator=(cVector&); 173 cVector& operator=(cVector*); 174 cVector& operator=(ml_data*); 175 inline bool operator==(cVector &a){return (size==a.size ? true : false);} 176 inline bool operator!=(cVector &a){return (size==a.size ? false : true);} 177 178 ml_data& operator()(int)const; 179 ml_data norm(void); 180 cVector& ones(void); 181 182 inline void set(ml_data *d){memcpy(p,d,sizeof(ml_data)*size);} 183 inline void set(ml_data a, ml_data b, ml_data c){p[0]=a;p[1]=b;p[2]=c;} 184 inline void clear(void){memset(p,0,sizeof(ml_data)*size);} 185 inline int getSize(void){return size;} 186 void resize(int); 187 188 //private: 189 int size; 190}; 191 192class cMatrix : public cBaseMath { 193 friend cMatrix& operator+(cMatrix&,cMatrix&); //matrix addition 194 friend cMatrix& operator-(cMatrix&,cMatrix&); //matrix subtraction 195 friend cMatrix& operator-(cMatrix&); // negate a matrix 196 friend cMatrix& operator*(cMatrix&,cMatrix&); //matrix multiplication 197 friend cMatrix& operator*(ml_data,cMatrix&); // scalar multiplication 198 friend cMatrix& operator*(cMatrix&,ml_data); // scalar multiplication 199 friend cVector& operator*(cMatrix&,cVector&); //matrix addition

PAGE 140

124 200 friend cMatrix& operator/(ml_data,cMatrix&); // scalar division 201 friend cMatrix& operator/=(cMatrix&,ml_data); // scalar division 202 friend cMatrix& operator+=(cMatrix&,cMatrix&); //matrix addition 203 friend cMatrix& operator-=(cMatrix&,cMatrix&); //matrix subtraction 204 friend cMatrix& operator*=(cMatrix&,cMatrix&); //matrix multiplication 205 friend cMatrix& skew(ml_data*,int=3); 206 friend cMatrix& eye(int); 207 friend cMatrix& zeros(int,int); 208 friend cMatrix& ones(int,int); 209 210 friend std::ostream &operator<<(std::ostream &s, cMatrix &mat); // print matrix 211 212 public: 213 cMatrix(int,int,char* =NULL,ml_data* =NULL);// constructor 214 cMatrix(void){r=c=0;name=NULL;p=NULL;} 215 ~cMatrix(); // destructor 216 cMatrix(const cMatrix &m);// copy constructor 217 218 cMatrix& eye(void); 219 cMatrix& ones(void); 220 cMatrix& diag(ml_data*); 221 cMatrix& inv(void); 222 cMatrix& trans(void); 223 cMatrix& skew(ml_data*); 224 ml_data trace(void); 225 void resize(int,int); 226 inline void reset(int a, int b, char *n =NULL){resize(a,b);setName(n);} 227 228 ml_data& operator()(int i, int j)const;// index 229 cMatrix& operator=(cMatrix&); 230 cMatrix& operator=(cMatrix*); 231 cMatrix& operator=(ml_data*); 232 cMatrix& fill(int,int,int,int,ml_data*); 233 inline bool operator==(cMatrix &a){return (r==a.r && c==a.c ? true : false);} 234 inline bool operator!=(cMatrix &a){return (r==a.r && c==a.c ? false : true);} 235 236 inline int getRow(void){return r;} 237 inline int getCol(void){return c;} 238 void setR(ml_data*,int); 239 void setC(ml_data*,int); 240 inline void set(cMatrix &m,int a, int b){set(m.p,a,b,m.r,m.c);} 241 inline void set(ml_data *d){memcpy(p,d,sizeof(ml_data)*r*c);} 242 inline void set(ml_data *d,int rr,int cc){set(d,rr,cc,r,c);} 243 void set(ml_data*,int,int,int,int); 244 inline void clear(void){memset(p,0,sizeof(ml_data)*r*c);} 245 246 //private: 247 int r,c; // index ranges 248}; 249 250 251class cQuaternion : public cBaseMath { 252 friend cQuaternion& operator*(cQuaternion&,cQuaternion&); 253 friend inline cQuaternion& inv(cQuaternion &q){return q.inv();}

PAGE 141

125 254 friend cQuaternion& slerp(cQuaternion&,cQuaternion&,int); 255 friend std::ostream &operator<<(std::ostream&, cQuaternion&); // print quaternion 256 //friend void operator=(ml_data*,cQuaternion&); 257 258 public: 259 cQuaternion(char*); 260 cQuaternion(const cQuaternion&); // copy constructor 261 262 // set data 263 inline void set(ml_data a, ml_data b, ml_data c, ml_data d){ p[0]=a; p[1]=b; p[2]=c; p[3]=d; } 264 inline void set(ml_data *a){set(a[0],a[1],a[2],a[3]);} 265 inline ml_data& operator[](int i){return p[i];} // 0 based 266 inline ml_data& operator()(int i){return p[i-1];} // 1 based 267 cQuaternion operator*=(cQuaternion&); 268 void operator=(cQuaternion&); 269 void operator=(cQuaternion*); 270 void operator=(ml_data*); 271 272 // get data 273 inline ml_data& x(void){return p[0];} 274 inline ml_data& y(void){return p[1];} 275 inline ml_data& z(void){return p[2];} 276 inline ml_data& r(void){return p[3];} 277 inline ml_data* getImg(void){return p;} 278 inline ml_data getReal(void){return p[3];} 279 ml_data* getAxis(void); 280 ml_data getAngle(void); 281 ml_data norm(void); 282 void copy(ml_data*); 283 284 // misc 285 cQuaternion& inv(void); 286 cQuaternion& conj(void); 287 void normalize(void); 288 void clear(void); 289 void error(cQuaternion&,cQuaternion&); 290 291 // conversions 292 void e2q(ml_data,ml_data,ml_data,int); 293 inline void e2q(ml_data *e,int t){e2q(e[0],e[1],e[2],t);} 294 ml_data* q2e(int=7); 295 ml_data* q2R(int); 296 297 private: 298 299}; 300 301 302#endif 303 304 305

PAGE 142

126cRK4 Header File This code was responsible for integrating the equations of motion using a Rung-Kutta integration method. Notice that it utilizes the cVectors from cMathlib shown above.1 #ifndef KEVINS_RUNGE_KUTTA 2 #define KEVINS_RUNGE_KUTTA 3 4 #include 5 #include "cMathlib.h" 6 7 8 class cRK4 { 9 public: 10cRK4(cVector& (*model)(cVector&,cVector&,cVector&), 11 ml_data dt, 12int x_size, 13int u_size 14); 15 16~cRK4(void); 17cVector& integrate(void); 18 19//protected: 20ml_data time,dt; 21cVector k1x,k2x,k3x,k4x,x,xx,y,u,dist; 22cVector& (*model)(cVector&,cVector&,cVector&); 23}; 24 25 26#endif 27 28

PAGE 143

127REFERENCES [1]Kang, W., Sparks, A., and Banda, S., “Coordinated Control of Multisatellite Systems,” Journal of Guidance, Control, and Dynamics, Vol. 24, No. 3, March-April 2001, pp. 360-368. [2]Mesbahi, M. and Hadaegh, F., “Formation Flying Control of Multiple Spacecraft via Graphs, Matrix Inequalities, and Switching,” Journal of Guidance, Control, and Dynamics, Vol. 24, No. 3, March-April 2001, pp. 369-377. [3]Chichka, D., “Satellite Clusters with Constant Apparent Distribution,” Journal of Guidance, Control, and Dynamics, Vol. 24, No. 1, January-February 2001, pp. 117122. [4]Kapila, V., Sparks, A. G., Buffington, J., and Yan, Q., “Spacecraft Formation Flying Dynamics and Control,” Journal of Guidance, Control, and Dynamics, Vol. 23, No. 3, May-June 2000, pp. 561-564. [5]Queiroz, M. S., Kapila, V., and Yan, Q., “Adaptive Nonlinear Control of Multiple Spacecraft Formation Flying,” Journal of Guidance, Control, and Dynamics, Vol. 23, No. 3, May-June 2000, pp. 385-390. [6]Ocampo, C., and Rosbrough, G., “Multiple-Spacecraft Orbit Transfer Problem: W eak-Booster and Strong-Booster Cases,” Journal of Guidance, Control, and Dynamics, Vol. 22, No. 6, November-December 1999, pp. 897-904. [7]Ocamp, C., and Rosborough, W., “Multiple-Spacecraft Orbit-Transfer Problem: The No-Booster Case,” Journal of Guidance, Control, and Dynamics, Vol. 22, No. 5, September-October 1999, pp. 650-657. [8]Wang, P., Hadaegh, F., and Lau, K., “Sychronized Formation Rotation and Attitude Control of Multiple Free-Flying Spacecraft,” Journal of Guidance, Control, and Dynamics, Vol. 22, No. 1, January-February 1999, pp. 28-35. [9]Sidi, M., Spacecraft Dynamics and Control Cambridge University Press, New York, 1997. [10]Agrawal, B., “Dynamic Characteristics of Liquid Motion in Partially Filled Tanks of a Spinning Spacecraft,” Journal of Guidance, Control, and Dynamics, Vol. 16, No. 4, July-August 1993, pp. 636-640.

PAGE 144

128[11]Genevaux, J. and Lu, D., “On Sloshing in a Container with Moving Wall,” Journal of Fluids and Structures, Vol. 14, 2000, pp. 275-278. [12]Anderson, J., Turan, O., and Semercigil, S., “Experiments to Control Sloshing in Cylindrical Containers,” Journal of Sound and Vibration, Vol. 240, No. 2, 2001, pp. 398-404. [13]Anderson, J., Turan, O., and Semercigil, S., “A Standing-wave-type Sloshing Absorber to Control Transient Oscillations,” Journal of Sound Vibrations, Vol. 232, No. 5, 2000, pp. 839-856. [14]Cho, S. and McClaroch, N., “Feedback Control of a Space Vehicle with Unactuated Fuel Slosh Dynamics,” AIAA Guidance, Navigation, and Control Conference, Denver, CO, August 14-17, 2000. [15]Kuang, J. and Leung, A. Y., “ Feedback for Attitude Control of Liquid-Filled Spacecraft,” Journal of Guidance, Control, and Dynamics, Vol. 24, No. 1, JanuaryFebruary 2001, pp. 46-55. [16]Grundelius, M. and Bernhardsson, B., “Control of Liquid Slosh in a Industrial Packing Machine,” IEEE Conference on Control Applications -Proceedings, Vol. 2, Kohala Coast, HI, 1999, pp. 1654-1659. [17]Wisniewski, R. and Blanke, M., “Fully Magnetic Attitude Control for Spacecraft Subject to Gravity Gradient,” Automatica, Vol. 35, 1999, pp. 1201-1214. [18]Chen, X., Steyn, W. H., Hodgart, S., and Hashida, Y., “Optimal Combined ReactionWheel Momentum Management for Earth-Pointing Satellites,” Journal of Guidance, Control, and Dynamics, Vol. 22, No. 4, July-August 1999, pp. 543-550. [19]Agrawal, B. N., “Dynamic Characteristics of Liquid Motion in Partially Filled Tanks of a Spinning Spacecraft,” Journal of Guidance, Control, and Dynamics, Vol. 16, No. 4, July-August 1993, pp. 636-640. [20]Hoskins, W., Cassady, R., Campbell, M., and Rayburn, C., “A Micro Pulsed Plasma Thruster (PPT) for the "Dawgstar" Spacecraft,” IEEE Aerospace Conference Proceedings, Vol. 4, Big Sky, MT, 2000, pp. 7-13. [21]Wilson, J., “Major New Thrust for MEMS Engines,” Aerospace America, February, 2003, pp.34-38. [22]Wie, B., Space Vehicle Dynamics and Control AIAA Education Series, AIAA Inc., Reston, 1998. H

PAGE 145

129[23]W. Singhose, K. Bohlke, and W. Seering, “Fuel-Efficient Shaped Command Profiles for Flexible Spacecraft,” AIAA Guidance, Navigation, and Control Conference. Baltimore, MD, 1995. [24]Ooten, E. and Singhose, W., “Command Generation for Flexible Systems using Numerator Dynamics and Sliding Mode Control”, Master’s Thesis, Georgia Institute of Technology, 2000. [25]Wie, B. and Liu, Q., “Classical and Robust Controller Redesign for the Hubble Space Telescope,” Journal of Guidance, Control, and Dynamics, Vol. 16, No. 6, November-December 1993, pp. 1069-1077. [26]Thorten, E., Thermal Structures for Aerospace Applications AIAA Education Series, AIAA Inc., Reston, 1996. [27]Walchko, K. and Mason, P. A. C., “Development of a Generic Hybrid Fuzzy Controller,” Artificial Neural Networks in Engineering, St. Louis, Ms, 1-4 Nov.1998. [28]Petroff, N., Mason, P., Walchko, K., “Numerical Stability Analysis of a Fuzzy Controller,” Artificial Neural Networks in Engineering, St. Louis, Ms, 1-4 Nov.1998. [29]Vreeberg, J., “Acceleration Measurements on Sloshsat FLEVO for Liquid Force Location Determination,” 4th ESA International Conference on Spaceflight Guidance, Navigation, and Control Systems, Noordwijk, Netherlands, 18-21 October 1999. [30]El-Sayad, M., Hanna, S., and Ibrahim, R., “Parametric Excitation of Nonlinear Elastic Systems Involving Hydrodynamic Sloshing Impact,” Nonlinear Dynamics, Vol. 18, 1999, pp. 25-50. [31]Hughes, P., Spacecraft Attitude Dynamics John Wiley & Sons, New York, 1986. [32]Greenburg, M., Advanced Engineering Mathematics Prentice Hall, New York, 1998, pp. 1033-1034. [33]Kreyszig, E., Advanced Engineering Mathematics John Wiley & Sons, New York, 1999, pp. 598-599. [34]Meirovitch, L., Analytical Methods in Vibrations Collier-MacMillian Limited, New Y ork, 1967. [35]Zill, D. and Cullen, M., Advanced Engineering Mathematics Jones and Bartlett Publishers International, New York, 2000, pp. 703-704. [36]McConnell, K., V ibration Testing, Theory and Practice, John Wiley and Sons, New Y ork, 1995. H

PAGE 146

130[37]Thompson, W., Theory of Vibrations with Applications Prentice Hall, New York, 1998. [38]Wertz, J., Spacecraft Attitude Determination and Control D. Reidel Publishing Co., Boston, 1978. [39]Meirovitch, L., Methods of Analytical Dynamics McGraw-Hill Book Co., New York, 1970. [40]Miller, A., Gray, G., and Massoleni, A., “Nonlinear Spacecraft Dynamics with a Flexible Appendage, Dampening, and Moving Internal Submasses,” Journal of Guidance, Control, and Dynamics, Vol. 24, No. 3, May-June 2001, pp. 605-615. [41]Wie, B., and Barba, P., “Quaternion Feedback for Spacecraft Large Angle Maneuvers,” Journal of Guidance, Control and Dynamics, Vol. 8, No. 3, May-June 1985, pp. 360-365. [42]Slotine, J. and Li, W., Applied Nonlinear Control Pearson Education, New York, 1990. [43]DeCarlo, R., Zak, S., and Mathews, G., “Variable Structure Control of Nonlinear Multivariable Systems: A Tutorial,” Proceeding of the IEEE, Vol. 76, No. 3, March 1988, pp. 212-232. [44]Zhou, D., Mu, C., and Xu, W., “Adaptive Sliding-Mode Guidance of a Homing Missile,” Journal of Guidance, Control, and Dynamics, Vol. 22, No. 4, July-August 1999, pp. 589-594. [45]Lewis, A. S. and Sinha, A., “Sliding Mode Control of Mechanical Systems with Bounded Disturbances via Output Feedback,” Journal of Guidance, Control, and Dynamics, Vol. 22, No. 2, March-April 1999, pp. 235-240. [46]Erbatur, K., Kaynak, O., Sabanovic, A., and Rudas, I., “Fuzzy Parameter Adaptation for a Sliding Mode Controller as Applied to the Control of an Articulated Arm,” Proceedings of the 1997 IEEE International Conference on Robotics and Automation, Albuquerque, NM, April 1997, pp. 817-822. [47]Hsu, F. and Fu, L., “Nonlinear Control of Robot Manipulators Using Adaptive Fuzzy Sliding Mode Control,” International Conference on Intelligent Robots and Systems, V ol. 1, Pittsburgh, PA, August 5 9, 1995, pp. 156 161. [48]Bekit, B., Whidborne, J., and Seneviratne, L., “Fuzzy Sliding Mode Control for a Robot Manipulator,” The 1997 IEEE International Symposium on Computational Intelligence in Robotics and Automation, Monterey, CA, July 10-11, 1997, pp. 320325.

PAGE 147

131[49]Kiriazov, P., Kreuzer, E., and Pinto, F. C., “Robust Feedback Stabilization of Underwater Robotic Vehicles,” Robotics and Autonomous Systems, Vol. 21, 1997, pp. 415423. [50]Song, G., Longman, R., and Mukherjee, R., “Integrated Sliding-Mode AdaptiveRobust Control,” IEE Proceedings of Control Theory Applications, Vol. 146, No. 4, July 1999, pp.341-347. [51]Conte, G. and Serrani, A., “Robust Control of a Remotely Operated Underwater Vehicle,” Automatica, Vol. 34, No. 2, 1998, pp. 193-198. [52]Alessandri, A., Hawkinson, T., Healey, A., and Veruggio, G., “Robust Model-Based Fault Diagnosis for Unmanned Underwater Vehicles Using Sliding Mode-Observers,” 1 1th International Symposium on Unmanned Untethered Submersible Technology, Autonomous Undersea Systems Institute, Durham, NH, Aug. 22-25, 1999, pp. 1-8. [53]Livneh, R. and Wie, B., “Asymmetric Body Spinning Motion with Energy Dissipation and Constant Body-Fixed Torques,” Journal of Guidance, Control, and Dynamics, Vol. 22, No. 2, March-April 1999, pp. 322-328. [54]Tsiotras, P., Haijun, S., and Hall, C., “Satellite Attitude Control and Power Tracking with Energy/Momentum Wheels,” Journal of Guidance, Control, and Dynamics, Vol. 24, No. 1, January-February 2001, pp. 23-35. [55]Collins, E. G. and Song, T., “Robust Estimation and Fault Detection of Uncertain Dynamic Systems,” Journal of Guidance, Control, and Dynamics, Vol. 23, No. 5, September-October 2000, pp. 857-864. [56]Fowell, R. A., Milford, R. I., and Yocum, J. F., “Spacecraft Spin Stabilization Using a T ransverse Wheel for Any Inertia Ratio,” Journal of Guidance, Control, and Dynamics, Vol. 22, No. 6, November-December 1999, pp. 768-775. [57]Kumar, K. K., and Nishita, K., “Robustness Analysis of Neural Networks with Application to System Identification,” Journal of Guidance, Control, and Dynamics, Vol. 22, No. 5, September-October 1999, pp. 695-701. [58]Shen, H. and Tsiotras, P., “Time-Optimal Control of Axisymmetric Rigid Spacecraft Using Two Controls,” Journal of Guidance, Control, and Dynamics, Vol. 22, No. 5, September-October 1999, pp. 682-694. [59]Crassidis, J. L., “Robust Control of Nonlinear Systems Using Model-Error Control Synthesis,” Journal of Guidance, Control, and Dynamics, Vol. 22, No. 4, July-August 1999, pp. 595-601. H

PAGE 148

132[60]Lo, S. and Chen, Y., “Smooth Sliding-Mode Control for Spacecraft Attitude Tracking Maneuvers,” Journal of Guidance, Control, and Dynamics, Vol. 18, No. 6, NovemberDecember 1995, pp. 1345-1348. [61]Boskovic, J. D., Li, S., and Mehra, R. K., “Robust Adaptive Variable Structure Control of Spacecraft Under Control Input Saturation,” Journal of Guidance, Control, and Dynamics, Vol. 24, No. 1, January-February 2001, pp. 14-22. [62]Boskovic, J. D., Li, S., and Mehra, R. K., “Robust Adaptive Variable Structure Control of Spacecraft Under Control Input Saturation,” Journal of Guidance, Control, and Dynamics, Vol. 24, No. 1, January-February 2001, pp. 14-22. [63]Crassidis, J. L., Vadali, S. R., and Markley, F. L., “Optimal Variable-Structure Control T racking of Spacecraft Maneuvers,” Journal of Guidance, Control, and Dynamics, V ol. 23, No. 3, May-June 2000, pp. 564-566. [64]Coleman, C. P. and Godbole, D., “A Comparison of Robustness: Fuzzy Logic, PID, & Sliding Mode Control,” Proceedings of the 3rd IEEE International Conference on Fuzzy Systems, pp. 1654-1560. [65]Yousefi-Koma, A. and Vukovick, G., “Vibration Suppression of Flexible Beams with Bonded Piezotransducers Using Wave-Absorbing Controllers,” Journal of Guidance, Control, and Dynamics, Vol. 23, No. 2, March-April 2000, pp. 347-354. [66]Reitz, R., Nonlinear Feedback Control of a Slewing Flexible Structure Master’s Thesis, State University of New York at Buffalo, 1992. [67]Reitz, R., Dynamics and Control of Thermally Induced Vibrations Ph.D. Dissertation, State University of New York at Buffalo, 1996. [68]Johnston, J. and Thorten E., “Thermally Induced Dynamics of Satellite Solar Panels,” Journal of Spacecraft and Rockets, Vol. 37, No. 5, September-October, 2000, pp. 604613. [69]Incropera, F. and DeWitt, D., Fundamentals of Heat and Mass Transfer, 4th Ed. John W iley & Sons, New York, 1996. [70]Brown, R. and Hwang, P., Introduction to Random Signals and Applied Kalman Filtering John Wiley & Sons, New York, 1997. [71]Tonque, B., Principles of Vibrations Oxford University Press, New York, 1996. [72]Capper, D., Introducing C++ for Scientists, Engineers, and Mathematicians Springer-Verlag, New York, 1994.

PAGE 149

133[73]Press, H., Teukolsky, S., Vetterling, W., and Flannery, B., Numerical Recipes in C, 2nd Ed. Cambridge University Press, New York, 1992. [74]Deitel, H. and Deitel, P., C++ How to Program, 2nd Ed., Prentice Hall, Upper Saddle River, New Jersey, 1998. [75]Josuttis, N., The C++ Standard Template Library, A Tutorial and Reference AddisonW esley, New York, 1999.

PAGE 150

134BIOGRAPHICAL SKETCH Kevin J. Walchko was born on 1 November 1972 in Sarasota, FL, and attended the University of Florida starting in the fall of 1991. After several years, a B.S. degree in mechanical engineering was earned in 1997. His college was paid for through the Montgomery G.I. Bill during his 7 years of service in the Florida Army National Guard as a Light Cavalry Scout. During his military career, Kevin participated in humanitarian relief for hurricane Andrew and Opal. Kevin next received a master’s degree in spring of 1999 for his work with NASA and satellite attitudes control. Most of the research Kevin conducted during this time was in controls and dynamics. Specifically, the research involved intelligent controls such as fuzzy logic and neural networks. Continuing his work with NASA, Kevin stayed at the University of Florida for a Ph.D. in mechanical engineering awarded in the spring of 2003. Kevin has also conducted extensive research in robotics and autonomous vehicles. He has been a main member of the autonomous submarine team, Subjugator, for several years. Working closely with several members of the Machine Intelligence Laboratory in

PAGE 151

135the Electrical and Computer Engineering Department at the University of Florida, he concurrently completed a master’s degree in electrical engineering during the summer of 2002 which involved low cost inertial navigation for autonomous vehicles. During his time in graduate school, Kevin also fell in love and finally married Nina C. Massie on 20 October 2001. Their first son, Michael Jacob Walchko, was born on 20 January 2003. Kevin’s future plans are to return to the military by joining the U.S. Air Force as a developmental engineer. Here he will conduct research into new technologies and supervise research projects for the Air Force.


Permanent Link: http://ufdc.ufl.edu/UFE0000818/00001

Material Information

Title: Robust Nonlinear Attitude Control With Disturbance Compensation
Physical Description: Mixed Material
Copyright Date: 2008

Record Information

Source Institution: University of Florida
Holding Location: University of Florida
Rights Management: All rights reserved by the source institution and holding location.
System ID: UFE0000818:00001

Permanent Link: http://ufdc.ufl.edu/UFE0000818/00001

Material Information

Title: Robust Nonlinear Attitude Control With Disturbance Compensation
Physical Description: Mixed Material
Copyright Date: 2008

Record Information

Source Institution: University of Florida
Holding Location: University of Florida
Rights Management: All rights reserved by the source institution and holding location.
System ID: UFE0000818:00001


This item has the following downloads:


Full Text











ROBUST NONLINEAR ATTITUDE CONTROL WITH DISTURBANCE COMPENSATION


By


KEVIN J. WALCHKO












A DISSERTATION PRESENTED TO THE GRADUATE SCHOOL
OF THE UNIVERSITY OF FLORIDA IN PARTIAL FULFILLMENT
OF THE REQUIREMENTS FOR THE DEGREE OF
DOCTOR OF PHILOSOPHY

UNIVERSITY OF FLORIDA


2003

























Copyright 2003

By

Kevin J Walchko















ACKNOWLEDGEMENTS

I would like to thank Dr. Paul Mason for his guidance of my education and research.

He has always made every effort to aid me in my work and life. Paul has held many roles

in my life such as mentor, academic advisor, friend, and best man at my wedding. He has

become more than an advisor over the years, and will forever be part of my family. I will

always be in his debt.

I would like to thank Dr. John Schueller who acted as my co-chair during this

dissertation. He provided much advice and guidance in the bureaucracy of University of

Florida. Schueller was always willing to help me with various paperwork and other such

matters pertaining to fellowships, pay checks, and graduation. He is a good friend and I

look forward to working with him in the future.

I would also like to thank my family who sacrificed much at times so this could be

possible. My wife Nina took great care of me and our new son Michael. She made sure I

had plenty of time to do my work and never complained about anything. She is the source

from which I draw all of my strength, ambition, and dedication.

Finally, I would like to thank the horde, our pack of cats. They provided much

amusement and stress relief through their wacky antics. I treasure each and every one of

them, and will never forget those who have passed away.















TABLE OF CONTENTS

page


A C K N O W L E D G E M E N T S .............................................................................................. iii

L IST O F FIG U R E S ............. ...................... .... .......................... .. .............. viii

LIST OF TABLES ......... ........................ ........ .......... ....... xi

L IST O F SY M B O L S ..... .... .................................................. .. ....... .............. xiii

A B ST R A C T ..............................................................................xv

1 IN TR O D U C TIO N .................. ...... ............................ .... ......... .. .............1

Form action Flying ................................... ..... .................. ................... 1
G general Problem s Controlling Single Satellites ........................................ ................2
External D isturbances ........................ .. .................... .. ....... .......... ...... .
Param etric U uncertainty ................................................. .............................. 3
L literature R review ................. ................................ .......... ....... ................ .4
Fuel Slosh ........................................ ..............
M odelling fuel slosh dynam ics .................................. ..................................... 4
C controlling fuel slosh.................. .......... ....................... ......... ... ..............
Do we need to worry about fuel slosh?..........................................................9
Missions involving fuel slosh disturbances ............................................. 12
Therm ally Induced Vibrations: Solar Snap ................................ ..................... 12
Controlling therm al vibrations...................................................................... 13
M missions w which involved solar snap ..................................... ............... 15
Param etric U uncertainty ............................................... .............................. 16
Sum m ary and Overview of the Thesis ........................................ ..... ............... 17

2 CONTROLLING THE UNKNOWN: FUEL SLOSH ..............................................19

In tro d u ctio n ......................................................................... 19
M odelling Fuel Slosh D ynam ics ........................................... ........... ............... 20
F uel Slosh M models ......................... ...... .............. ..................... .... .. ..23
Sim ulation M odel ............................................................. 23
C controller M odel ................................................................. .. ..... 24









3 D YN AM ICS OF SOLAR ARRAY S ................................................. .....................26

In tro d u ctio n ........................................................................................................... 2 6
Therm ally Induced Vibrations ............................................................................. 26
H e a t F lu x ............................................................................................................... 2 7
Sum m ing F forces ............. ........ ............................................ ........ ..... .......2 8
Nonhomogeneous Boundary Conditions ................................... ............... 33
Separation of V ariables ............................ ........... .............. ... .. ............... 34
Solving the position function F(x) ....................................... ............... 35
Solving the tim e function G (t) ........................................ ...... ............... 38
Su m m ary of resu lts ..................................................... .. ........ ........ ..... ... ... 3 8
Therm ally Induced Dynamics .................................. .....................................39
M odal Equations .......................... ........ .... .. ..... .. ............39
D isturbance Torques ........................ ........................ ......... ........... 40
Thickness Tem perature Gradient ........................................ ...................... 42
C on clu sion ............................................................................44

4 SATELLITE ATTITUDE DYNAM ICS ........................................... .....................46

Reference Fram es ................................... .. .......... .. ............46
E arth C entered In ertial ........................................ ............................................4 6
O rb ital F ram e ...............................................................4 7
Body Fixed Frame ....................................... .......... ................47
Fun With Vectors and Rotating Coordinate Systems ...............................................47
Spacecraft Attitude ................................................ ..............49
Q u aternion s ..................................................... ................. 50
Rotations of Rigid Bodies in Space. ............................... .. ....................... 51
Spacecraft Equations of M otion ......................... .................... .. ............ ... 51
A attitude K inem atics ....................................... ................... ......... 52
Spacecraft D ynam ics ....................... .. ....................... .... .. ........... 52
Sim ple Spacecraft ................................. ... .. ........ ............ 52
Angular Velocity .................................... .. .. ..... .. ............56
Internal D disturbances ............. .... ......... .......................... ............ ........ ......56
Spacecraft with Reaction or Momentum Wheels ..............................................56
F u el S lo sh ...................................................................................................... 6 0
Environm ental D isturbances in Space ............................................. ............... 61
Gravity Gradient ................................. .. .. .. ...... .. ............62
S o la r S n a p ....................................................... ................ 6 2
C on clu sion ............................................................................62

5 ATTITUDE CONTROL OF SPACECRAFT .................................... ...............64

Satellite C control H ardw are .............................................................. .....................64
Control Com putational H ardw are ........................................ ...... ............... 64
C control A ctu ators ............... .......................................... ...... .... ......... ........ 65
Reaction and m om entum wheels .......................................... ............... 65









M agnetic torque rods ............................................... ............................ 66
G a s jets ................................................................................ 6 7
Proportional D erivative Control .......................................... ........................... 67
Sliding Mode Control Theory and Background .................................. ...............69
W hat is Sliding M ode? ........................................................... .........................69
T h e o ry .......................................................................... 6 9
Sim ple E xam ple .................. ................... .............. ........ .. ...... .... 72
Other Sliding Mode Designs for Spacecraft .................................... ...............73
Satellite Controller Design For This Work ..................................... ............... 75
Controller D erivation ......................... ........... ........ ............ 75
Fuel Slosh D isturbance ............................................... .............................. 77
Stability of C ontroller .................... .. ........................ ......... ........... 78
Saturation ..................... .......................80
Com paring The Tw o Types of Controller ........................................ .....................81
C o n c lu sio n ................................................................. 8 1

6 SIMULATION AND RESULTS ............................................................................83

Simulation Dynamics and Parameters ............................................. ............... 83
S a te llite ......................................................................... 8 3
Satellite D ynam ics ....................... .................. ................... .. ...... 85
F u e l T a n k ............................................................................................................... 8 6
Reaction Wheels ................................. ..... ........ ..............86
S o la r P a n e ls ..................................................................................................... 8 7
C o n tro l S y ste m ............................................................................................... 8 7
F finding an O ptim al P D G ain ........................................................... .....................88
F u e l S lo sh ............................................................................................................. 9 3
IA E ...............................................................................9 4
IT A E ..............................................................................9 5
S ettlin g T im e ................................................................................................... 9 5
C M ................. ....... ... ...................................9 6
Torque and M om entum Lim its ....................................................... 96
Qualitative Analysis of Results .......... ... ...... ............. .....98
S o la r S n a p ..............................................................................9 8
Perform ance .............................................................98
Qualitative Analysis of Results ....... .............................99

7 C O N C L U SIO N ....................................................... 101

A ATTITUDE REPRESENTATIONS AND ROTATION MATRICES ....................104

F ix ed A ngle R rotation s ......................................................................................... 104
E uler A ngles ........................................ 105
Q u atern io n s ...............................................................10 7
Q uaternion A lgebra ................................................ ............... 108
Rotations of Rigid Bodies in Space. ........................... ..............110



vi









A attitude Errors in Q uaternions .............. ........................................... .............. 111
Sum m ary of Q uaternions ................................................. ....................... 112

B M O D A L E Q U A TIO N S .................................................................... ...................113

M option of M ulti-Degree of Freedom System ................................. .................... 113
Dam opening ................................. ......................... ..... ..... ........ 115
Sim ple E xam ple .............................................................................................. 116

C SIM U L A T IO N C O D E ............................................................................................... 117

Simulation Code .................. ..... ................... 117
M a tla b ................................................................1 1 7
C Programming Language .....................................................................117
C++ Programming Language ........... .......................................118
cM athlib H leader File ........... .................................... .......... .. ............ 119
cRK 4 H leader File ................................... .. .......... .. ............126

REFEREN CES .................................. .. ... ........ .............. 127

B IO G R A PH IC A L SK E T C H ............................................ ........................................ 134
































vii















LIST OF FIGURES


Figure age

1-1 Finite element analysis of a fuel tank inside of a satellite. ......................................6

1-2 PP T for E O -1 ............................................................................ 10

1-3 D iagram of a PP T ........... .................................. ... ................... 10

1-4 NEAR spacecraft orbiting the asteroid Eros. ............. ......................................... 11

1-5 Hubble Space Telescope diagram .............................................. ............... 13

2-1 Traditional fuel slosh pendulum model where M is the total fuel mass center, mp is
the pendulum mass, ms is the stationary mass fraction, L is the pendulum arm
length, and g is the local acceleration due to some external force. .....................22

2-2 Propagation of fuel slosh in an elliptical tank over time. These images demonstrate
the difficulty in modelling the fluid movement with simple representations. ......22

2-3 Diagram of a single sloshing mode, with its mass attached to the fuel tank walls by
springs and dampers aligned along the x, y, and z axes. .....................................23

2-4 Plot of the effect of multiple modes in the velocity of a simply supported beam. This
is representative of the fuel slosh disturbance in flight data. .............................24

3-1 The Soviet satellite Sputnik which orbited the Earth in October 1957. ................27

3-2 Diagram of the penumbra and umbra regions which block a satellite from sunlight.
....................................................................................................... 2 8

3-3 Another view showing the path of a satellite and how its orbit does effect the heat
flux it is exposed to. ..................... .................. ................... ........ 28

3-4 Cross section of solar panel with coordinate systems superimposed. Note that the
S's in the figures represent 's and M's are moments. ........................................30

3-5 Free body diagram of a section of the beam of size dx where M(x,t) is a moment,
V(x,t) is a shear force, f(x) is a force per unit length, and u(x,t) is the height of the
m ass elem ent at position x at tim e t. ........................................ ............... 32









3-6 Plot of [3-49]. Intersections of the two lines are solutions of the characteristic
equation ...........................................................................36

3-7 A cross section of the beam in which a heat flux is applied to one side. Notice the
beam starts off at a uniform temperature. .................................. ............... 42

3-8 Plot of the temperature change of a beam. Notice that the temperature difference
between the two sides remains constant eventually. ...........................................42

3-9 Hubble Space Telescope which shows the solar panels bent due to a constant
tem perature gradient. ..................... .. .... ................ ......................... 43

3-10 Another picture of the Hubble Space Telescope showing deflection of the solar
panels during maintenance by the space shuttle. ............ ................................... 44

4-1 The three frames of reference commonly used in spacecraft dynamics. ...............47

4-2 A reference frame in motion relative to an inertially fixed frame of reference. ....48

4-3 A rigid body satellite. ..................... .. .... ................ ......................... 53

4-4 A rigid body spacecraft with a single reaction wheel. ........................................57

4-5 A plot of common environmental disturbances a spacecraft is subjected too. ......60

5-1 Two types of reaction wheels produced by Ithaco. ............................................66

5-2 Torque rods produced by Ithaco. ........................................ ....................... 67

5-3 Spacecraft engine produced by Boeing. ..................................... ............... 68

5-4 Two possible solutions for the double integrator, where (a) is and (b) is These plots
are typically referred to as phase portraits. They plot position and velocity in this
example, but later (when we are doing controls with full state feed back) the
portraits will represent error and error velocity. ............ .................................... 70

5-5 Switching between the two controllers results in driving the system to zero in a stable
m an o r. .......................................................... ................ 7 1

5-6 Graphical representation of how effects the sliding surface's orientation in the phase
plane. The dotted lines represent the sliding surface (s = 0). .............................72

6-1 EO-1 satellite with solar panel deployed. ................................... ............... 83

6-2 Diagram of satellite's internal components. .................................. .................84









6-3 Composition of a typical simple solar panel. ................................. .................86

6-4 Cost function surface plot for the optimal gain of the PD controller. ...................89

6-5 ITAE results from the optimal PD search. ................................. .................90

6-6 IAE results from the optimal PD search. .................................... .................90

6-7 Settling time results from optimal PD search ....................................................... 91

6-8 Control momentum or fuel used searching for the optimal PD controller. ...........92

6-9 A depiction of the step maneuver showing the orientation of the satellite at start and
end. The angular distance and axis of rotation are also shown for better
understanding of the term s. .............................................................................. 94

6-10 Quaternion error and control effort of the PD controller rotating from -30 to 30.
...................................................................................................... 9 6

6-11 Quaternion error and control of the SM controller rotating from -30 to 30. .......97

6-12 Quaternion error and control of the SMCNF controller rotating from -30 to 30.
...................................................................................................... 9 7

6-13 Arcsecond error of the PD controller during solar snap. ............... ............ 99

6-14 Arcsecond error of the SM controller during solar snap. ...........................100

6-15 Arcsecond error of the SMCNF controller during solar snap. ........................100

A-i Body reference frame attached to a rigid body. The x-axis points out the front of the
v eh icle ........................................................................ 10 5

B-l Simple multi-degree of freedom system composed of mass and spring elements.
.................................................................................................... 1 1 3

C-1 UML diagram of Mathlib class which contains C++ objects for vectors, matricies,
and quaternion used in this work. .............................. .... ....................... 119















LIST OF TABLES


Table page

1-1 Spacecraft with thermally induced disturbances. ................................................15

3-1 Solar panel thermal material properties. ..................................... ............... 32

3-2 First four eigen modes of the beam ............................... ..................35

5-1 Comparison of common active attitude control hardware. ....................................65

5-2 Sum m ary of the switching law ........................................ ......................... 70

6-1 Technical Specifications for the EO-1 spacecraft ............................................84

6-2 EO -1 fuel slosh param eters. ............................................ ............................ 86

6-3 E O -1 w heel specifications. .......................................................... .....................86

6-4 EO-1 solar panel param eters. .............................................................................87

6-5 Initial and final euler angles where roll, pitch, and yaw are equal. ......................93

6-6 IAE results with values having units of radians ............................... ...............94

6-7 ITAE results with values having units of radian-seconds2. ............ .................95

6-8 Settling time results with values having units of seconds. ....................................96

6-9 CM results with values having units of Nms. ................................. ...............96

6-10 Maximum wheel momentum achieved where units are in Nms. .......................98

6-11 Solar snap results of the three controllers. ...................................................99

A-i Properties of a rotation m atrix. ............ ...... ............................................105

A-2 Comparison of the two major types of rotations, sequential and first and third axes
ro tatio n ........................................................................ 1 0 6










A -3 Quatem ion Algebra Sum m ary ........................................ ......................... 109





























































xii

















LIST OF SYMBOLS



Thermally Induced Vibrations Symbols

a Absorptivity

CT Coefficient of thermal expansion.

n Natural frequency of mode n.

P Density of a material.

K Thermal diffusivity.
A Cross sectional area of the solar panel.
c Heat capacity.
E Young's Modulus
F(x) Position function from the separation of variables method.
Fn Position function for mode n.
G(t) Time function from the separation of variables method.
h Thickness of the solar panel.
k Thermal conductivity.
MT(t) Thermal moment.
qo Total heat flux on surface.

qs Heat flux from sun.

qe Earth emitted heat flux.

qa Earth reflected heat flux.

u(x,t) Displacement of solar panel.
V(x,t) Shear force.




Satellite Dynamics and Fuel Slosh Symbols

OB/I Angular velocity of satellite system relative to inertial frame.

oB/O Angular velocity of satellite system relative to orbital frame.

0O/1 Angular velocity of orbital frame relative to inertial frame.











Satellite Dynamics and Fuel Slosh Symbols


C Rotation matrix.
c First moment of inertia.

Cfn Damper coefficient for sloshing mode n.
F External forces and torques.
h Angular momentum.
J Second moment of inertia.
kfn Spring coefficient for sloshing mode n.
L Length of pendulum.
M Mass matrix of satellite system.
mp Pendulum mass.
ms Stationary mass fraction.
p Linear momentum.
q Quaternion
The imaginary part of the quaternion, 1-3
q
Rc Distance of satellite from center of Earth.
T Kinetic energy of satellite system.
V Composite linear and angular velocity of satellite system.


Controls Symbols

x Weighting term on the sliding surface.

Kd Derivative controller gain.

Kp Proportional controller gain.
PD Proportional Derivative controller
s Sliding surface.
SM Sliding Mode controller.
SMCNF Sliding Mode controller with a colored noise filter
Td Gain ratio between the proportional and derivative gains.















Abstract of Dissertation Presented to the Graduate School
of the University of Florida in Partial Fulfillment of the
Requirements for the Degree Doctor of Philosophy

ROBUST NONLINEAR ATTITUDE CONTROLWITH DISTURBANCE
COMPENSATION

By

Kevin J Walchko

May 2003


Chair: Paul Mason
Cochair: John Schueller
Major Department: Mechanical and Aerospace Engineering

Attitude control of small spacecraft is a particularly important component for many

missions in the space program: Hubble Space Telescope for observing the cosmos, GPS

satellites for navigation, SeaWiFS for studying phytoplankton concentrations in the ocean,

etc. Typically designers use proportional derivative control because it is simple to

understand and implement. However this method lacks robustness in the presence of

disturbances and uncertainties. Thus to improve the fidelity of this simulation, two

disturbances were included, fuel slosh and solar snap.

Fuel slosh is the unwanted movement of fuel inside of a fuel tank. The fuel slosh

model used for the satellite represents each sloshing mode as a mass-spring-damper. The

mass represents the wave of fuel that propagates across the tank, the damper represents the

baffling that hinders the movement, and the spring represents the force imparted to the

spacecraft when the wave impacts the tank wall. This formulation makes the incorporation









of multiple modes of interest simple, which is an advance over the typical one sloshing

mode, pendulum model.

Thermally induced vibrations, or solar snap, occur as a satellite transitions from the

day-to-night or night-to-day side of a planet. During this transition, there is a sudden

change in the amount of heat flux to the solar panels and vibrations occur. Few authors

have looked at the effects of solar snap. The disturbance dynamics were based on the work

by Earl Thorten. The simulated effects compared favorably with real flight data taken

from satellites that have encountered solar snap.

A robust sliding mode controller was developed and compared to a more traditional

proportional derivative controller. The controllers were evaluated in the presence of fuel

slosh and solar snap. The optimized baseline proportional derivative controller used in this

work showed little effort was needed to obtain better performance using sliding mode. In

addition, a colored noise filter was developed to compensate for the fuel sloshing

disturbance and incorporated into the sliding mode controller for greater performance

increase at the expense of requiring a little more control effort.















CHAPTER 1
INTRODUCTION

This chapter provides the motivation and a literature review for this work. Additional

literature will be presented in subsequent chapters as appropriate.

Formation Flying

As scientists aim for higher goals of discovery, satellites prove to be an increasingly

important tool. Bigger, more powerful satellites can search the cosmos more easily for

distant planets, map the Earth, or even try to discover the origins of the universe itself.

However, when a satellite reaches a certain size it effectively can no longer be launched

into space.

The only answer to bigger, more powerful satellites is either construction in space or

formation flying. Construction is space is difficult and costly, but perhaps one day when

the international space station (or its successor) is finished this will be a viable solution.

Formation flying is the only near term solution. Thus, the idea of not one big satellite but

rather many little satellites, each acting as an individual element of a much larger sensor

array, is actively being pursued by National Aeronautical and Space Administration

(NASA) and other agencies. These missions greatly improve the sensing capabilities

available to NASA, while keeping the associated cost and risk down.

One of the main problems with formations and individual spacecraft in orbit is the

amount of fuel consumed. Small differences in altitude between satellites will result in

different velocities and orbits. This will require some satellites to expend more fuel in









order to continually maintain formation. In addition, if they are in different inclinations,

they are in different orbital planes (planes that cut through the center of the Earth) and thus

will naturally move farther apart or closer together depending on where they are in the

orbit. Therefore fuel is required for formation maintenance, which inevitably leaves fuel

tanks partially filled. Partially filled tanks can lead to fuel slosh which can degrade the

pointing accuracy and therefore affect the science mission. Fuel slosh compensation is not

an easy task and unfortunately is often neglected.

Satellites do not typically change their orbits unsupervised, but with formations there

will have to be some level of autonomy where the satellites must make corrections to their

orbits. Satellites will be required to change orbit to close gaps in formations or to avoid

collisions from orbital debris. Gaps in the formation will occur as satellites run out of fuel

or have failures. The control solution must be robust to account for uncertainty, produce

stable motion, be easy to understand, and easy to implement. There has been much

research on satellite formation control [1-8]; however these works assume satellites with

robust individual control systems.

General Problems Controlling Single Satellites

To make a formation more stable, we must improve individual satellite control. This

section will provide a brief overview of some of the problems with satellite control. Later

greater detail will be presented on two specific problems: fuel slosh and solar snap.

External Disturbances

NASA engineers must contend with many types of disturbances when controlling a

spacecraft. Low altitude satellites must fight gravitational and aerodynamic torque which

can deteriorate their attitude accuracy and over time degrade their orbit. These









gravitational torques can be present on large satellites where gravity has more of an effect

on one end of a satellite than the other.

In deep space or geosynchronous orbits, solar based disturbance becomes the

dominate environmental disturbance since there is no friction or dampening. Also in these

orbits, the gravity effects become negligible compared to effects produced by the sun. One

of the major solar disturbances is solar snap which is thermally induced vibrations.

Vibrations due to thermal changes have been known for a long time, but here on Earth

there are few situations where these vibrations are of any interest.

Another disturbance due to the sun is called radiation torque. Photons emitted from the

sun slam into the satellite (a force) which produces a torque about the center of mass of the

craft. This radiation torque grows and wanes with solar flare activity, but none the less

pushes and prods satellites until they slip behind the dark side of a planet or moon.

Parametric Uncertainty

Control systems are typically designed using mathematical models of the system to be

controlled. Thus there is a direct correlation between the performance of a standard

classical controller and the accuracy of the mathematical model of the plant.

Unfortunately in satellite control, the inertial matrix is only known within 5% of its in

orbit value. This model inaccuracy must be addressed in order to achieve the higher level

of performance expected from the missions. Thus engineers are currently examining new

control schemes which help to combat this problem and improve performance in the

presence of modelling error.

Another problem that effects the inertia matrix is fuel slosh and fuel expenditure. Fuel

can compose 40% of the total mass of a satellite prior to launch [9]. During orbit insertion,









some of this fuel will be burned to perform various Av maneuvers to achieve a desired

orbit. However, some of the fuel will be kept so that later orbital corrections can be made

during the life span of the spacecraft. This remaining fuel poses a problem as it moves

around inside of the satellite and thus changes the inertial matrix on orbit.

Literature Review

This section will provide a review of spacecraft that have had problems with fuel slosh

and solar snap. The controllers developed to compensate for these disturbances will also

be discussed. An analytical description of these disturbances will be provided in a later

chapter.

Fuel Slosh

Modelling fuel slosh dynamics

Accurate analysis of the fuel requires the use of the Naiver-Stokes equations, which

are extremely computationally intensive. Since most satellites have a limited amount of

processor power assigned to control, this is not practical for controller implementation.

Therefore the most used models are the simplest ones such as the pendulum, mass-spring-

dampener, and moving fuel tank walls [11]. For small motions, the objective is to

reproduce the sloshing mode resonance frequencies. Typically a mass is attached to the

wall of the fuel tank by a spring. The movement of this mass accounts for the free surface

wave motion and the influence on the spacecraft's dynamics is represented by the force

the spring exerts on the wall of the fuel tank. Dampeners are commonly introduced to

represent the viscous dampening of the fuel. There is typically one set of mass-spring-

dampeners for each sloshing mode that is being represented. Swirl effects of the fuel can

also be included as reaction wheels with torsional springs attached. These models however









are only accurate in representing the sloshing modes when small linear or angular motions

are performed. For accuracy during larger motions, better models must be utilized.

However these models are not ideal for control purposes. The typical pendulum model

used for fuel slosh is derived from work done with rockets. Typically in this scenario, the

rocket's main engine is firing and produces a local acceleration that pools the fuel at one

end of the tank. The amount of fuel that remains pinned to this end of the tank is the

stationary fuel mass fraction. The rest of the fuel is able to oscillate about the center line of

the fuel tank in a pendulum-like motion. These oscillations are typically small since the

fuel has only a small momentum compared to the forces produced by the local gravity.

Unfortunately this model is not as useful in the situation of satellite attitude control,

since there is no external force produced by a thruster that is constantly firing throughout

the life of the vehicle. Thus a model is needed in which the fuel has more freedom of

movement, can exhibit a larger range of movement, and produce accurate results.

The work developed by Agrawal [10] for use on the INTELSAT IV, dual spin

spacecraft, is a boundary layer model of the fuel for a spinning spacecraft. The author

reached the conclusion that the analytical model gives the most accurate prediction of the

fuel dynamics. The finite element model consisted of 81 nodes which represented a

spherical fuel tank, shown in Figure 1-1. This fuel tank had a radius of 0.42 m, was 1.31 m

from center of mass of the craft, and spun at 30 rpm. The solution for the fluid motion is

obtained by solving three equations: an inviscous fluid problem, a boundary layer

problem, and a viscous correction problem. The inviscid and viscous correction solutions

are obtained by using finite element methods, while the boundary layer problem is solved

analytically. The numerical results showed two slosh modes which could be accurately















-- ----J


Il.: .




\ L




Figure 1-1. Finite element analysis of a fuel tank inside of a satellite.

modeled by the pendulum. However there were also first azimuth and elevation modes,

lower modes to be inertial modes and higher modes to be a combination of inertial and

slosh modes. These other modes would not be accurately modeled by a pendulum.

Controlling fuel slosh

Fuel slosh is an example of under-actuated control. The objective is to control both the

rigid body dynamics of the spacecraft and the fluid dynamics of the fuel only having

access to effectors that act directly on the spacecraft. The fuel's dynamics are controlled

through the coupling of the two. One of the main problems with controlling this type of

disturbance is that one can not measure the position, orientation, etc. of fuel in a satellite.

One can only measure the effects the fuel slosh on the total system. The states and

parameters of the fuel must be estimated, and a controller must try to account for the fuel's

influence. Therefore it is difficult to control what can not be measured directly. Thus

several types of passive methods are used to control slosh such as baffles [12], slosh









absorbers [13] to dissipate energy through sloshing, and large fuel tanks are broken up into

smaller ones. However these steps add complexity and weight to a satellite's design. The

added complexity increases the construction time, the cost, and the chances of mechanical

failure of the satellite. The added weight also increases the cost of launching the satellite

since it costs approximately $10,000 per pound to put things into space.

Fuel slosh directly affects the performance of the attitude control system, in both

pointing and tracking maneuvers. One of the simplest ways to handle fuel slosh is to

identify the sloshing frequencies and eliminate the influence on the feedback with a notch

filter. Unfortunately this denies the satellite the ability to operate in certain frequency

ranges since the filter attenuates the slosh frequencies in the feedback. Hence, normal

feedback without sloshing effects will also be attenuated. Also as fuel is expended and the

sloshing frequency changes, the filter would need to be updated to notch out the new

sloshing frequencies again.

Cho and McClaroch [14] developed a controller which utilized a common 2D

pendulum model. Their equations only dealt with motion in a plane, but "a significant

extension of what has been done previously, since it simultaneously controls both the rigid

vehicle motion and the fuel slosh dynamics." One of the major problems they set out to

solve, which others were unsuccessful in solving, was the transverse motion of the craft.

Their method utilized a Lyapunov function approach, and was able to control the system.

They utilized additional control effectors to counter this movement (i.e., gas jet thrusters)

while also controlling the pitch and suppressing the slosh dynamics.

Kuang and Leung [15] use a simple model of a liquid-filled spacecraft subject to

external disturbances, and control it with a state feedback H, controller. The spacecraft's









elliptical fuel tank is assumed to be completely filled with an ideal liquid in uniform

vortex motion. Because of these assumptions, the slug of liquid can be represented by a

finite number of variables. These variables can then be solved using Helmholtz equations.

For the simulations, the liquid slug is located at the center of mass of the vehicle and

surrounded by a viscous layer. An H. attitude controller was then developed using the

theory of Hardy-infinity optimal control theory. The controller was capable of attenuating

the vortices in the fuel tank. The pointing accuracy remained constant regardless of small

changes in the viscosity. The authors plan to next look at partially filled fuel tanks, which

is a more difficult problem.

Thus far we have discussed more complex models to model liquids and said that in

order to provide the level of control that a designer wants, one must use a complex model.

An alternate example is Grundelius and Bernhardsson [16], who developed a minimum-

time optimal controller for moving open rectangular fluid containers down an assembly

line. Proper control of the containers is important. If the accelerations applied to the

containers are too high, then the contents of the containers may spill out of the open top.

This is obviously a loss to the company in material and results in a more expensive

process. Also, if the top of the container is to be eventually sealed, the contents may spill

on to an area where glue may be applied, thus making for a bad seal. This could result in

customers disliking the product because of "cheap" design. When the speed of the system

began to exceed two meters per second, fluid began to slosh out of the containers and the

controller was not able to achieve the desired result. The authors incorporated a minimum-

energy cost function into the design. Using only a simple linear model for the position of

the surface of the liquid and the container, a controller that used a quadratic penalty on the









slosh and a quadratic penalty on the terminal state was able to produce the desired results

without a complex hydro-dynamical model.

Do we need to worry about fuel slosh?

Fuel slosh is associated with liquid fuel which is used by thrusters. These types of

thrusters are efficient and can provide a significant amount of control authority. Slosh

forces can be eliminated via the use of alternate actuators. Let us take a look at some of the

various types of control actuators available to satellites.

Wisniewski and Blanke [17] developed a three axis magnotorquer controller for

satellites subject to a gravity gradient while in a polar orbit about the Earth.

Magnotorquers are an attractive form of control for small inexpensive satellites since they

use the Earth's magnetic field to control the orientation of the spacecraft. However

magnotorquers produce very small torques, and their principle of operation is nonlinear

and difficult to use. They are difficult due to the fact that control torques can only be

generated perpendicular to the Earth's magnetic field. Also the size of the spacecraft might

prohibit the use of magnotorquers since they produce low torque.

Chen et al. [18] developed an optimal controller which utilizes both thrusters and

magnotorquers to unload extra momentum that has accumulated and saturated the reaction

wheels of a spacecraft. Their method separates the required torques for the thrusters and

magnotorquers using a simple optimal algorithm. Simulations conducted by Chen showed

a 20% fuel savings using this method.

Liquid [19] apogee motors are common in geosynchronous spacecraft design. This

results in liquid fuel constituting almost half of the vehicles' mass. In space, liquid fuel

motion greatly influences the attitude stability and control. The level and nature of these









influences are dependent on the geometry and location of the fuel tanks, the ratio of

spacecraft mass to fuel mass, the fill level of the fuel tanks, and the physical properties of

the fuel itself. They also depend on the flexible structure of the spacecraft and the behavior

of the estimation and control schemes that are present. Unfortunately, accurate prediction

of the fuel dynamics is a difficult problem due to the complexity of the hydrodynamical

equations of motion. LEASAT, which was launched in 1984, experienced attitude

instability during the pre-apogee injection phase which immediately followed the

activation of despin control. The instability was the result of the interaction between the

lateral sloshing modes and the attitude control.

Pulsed plasma thrusters (PPT) are becoming a more popular control solution for small

spacecraft [20]. This fact is due to the small electrical consumption and the use of a solid

teflon fuel. The fuel is attractive because it is stable, compact, and does not provide

complications due to sloshing liquids. EO-1 provided a successful demonstration of an

attitude control system using these PPT's shown in Figure 1-2 and Figure 1-3. PPT's are

capable of very high pulse width modulation frequencies (PPT's are only capable of on-

off or bang-bang type of control). Although PPT's have some very favorable properties,

they still have the liability of using a non-renewable fuel and are only capable of operating







F





Figure 1-2. PPT for EO-1. Figure 1-3. Diagram of a PPT.





























Figure 1-4. NEAR spacecraft orbiting the asteroid Eros.


on small spacecraft.

Another technology that shows promise for satellites ranging from 100 grams to 500

kg is MEMS engines [21]. MEMS (micro-electrical-mechanical systems) is an exciting

new technology field which currently is beginning to make its way out of the R&D labs

and into commercial products. MEMS are microscopic mechanical-electrical systems that

are developed using a combination of micromachining and standard integrated circuit

design. Some of the common MEMS components on the market are accelerometers,

gyros, and GPS receivers. Here MEMS rocket engines are developed for space

applications. The MEMS engines are the size of a human hair and have a thrust-to-weight

ratio hundreds of times greater than the space shuttle.

Although there are many good alternatives to liquid thrusters for small spacecraft,

liquid fueled thrusters still provide the best source of high thrust capability for both small









and large spacecraft. Therefore fuel slosh will always be present in these types of systems

and thus affect high precession pointing and tracking accuracy.

Missions involving fuel slosh disturbances

We now understand why fuel slosh exists. Next let us look at missions where slosh has

been observed. The Near Earth Asteroid Rendezvous (NEAR) spacecraft was sent to orbit

and observe the Eros asteroid on December 23, 1998. The spacecraft, shown in Figure 1-4,

had to fly by the asteroid following an unsuccessful firing of its main engine a few days

earlier. A subsequent successful firing of the engine put NEAR on course to rendezvous

with Eros to begin its planned year long orbital mission starting in mid-February 2000.

The reason the main engine firing did not work is the on board control system terminated

the bum when sensors registered higher than expected lateral accelerations. These

accelerations were later determined to be caused by the spacecraft's fuel. The NEAR's

mission had to be delayed for over a year because the control system was incapable of

handling the coupled dynamics.

Situations like NEAR do not just happen around asteroids, but with the advent of the

International Space Station similar problems could arise. As resupply ships travel to and

from the station, similar problems could occur and a ship could crash into the station. This

scenario is similar to the accident that occurred with MIR and one of its resupply vessels.

Although that incident was not caused by fuel slosh, it still does represent a potential

problem.

Thermally Induced Vibrations: Solar Snap

Thermally induced vibration or solar snap occurs when dynamics are excited by a

rapid temperature change. This type of vibration has been known for some time but rarely






























Figure 1-5. Hubble Space Telescope diagram.


appears in applications on Earth. However in the harsh, irradiated vacuum of space,

without the presence of an atmosphere to dampen the vibrations, this disturbance can

adversely affect the performance of high precision instruments and equipment. This

section, like the proceeding one, will present missions that where adversely affected by

solar snap and control systems that were designed to compensate for it.

Controlling thermal vibrations

The Hubble Space Telescope (HST), shown in Figure 1-5, was launched on 24 April

1990 with the hopes of making numerous new scientific discoveries. The HST was a 13

ton, free-flying spacecraft with a pointing precision requirement of 0.007 arc-sec over a 24

hour period (which is the most stringent requirement ever imposed). During its initial

checkout period, a pointing "jitter" was discover. The jitter was produced by a thermally

induced vibrations imitating from the solar panel array. Later in December 1993 the solar









arrays were replaced by the Space Shuttle Endeavour in an attempt to reduce the jitter. The

solar panels exhibited an end-to-end bending oscillation when the spacecraft transitioned

between sunlight and shadow. They also had a sideways oscillation on the day side of the

Earth. At its worst, the 20 ft. solar arrays would deflect as much as 3 ft. This would have

the effect of rendering any long time exposures of over 25 min. pointless. The original

controller was a digital PID controller with an FIR filter in the rate path to attenuate high-

frequency, main-body bending modes. This design was incapable of handling the induced

disturbances from the solar arrays, which caused a peak pointing jitter of 0.1 arc-sec.

([22], p. 583)

Wie and Liu [25] redesigned the attitude control system using an H. controller. This

controller was designed to attenuate the two major modes from the solar panels which

operated at 0.12 and 0.66 Hz. Through a trial-and-error process, they were able to tune the

weights of the controller and properly achieve the mission requirements.

Singhose et al. [23], while working at M.I.T., designed a method to minimize

vibrations in flexible spacecraft which utilized command profiles. The commands are

designed by the following:

1. Formulate the equations of motion for a system subject to a sequence of impulses.

2. The equations are solved for the sequence of impulses that results in the smallest
amount of residual vibrations in the system.

3. The sequence is convolved with the desired input to generate the command sequence
which is issued to the real system.

The command sequence is a pulse width modulation (i.e., on-off sequence) of commands

for the spacecraft's thrusters. Later at Georgia Tech., Ooten and Singhose [24] extended

the work to sliding-mode control. The drawback for this work is, although flexible space-










craft are mentioned, all simulations involve only a mass-spring-damper system. Also,

solving the equations of motion (step 2) for a nonlinear, flexible satellite would be very

difficult.

Missions which involved solar snap

Originally spacecraft were very simple, but as technology progressed into the 1960s

more appendages were being put on. These structures offered new challenges in

controlling flexible structures and vibrations. Below is a table of spacecraft that had their

performance compromised by thermal vibrations. These satellites do not include military

and Soviet spacecraft.

Table 1-1: Spacecraft with thermally induced disturbances.

Spacecraft Dale Description
Ulysses 1990- Spin-stabilized solar probe experiences unexpected nutation due to
1991 thermally induced vibrations of beryllium-copper axial boom
Upper Atmosphere 1991 "Thermal snap" at orbital night/day transition attributed to rapid
Research Satellite heating of large single solar array
Hubble Space Tele- 1990 Jitter phenomena attributed to thermally induced vibrations of
scope FRUSA-type solar arrays
LANDSAT-4/5 1980s Thermal elastic shock experienced during orbital night/day transition
due to disturbance of large single solar array.
Communications Tech- 1978 Three axis stabilized satellite experienced "thermal elastic shock"
nology Satellite during orbital night/day transition
Voyager 1977 Low frequency oscillations of PRA booms during LEO operations
attributed to "thermal flutter"
Apollo 15 1971 Thermally induced vibrations of twin Bi-STEM booms filmed during
64th lunar orbit by astronauts.
Explorer 45 1971 Unexpected nutation of spin stabilized satellite attributed to thermal
bending of four radial booms.
NRL 161, 163, and 164 1969 Thermally induced vibrations of STEM booms used for gravity-gra-
dient stabilization leads to large satellite motions.
OGO-IV 1967 NASA Orbiting Geophysical Observatory Satellite STEM booms
experience strong thermally induced vibrations during orbital night/
day transition.
OV-10 1966 US Air Force gravity-gradient satellite using STEM booms equipped
with tip masses experiences thermal flutter and complete inversion in
orientation.









Table 1-1: Continued
Spacecraft Date Description
Explorer XX 1964 NASA spin-stabilized satellite experiences spin decay due to interac-
tion between STEM boom thermal bending and solar radiation pres-
sure.
1964-83D and 1963- 1964/ Dynamic thermal bending of STEM booms on APL gravity-gradient
22A 1963 satellites.
Aloutte 1 1962 Canadian spin-stabilized satellite experiences spin-rate decay due to
thermal structural response of STEM booms similar to Explorer XX.
a. Information contained in this table comes from Thorten [26] on page 345.

"In spite of the practical importance of thermally induced vibrations to the successful

operation of spacecraft, recurrent vibration problems from the 1960s to the present day

suggest that the phenomena is either underemphasized or not well understood in

spacecraft design" (Thorten [26], p. 344)

Parametric Uncertainty

A very simple way to model the effects of fuel and thermal vibrations is to represent

them as external disturbances (i.e., deterministic noise). Depending on the mass fuel ratio,

fuel slosh affects the mass distribution of the vehicle which affects the inertia matrix.

Since control systems are model based, this will directly affect the level of control. If a

control system could be designed that provided generic control, or model independent

control, then the effects of these disturbances could be negated.

This idea of generic control was the subject of work done by Walchko and Mason [27]

utilizing fuzzy logic. The idea was that fuzzy logic is rule based and not model based

control. Therefore fuzzy logic would be capable of providing better control in the presence

of uncertainty. Although several different configurations of fuzzy logic were developed

(PD, PID, and sliding mode) none produced the desired results when there was a

significant amount of uncertainty. In fact it was shown that the two best were sliding mode









and fuzzy sliding mode. However the non-fuzzy sliding mode had the added advantage

that stability could be shown for the controller, whereas fuzzy logic lacks any general

method to prove stability

Petroff et al. [28] attempted to address this problem of fuzzy logic stability. They

looked at simple linear single input single output systems and were able to develop a

stability analysis. However, this method did not naturally extend to multi-input multi-

output systems and definitely not to nonlinear systems. Thus although fuzzy logic

appeared to be attractive as a model independent control architecture, it was extremely

weak with respect to stability analysis. This fatal flaw ultimately disqualified it for our

applications.

Summary and Overview of the Thesis

What is needed are controllers that are easy to adapt from one satellite to another but

still retain their robustness and efficiency. Problems due to simple control systems being

implemented on satellites are numerous. This usually results in a redesign of the control

system which is costly and a waste of time and resources. Now couple these problems

with the difficulties of formation flying and things get very risky.

This work will primarily focus on fuel slosh and solar snap. Both of these phenomena

are interesting in that they have rarely been observed. Rather, engineers have examined

the data from satellites that have experienced unexplained problems and deduced these

were the causes. Neither of these problems is generally looked at when a control system is

developed, but rather they are dealt with only when they appear. This attitude proves to be

a problem when applied to formation control. The possibility of one satellite going









unstable due to these phenomena then introduces the possibility of that satellite affecting

the entire formation, thus endangering the mission.

Chapter 2 will provide a more in depth discussion of fuel slosh. Various models will be

described that are used to model fuel slosh for simulation purposes and models that are

used in controllers.

Chapter 3 will cover the interactions of solar panels with satellite dynamics. Thermally

induced vibrations are described by the solution to their partial differential equations and

modal equations.

Chapter 4 will cover satellite attitude dynamics and control. This chapter will

introduce the reader to the fundamental dynamics of satellites.

Chapter 5 will introduce the basics of sliding mode control and other attitude control

algorithms and control hardware. The specific controller for this work will also be derived

and discussed.

Chapter 6 will provide an overview of the simulations. All parameters, equations, and

assumptions will be presented here. Also all results will be discussed in this chapter.

Chapter 8 will provide conclusion to the work done.















CHAPTER 2
CONTROLLING THE UNKNOWN: FUEL SLOSH

The introduction covered much of the problems associated with fuel slosh and

highlighted missions where this phenomena occurred. This chapter will expand on the

introduction by focusing more on the fuel slosh disturbance. Two models were developed,

one for use in the simulation to provide the disturbance torques on the satellite and another

model to help the controller account for the fuel slosh disturbances.

Introduction

The sloshing behavior of liquids and its modeling has been studied in many different

areas: fuel slosh in aerospace applications, movement of fluid filled containers in

industrial/manufacturing applications, earth quakes, vehicle and ship dynamics.

Unfortunately the choice of a model is not a simple one. Most analytic models that try to

accurately represent the dynamics are three dimensional partial differential equations.

These models are heavily dependant on boundary conditions and computationally

expensive to solve. Thus they are typically not used in controller design. However, the

nonlinear dynamics that appear from these complex models are important when designers

wish to move large amounts of fluid rapidly.

The aerospace industry has looked at controlling the effects of fuel movement in

aircraft fuel tanks for years. Most of the ideas (i.e., baffled fuel tanks, breaking one big

tank up into many smaller ones, etc.) can be traced back to standard aircraft design.









The interaction of the fuel and the spacecraft is more complex in three-axis stabilized

systems than spin-stabilized systems [29]. This is due to the fact that this system has

constraints that reduce the number of degrees of freedom, simplify the model structure,

and is numerically easier to compute. However, spacecraft will typically operate in several

modes of operation over their life span and thus it is important to be able to control fuel

slosh in any mode of operation.

Modelling Fuel Slosh Dynamics

Models for fuel slosh have not changed much since Abramson started in 1961 [9].

Models over the years have been: (single and multi) mass-spring-damper, pendulum,

liquid slug, and CFD/FEA models. Typically, the controller uses a model to account for

the fuel present inside of the spacecraft. Unfortunately, the models listed are too simple or

based on a static representation of the fuel and do not attempt to modify the parameters of

the model during run time.

The coordinate system for the fuel is typically the body coordinate system [29]. This

has the result of making the spacecraft's motions appear as a body-force field in the

Naiver-Stokes equations. The boundary conditions at the wall also become stationary with

only the free surface remaining transient while the viscous forces provide a dampening of

the liquid. The fuel motion is influenced by four classes of forces: gravity, inertia, viscous,

and capillary forces.

The early work with slosh dynamics and spacecraft centered around liquid filled

rockets and their stability [29]. This topic is still of interest and the assumptions [9] for this

problem typically are:

SThe fuel has small displacements, velocities, and slopes of the liquid-free surface.









* The fuel tank is rigid.

* The fuel is assumed to be a nonviscous liquid.

* The fuel is also modelled as an incompressible, homogenous fluid.

Propellant slosh is the induced motion of liquid due to acceleration of the container.

This motion along with its reactive forces can deteriorate the pointing performance. This

is especially true for systems with large fuel/weight ratio. Due to the nonlinear dynamics

of sloshing, it is not always possible to compensate for its affects using a standard attitude

controller.

Typically a pendulum system, shown in Figure 2-1, is used to model a single sloshing

mode. Since slosh dynamics typically display several modes, several pendulums should be

included in the model of the system. Thus the more dominate modes are modelled, the

more realistic the simulation.

A problem with this type of model is determining the values of the pendulum lengths,

masses, springs, and dampeners. One way to do this is to look at actual flight data where

you can determine the modal frequencies. Assuming the amount of mass in the slosh

mode, the value of the spring constant can be solved. The damper value is picked so that

the system's oscillations die out after an appropriate amount of time.

The main problem with liquid slosh dynamics is estimating the hydrodynamic

pressure distribution, forces, and moments. One reason that this is so difficult is the

dynamic boundary conditions at the free surface varies with time in a manor not known a

priori, as depicted in Figure 2-2. Hydrodynamic pressure in many rigid tanks is composed

of two parts. The first part is the fluid that moves with the tank in unison. The second

component is the sloshing at the free surface. This component is typically modelled as a












L









Figure 2-1. Traditional fuel slosh pendulum model where
M is the total fuel mass center, mp is the pendulum mass,
ms is the stationary mass fraction, L is the pendulum arm
length, and g is the local acceleration due to some exter-
nal force.


Figure 2-2. Propagation of fuel slosh
images demonstrate the difficulty in
simple representations.


in an elliptical tank over time. These
modelling the fluid movement with


mass-spring-damper or pendulum. Suggested values for the pendulum shown in Figure 2-

1 are given by EL-Sayad et al. [30] as


L = coth 1.84h) mp
1.84 \ RI)


[2-1]


where L is the length of the pendulum, R is the radius of the tank, and mp is the mass at the


end of the pendulum.


~k~E~q~,


totally tanh( 1.84









Fuel Slosh Models

Two different models for the fuel will be utilized in this work. The first model will be

included in the satellite dynamics and provide the proper disturbance torques on the

system. The second model will be computationally simpler, but still contain the important

aspects of the fuel's disturbance on the system.

Simulation Model

The model utilized is based on the mass-spring-damper models by Sidi [9] and Hughes

[31] shown in Figure 2-3. Each sloshing mode will be composed of a mass, a spring, and a

damper. The mass represents the mass of fuel in a specific mode. The spring represents the

force exerted on the tank wall by a sloshing mode. The damper represents the baffling in

the fuel tank that dissipates the sloshing movement.

The actual dynamics of this disturbance will be derived later when the dynamics of the

satellite are covered. This is due to the fact that both the dynamics of the spacecraft and

fuel are tightly coupled and it is easier to present them all at once.


















Figure 2-3. Diagram of a single sloshing mode, with its mass attached to the
fuel tank walls by springs and dampers aligned along the x, y, and z axes.

















S107
Mi 0de 1





6-


10 100 o
FrPieqcy (Hz)

Figure 2-4. Plot of the effect of multiple modes in the velocity of a sim-
ply supported beam. This is representative of the fuel slosh disturbance
in flight data.


Controller Model

The model used in the controller needs to be capable of representing several sloshing

modes, but can not be computationally expensive. Figure 2-4 provides a frequency

domain visualization of the proposed model. The model represents the sum of the slosh

modes as a band pass filter. This will allow us to incorporate the effects of the multiple

sloshing modes without having to specifically identify their individual characteristics.

A simple band pass filter is the combination of a low pass and high pass filter in series.

The standard transfer function for this filter is as follows:

as 1 as
[2-2]
as+ lbs+ 1 abs2 + (a+b)s+ 1

1 1
a = b = [2-3]
freqhighpass freqlowpass








where the transfer function on the left is a high-pass filter and the one on the right is a low-

pass filter. The variables a and b are the inverse of the high and low pass cutoff frequen-

cies. Putting [2-2] into state space form yields


ab ab + u [2-4]
x 1 0 -x


y = o0 x [2-5]


The interesting thing to note about 2-5 is the output given by the model is velocity and

not position. This is due to the derivative term (i.e., s in the numerator) in 2-3. This is a

reasonable characteristic, since the fuel slosh effects show up in the rate terms and not the

position terms.















CHAPTER 3
DYNAMICS OF SOLAR ARRAYS

As discussed in the introduction, thermally induced vibrations are a source of

problems for satellite attitude control. This chapter will provide the theoretical

development of thermally induce vibrations in long, thin beams. These results will finally

be used to generate the solar snap disturbance in the simulations.

Introduction

Research on thermally induced vibrations began back in the mid-1950's, before

Sputnik orbited the Earth in October 1957. However it was not until the 1960's that

applications for this new idea appeared. Initially spacecraft appendages (i.e., solar panels,

communication arrays, booms, etc.) where thought to be simple systems like those used

by Sputnik (Figure 3-1), but slowly it became apparent that they were much more

complicated. Strange and unexplained spacecraft behavior was now being attributed to

this new nonlinear behavior being observed.

Thermally Induced Vibrations

First, this section will derive and solve the partial differential equations associated

with solar snap. Then differential equations will be obtained from the PDEs which are

utilized in the simulations.

The derivation that follows will parallel work done by Thorten [26] and Rietz [67].

However, there was great difficulty in actually following their work. Fortunately, all key

equations were independently derived. The results in this work were also similar to

























Figure 3-1. The Soviet satellite Sputnik which orbited the Earth in Octo-
ber 1957.

derivations performed on cantilever beams excited by sources other than thermal

vibrations by Greenburg [32], Kreyszig [33], Meirovitch [34], and Zill and Cullen [35].

Heat Flux

Thermal vibrations are excited by a sudden change in heat flux to the object. In this

work, heat flux changes greatly during the day-to-night and night-to-day transitions. Here

the satellite passes into and out of the Earth's shadow.

The Earth's shadow is actually composed of two regions shown in Figure 3-2. The

penumbra is a region of partial shadow and the umbra is a region of full shadow. When a

satellite is in low Earth orbit (LEO), the penumbra is small and the satellite essentially

transitions from the day side to the night side (i.e., umbra). This sudden change from day-

to-night or night-to-day results in larger thermal vibrations. This is due to the sudden

change in heat flux, much like turning on or off a light switch.

The amount of heat flux (q) seen by a satellite is composed of solar flux from the sun

(qs), Earth emitted radiation flux1 (q,), and Earth's reflected heat flux (qa) from the sun.









solar vector
S1 1 1 1 1



rr






:i I
Figure 3-2. Diagram of the penumbra and Figure 3-3. Another view showing the
umbra regions which block a satellite path of a satellite and how its orbit does
from sunlight, effect the heat flux it is exposed to.

This amount changes with the satellite's position in orbit and its orientation relative to the

solar flux vector (see Figure 3-3).

q = qs+qe +a [3-1]

q, = 1350acos(p) [3-2]

where as is the surface absorptivity for solar radiation and Vp is the angle between the sur-

face normal and the solar flux vector.

Summing Forces

The solar panel is modelled as an Euler-Bernoulli cantilever beam subjected to a

thermal moment. The beam, a cross section is shown in Figure 3-4, will have the

following assumptions:

* Plane cross sections before bending will remain planes after bending.


1. The radiation referred to here is black body radiation.









* Lateral displacements in v and w are due to bending only. Deformations from shearing
and changes in transverse beam dimensions are negligibly small.

* The beam is a Euler-Bernoulli cantilever beam.

* The system is uncoupled, meaning that changes in beam orientation do not effect the
amount of heat flux entering the system.

* The temperature gradient is only one dimensional (i.e., through the thickness of the
beam).

The moments on the system are defined as


Mz = foydA [3-3]
A


MT = fEaTydA [3-4]
A

Linear strain displacement is as follows:

2 2
Ou duo v W [3
x -x dx ax 2 2[3-5


Hooke's uniaxial law is

ox = EEx-EaAT [3-6]

Then substituting [3-5] into [3-6] yields:

2 2
du o Ov 8 w
x = E y-2- EaAT [3-7]
x = dx 9x2 x2

Now substituting [3-7] into [3-3] produces

S2 2 2 2
du o v O wL V E W
S = -y z EaTydA = -El EM [3-8]
z dx Tx2 yx2 9x2 yz9x2 T


where















Sxz




Mz z x



Mx

Figure 3-4. Cross section of solar panel with coordinate systems superim-
posed. Note that the S's in the figures represent o's and M's are moments.


I= fy2dA and Iy = fy2dA [3-9]
A A

Next, summing the forces in they direction in Figure 3-5 results in


Fy = -f(x)- V(x, t) + V(x, t) + V(x, t)dx + mu(x, t) = 0 [3-10]


Note that from here on out, the terms in parenthesis will be dropped unless it is felt

that they are necessary for clarity.

2
8 u aV
pA dx + dx = fdx [3-11]
Oat ax

2
a u OV
pA + =f [3-12]
at2 Ox

Summing the moments about the center of mass of the beam yields:


M= M+ dx+ Vdx+- V+ dx dx -M = 0 [3-13]
Sx 2 2\ x









M+V = 0 [3-14]
ax


V [3-15]
ax

Substituting [3-8] into [3-5] and ignoring the middle term (since we are only looking

at motion in the u direction and not the w direction).


al__E- + M) = V [3-16]
ax\ ax2 T) =

Now taking [3-12] and substituting [3-16] gives us:

2 2
d ( ( a u f
pA +-- EI + = [3-17]
dt2 Ox[Ox Ox2 T))


2 4 2
du u MTr
pA- +EIf + = [3-18]
dt2 aX4 X2

This equation can be further simplified by assuming that no uniform load (f(x)) on the

solar panel exists and realizing the thermal moment is not a function of x, thus its second

derivative is zero.

f= 0 MT = 0 34 = E[319]
pA

Incorporating these, [3-18] becomes

utt + (34Uxxxx = 0 [3-20]

where a subscript ofx or t refers to differentiation with respect to that variable. The bound-

ary conditions of the system are as follows:

u(x, t) = 0 [3-21]


u (O, t) = 0


[3-22]










f(x)

i iI


Figure 3-5. Free body diagram of a section of the beam of
size dx where M(x,t) is a moment, V(x,t) is a shear force,
f(x) is a force per unit length, and u(x,t) is the height of the
mass element at position x at time t.


Mz(L) = -Eluxx-MT = 0

V(L, t) = Eluxxx + M = 0


Table 3-1: Solar panel thermal material properties.
Variable Description Units
k thermal conductivity W
mK
q heat flux W
m2K

p density kg
m3
c heat capacity j
kgK
K thermal diffusivity m2
s
a coefficient of thermal expansion 1
K


[3-23]

[3-24]









Nonhomogeneous Boundary Conditions

Unfortunately the boundary conditions for the system are nonhomogeneous and in

order to solve the equations this must be corrected. Getting rid of the nonhomogeneous

boundary conditions will be accomplished by changing the dependent variable.

u(x, t) = v(x, t) + p(x) [3-25]

Substituting [3-25] into [3-23] yields:

Mr
Ux = + xx [3-26]
+x xx El

The second derivative of u must be equal to zero to make the system homogenous, thus p

must negate the right hand side of the equation.


Mxx [3-27]
Pxx

Now we must solve for Vp(x), thus integrating [3-27] twice results in an equation with

two constants that can be solved by looking at the initial conditions.


= -MT2 + C X + C2 [3-28]
2EI

p(0) = = 0= C2 x(0) = 0 = C1 [3-29]

MT
(x) = 2 [3-30]
V (x) = 2Ex

Thus [3-25] can be rewritten as

MTX2
u(x, t) = v(x, t) [3-31]
2El

Our new beam equation and boundary conditions using [3-31] are shown below. Notice

how this is now a simpler homogeneous equation to solve.









vtt + 34vxxxx = 0 [3-32]

v(0, t) = 0 [3-33]

Vx(O, t) = 0 [3-34]

vxx(L, t) = 0 [3-35]

xxx(L, t) = 0 [3-36]

Separation of Variables

Now [3-32] with its homogeneous boundary conditions can be solved via separation of

variables. We will assume that the function v(x, t) can be broken up into two separate

functions. Each of these two functions are dependent on only one variable, either position

(x) or time (t).

v(x, t) = F(x)G(t) [3-37]

Substituting [3-37] into [3-32] gives us an equation which can be separated.

FxxxxG + 4FGt = 0 [3-38]

G Fxx
S xxxx = (02 [3-39]
G [4F


p4 [3-40]
pA

Since each side of the equation is dependent on a single variable, it must be equal to a

constant. We will choose this constant to be co2, which is squared to simplify the

derivation. Each of the two sides can now be solved, which gives us one differential

equation in x and one differential equation in t.

Gtt + 2G = 0 [3-41]








Fxxxx- (024F = 0 [3-42]

Solving the position function F(x)

We will begin by solving the F(x) equation first. Assuming that the solution to this

differential equation is an exponential, the following solution is found.

F(x) = B cosh(pnx) + B2sinh(p 3x) + B3cos(nx) + B4sin(3nx) [3-43]

F(0) = B +B3 = 0 [3-44]

Fx(O) = B2 +B 4 0 [3-45]

Fxx(L) = B(cosh(P3nL) + cos(P3nL)) +B2(sinh (P3L) + sin([3nL)) = 0 [3-46]

Fxxx(L) = B(sinh( 3L) sin([3nL)) +B2(cosh(P3nL) + cos(P3nL)) = 0 [3-47]

Combining the last two equations into a matrix give us:


cosh(3PL) + cos(3nL) sinh(PnL) + sin(InL) 0 [3-48]
sinh(PnL) -sin(PnL) cosh (3L) + cos (iL) B

The determinate of this matrix gives us the characteristic equation from which we can

solve for the values of P3,L that result in zero. The characteristic equation is shown below

and a plot (shown in Figure 3-6) which can be used to graphically solve the equation.

cosh (PL) cos (PL) + 1 = 0 [3-49]



Table 3-2: First four eigen modes of the beam
n 1 2 3 4
[3nL 1.875104 4.694091 7.854757 10.995541

It is now possible to solve the equations in the matrix as a ratio B2/B1. Thus it is

possible to solve for 3 of the constants in terms of the fourth.










F,(x) = cosh(3px)- cos(3px)+ ,(sinh(Px)- sin(3,x)) [3-50]


(P3,L)2 4 pA r 2 B2 sin(3L)- sinh(3,L)
SK n a [3-51]
n L2 p E" El n B1 cos(PL)+ cosh(P3,L)


This represents the n'' spatial solution or mode shape. In order to calculate the

complete response of the beam, a summation over all modes must be done.

00
F(x) = aFn(x) [3-52]
n= 1


The coefficient an is found by the orthogonality property of the modes. This coefficient is

what defines the amplitude of the vibrations while F,(x) only defines the shape of the

vibrating modes. Now because it is assumed that the eigenvectors form an orthonormal

basis, the dot product of any one eigenvector with another is zero. For example, if we take

the above equation and dot both sides with the first eigenvector only the first eigenvector






0.4 -


0.2 -




0.4




0 Z 4 6 8 10 12 14 16 18 20


Figure 3-6. Plot of [3-49]. Intersections of the two lines are solutions of
the characteristic equation.









and its coefficient will remain. We are now free to assume that the solution to the differen-

tial equation is a form of the first eigenvector and solve for the unknown coefficient.

Mt(0) -x2
S2EI F,(x)dx
an = L2 [3-53]
foF2(x)dx


where the thermal moment (Mt(t)) comes from Johnson and Thorten which is similar to

Reitz [66] shown below respectively.

Eawh2
MT(t) = h2AT [3-54]

1n272Kt
48EIqoc J4 0 h2
MT(t) = --- 4 e [3-55]


The AT in Thorten's equation is the temperature difference through the thickness of the

solar panel. An equation for AT, assuming the solar panel can be modelled as a thin plate,

will be derived later in the chapter. The temperature difference equation is dependent on

the physical shape of the solar panel, while the thermal moment equation is more general.

A quick numerical example, to calculate the coefficients for the first mode, set P3,L to

1.875104 and L to a position on the beam (w).

(1.875104)2 El 3L 1.875104
01 = 2 L w

sin(1.875104) sinh(1.875104)
cos(1.875104)+ cosh(1.875104)

Then solve for the coefficient an and finally for F(x). Although the derivations show that

the summation symbol goes from zero to infinity, in actuality you would only sum over a

finite number of modes.








Solving the time function G(t)

The solution of the time part of the separation of variables method is a little

complicated. The change of the dependent variable is what complicates finding the result.

Our solution for the time variable has the form shown below.

Gtt(t) + 02G(t) = -MTt [3-57]

The solution of [3-57] will be in the form shown below, where the term Gp is the

particular solution which can be determined if we have an equation for the thermal

moment. The general form of the solution is shown below.

G(t) = Asin(nt) + Bcos(ont) + Gp [3-58]

Summary of results

We have taken a long, winding road in order to solve the equations of motion for the

beam. The final solution of the beam, which describes the motion in one direction, will

need to be numerically solved is shown below.

M,(t)x2
u(x, t) = v(x, t) + p(x) = F(x)G(t) [3-59]
2El

Eawh2
M(t) = wh2AT [3-60]


F(x) = a,[cosh(Px)- cos(P3x)+ c,(sinh(P x) sin(P3,x))] [3-61]


(IPL)2 ELI InL sin (3nL)- sinh(P3,L)
Sn P = cos(3 +cosh(3nL [3-62]
L2 N P L cos(PnL) + cosh(PnL)


G(t) = A sin(ct) +Bcos(cmt) + Gp


[3-63]









2Mt(0) x2
an 2Ej X [3-64]
F F(x) dx


Again, summation occurs only over the number of modes of interest.

Thermally Induced Dynamics

Modal Equations

The preceding derivation results in the path of the beam at some point x over time.

However for the simulation modal equations are developed which are ordinary differential

equations (ODE). This section follows the ideas in McConnell [36] and Thompson [37]

and will cover the formulation of the modal equations. Further explanation of the theory

pertaining to modal equations can be found in Appendix B or Tongue [71].

Each mode of the system will be a second order equation as shown below.

qn + C,q, + Coq,, = Q, [3-65]

The subscript n denotes that the equation is for the n'1' mode. The q's are the generalized

coordinates of the mode and Q is the generalized force applied to that mode. Note that in

this formulation, the generalized force for the second mode will only effect the second

mode and have no effect on the first mode or any other mode. This idea is due to orthogo-

nality and as in the previous derivation, the modes are assumed to be orthogonal and do

not interact with each other. The relationship between our physical parameter of distance x

and the generalized coordinate q is as follows:


x = Fq, [3-66]
n









where F, is our position function1 for mode n from the previous derivation, note that

again it's dependence on x is not shown for simplicity.

L
Q = M P(x)Fndx [3-67]
0

where P(x) is some external distributed load per unit length. The mass, spring, and damp-

ener matrices are shown below:

L
M, = fpAF2dx [3-68]


L 2 2
Kn = fEl 2 dx [3-69]


C,n =a + 3,2 [3-70]

The dampening matrix is assumed to be the same magnitude as the mass and spring

coefficient matrices. Thus the dampening matrix is just a scaled version of the two where

a and P3 are scalar weighting coefficients that are tuned until the desired amount of

dampening is present in the system. The natural frequencies ((cn) are defined by [3-51]

and are shown again below:


(IPL)2 EI
mcn L [3-71]
n L2 N p

Disturbance Torques

In order to include the solar snap disturbance torques in a simulation, not all of the

preceding steps need to be performed. This section will follow along the work of Johnson

and Thorten [68]. Instead of solving [3-25], a differential equation is developed which

1. Thompson refers to F, as a normal mode and uses the variable %,, see page 438 in [37]









draws on the material already presented. This solution is easier to understand and integrate

into a computer simulation.

The equations Johnson and Thorten give for the coupled satellite solar panel system

using a generalized form of Lagranian equations is as follows:

L
Isat0 +fpA(Ro + x)u(x, t)dx = 0 [3-72]
0

An expression for the disturbance torque can be obtained by moving terms associated with

movement of the solar panel to the right side of the equation.The disturbance torque (g,,)

is a composite of two torques. The first torque (gqs) is generated from the quasi-static

terms, or the nonhomogeneous boundary conditions. The second torque (g,) is derived

from the vibrational aspects of the beam.

gss = gqs + v [3-73]

This can be seen by substituting the second derivative of [3-25] into the disturbance

torque equation obtained from [3-72].

L
gss = -fpA(Ro + x)u(x, t)dx [3-74]
0

d2 2
u(x, t) = u(x, t) = -(v(x, t) + p(x, t)) [3-75]


L d2 pA2ha RL3 3-76]
gqs = -fpA(Ro+ x) dd + -A [3-76]
0

L 2 N
= -pA(R,+x)2v(x, t)dx = pA [(Ro+x) (x)dx qn [3-77]
0 n =












T
T qsL
2k


T(L,t)

t


Figure 3-7. A cross section of the beam in Figure 3-8. Plot of the temperature
which a heat flux is applied to one side. change of a beam. Notice that the tem-
Notice the beam starts off at a uniform tem- perature difference between the two
perature. sides remains constant eventually.

where is the shape function (which use to be F, but is changed to match Johnson and

Thorten's work) and q, is the generalized coordinate for mode n, and AT is the tempera-

ture difference through the thickness (h) of the solar panel.

Thickness Temperature Gradient

The equation for the temperature in a thin plate of thickness h can be found in

Thorten's book [26] on pg. 86 or Incropera and DeWitt [69]. A diagram of the solar panel

cross section is shown in Figure 3-7, the equation for the temperature change through the

thickness is given below.

o -242Kt
qh xKt l(y + 1 2 1 2 (-1)" ( 1 -7
T(Y, = + 6 cos n2 + 2 eh [3-78]
k h2 2 h 2 6 2 h 2


The temperature difference through the thickness of the beam can be found by subtracting

the temperature at the heated surface from the temperature at the insulated surface at time

t. After subtracting the two equations and some algebra, the temperature gradient at time t,

steady state at time oo, and its acceleration becomes:


q(O,t)

T(x,O) = T
Sy

Sq(L,t)=0
SL ---











ko 2 7U2 o n2-
qfh 1 2 -2 h
AT(t) -I-+- -2 [3-79]
Sn = 1, 3, 5, ... -


qoh
AT(oo) =q [3-80]
2k


AT(t) = [3-81]
dt2

Notice that Figure 3-8 and [3-80] show the temperature gradient of the solar panel will

eventually reach a steady state. The two sides of the beam never simultaneously achieve

the same temperature while there is a heat flux constantly applied. This means there is

always a moment on the solar panels, and thus always a deflection of the solar panels.

Examples of this are shown in Figure 3-9 and Figure 3-10 where the solar panels of

Hubble are deflected due to the constant temperature gradient. Hubble's solar panels are

more complex than what is represented in this work. Hubble's solar panels incorporate bi-

stems, solar blankets, and spreader bars which are very thin, light weigh structures that are

highly susceptible to solar snap. Hubble's solar panels have been replaced several times

















Figure 3-9. Hubble Space Telescope which shows the solar panels bent
due to a constant temperature gradient.



























Figure 3-10. Another picture of the Hubble Space Telescope showing
deflection of the solar panels during maintenance by the space shuttle.

over the years, and the current ones are smaller and stiffer than the original design. This is

an attempt to reduce the influence of solar snap on Hubble.

Even though there is a constant temperature gradient through the panels, this does not

mean there is constant vibrations. The vibrations are driven by a change in the heat flux,

and thus a change in the temperature gradient. These changes primarily occur during the

day-to-night and night-to-day transitions of the satellite around the Earth.

Conclusion

Thermally induced vibrations or solar snap is a disturbance which occurs when a

spacecraft transitions into and out of the Earth's shadow. A temperature gradient occurs in

the long thin solar panel which in turn creates a thermal moment. This moment causes the

panel to vibrate until the transients die out. During this time frame, the scientific mission

has to be temporarily halted and the spacecraft put into a stable orientation.






45


The theory and math behind this phenomenon is not a mystery and the aim of this

chapter was to provide the reader with a complete derivation of this problem for a simple

solar panel. Here the assumption was made that the temperature gradient only occurs

though the thickness of the beam which reduced the problem to one dimension.















CHAPTER 4
SATELLITE ATTITUDE DYNAMICS

Satellite attitude control can be difficult due to the various types of disturbances and

limited control capabilities of satellites. In order to develop a controller, a good

understanding of the dynamics of the system is needed for proper simulations to be

conducted. This chapter will give an overview of the dynamics of rigid bodies in rotating

reference frames.

Reference Frames

Several different coordinate systems will be used to develop the equations of motion

for a satellite. This section will cover the reference frames used in this work and are shown

in Figure 4-1.

Earth Centered Inertial

The first reference frame is the Earth centered inertial (ECI) frame. This frame is a

non-rotating frame relative to the fixed stars1. Another way to say this is, the ECI frame is

an inertially fixed frame where Newton's laws are valid. This frame has its origin located

at the center of the Earth with its z-axis along the Earth's mean axis of rotation. The y-axis

and x-axis are pointing to some convenient set of stars. Note that the Earth is allowed to

rotate in this fixed frame.






1. Actually the stars move, and thus our ECI frame moves with them. This movement is small
enough to be ignored in most cases.















Z {Inertial}






Earh l


X

,- $ Reaction
"- Wheel Z


{ XX
Y
{LVLH}

{Body}


Figure 4-1. The three frames of reference commonly used in spacecraft
dynamics.



Orbital Frame

The orbital reference frame or local vertical, local horizontal (LVLH) frame is the

reference frame from which all of our spacecraft's angle will be defined. The orbital

frame's origin is located at the center of mass of the spacecraft. The z-axis points towards

the center of the Earth. The y-axis is perpendicular to the orbital plane. The x-axis points

in the general direction of travel. When the orbit is circular, then the velocity vector and

the x-axis are co-linear.

Body Fixed Frame

The body reference frame is located with its origin located at the center of mass of the

spacecraft. The x, y, and z axis are located along the x, y, and z principle axis of the

spacecraft.

Fun With Vectors and Rotating Coordinate Systems

This section will present some simple vector identities and properties which maybe

useful to the reader. Although all of the important steps are shown in the following


x


f


I


\x









derivations for spacecraft dynamics, at times some algebraic simplifications are not

shown.



a-(b+c) = ab+a-c axb = -bxa

(axb)-c = a (bxc) a-b = b'a
ax(bxc) = b(a c)-c(a b) a-b = aTb


Conventions used:

* All vectors and matrices are bold face.

* 1 and 0 represent the identity matrix and the zero matrix respectively.

The dynamics for a spacecraft are written with respect to several frames of reference.

Two important reference frames will be the ECI inertially fixed frame and a body fixed

frame which is allowed to rotate in space. Differentiating a vector relative to a rotating

reference frame is given by:

0A1 dAB
-a = d + oxAB [4-1]
at dt






Y p











Figure 4-2. A reference frame in motion relative to an
inertially fixed frame of reference.









where (o is the relative angular rotation rate between the inertially fixed frame I and the

rotating body frame B.

Now for the point located in the rotating system at position r, its velocity is:

dr1 dR, dr' dr's
V + = R+- +coxr'B = R+v+coxr'B [4-2]
S dt dt dt dt

r = R + r' [4-3]

where R is the velocity of the origin of the rotating system relative to the fixed system and

vB is the velocity of the point in the rotating system. A similar derivation can be shown for

the acceleration of the point which results in:


a, = aB+R+2c xvB+c x(coxr')+ Ixr' [4-4]

When dealing with rigid bodies, with the rotating frame representing a body fixed

frame, the point does not move within this frame. This is due to the fact that is fixed in

the body. Thus vB = a = 0, and the proceeding equations reduce to:


v = R+ x r' ai = R+cox(coxr')+ Cxr' [4-5]

where the velocity and acceleration are functions of the rotation rate of the rotating frame

and its velocity and acceleration.

Spacecraft Attitude

Attitude refers to the orientation a spacecraft occupies in 3D space. Although there are

several different ways to represent this (i.e., euler angles, Gibbs vector, fixed angles, Euler

symmetric parameters, etc.), typically space applications utilize quaternions. This section

does not attempt to provide the extensive understanding needed to employ quaternions but









rather a simple introduction. Further information can be found in Wertz [37] or Crane and

Duffy [42], or Appendix A.

Quaternions

Quaternions were invented by William Rowan Hamilton in 1843. Prior to his

discovery, it was believed impossible that any algebra could violate the laws of

commutativity for multiplication. His work introduced the idea of hyper-complex

numbers. Here real numbers can be thought of as hyper-complex numbers with a rank of

1, ordinary complex numbers with a rank of 2, and quaternions with a rank of 4.

Hamilton's crucial rule that made this possible:

i2 = j2 = k2 = ijk = -1 [4-6]

Hamilton supposedly developed this rule while on his way to a party. When he

realized what the solution was, he took out his pocket knife and carved the answer into a

wooden bridge. This rule would forever change mathematics as was known at the time.

Now mathematicians could look at algebra where communitivity did not work. This is

where Gibbs and others developed algebra of vector spaces, and quickly eclipsed

Hamilton's work until recently.

Quaternions, also known as Euler symmetric parameters, are more mathematically

efficient ways to compute rotations of rigid and non-rigid body systems than traditional

methods involving standard rotational matrices or Euler angles. Quaternions have the

advantage of few trigonometric functions needed to compute attitude. Also, there exists a

product rule for successive rotations that greatly simplifies the math, thus reducing

processor computation time. Quaternions also hold the advantage of being able to

interpolate between two quaternions (through a technique called spherical linear









interpolation or slerp) without the danger of singularities, maintaining a constant velocity,

and minimum distance travelled between points.

The major disadvantage of quaternions is the lack of intuitive physical meaning. Most

people would understand where a point was in cartesian space if they were given [1 2 3].

However, few would comprehend where a point was if given the quaternion [1 2 3 4].

Rotations of Rigid Bodies in Space.

Quaternions are able to represent a rotation of a rigid body in space. To perform a

rotation (|) of a rigid body about an arbitrary moving/fixed axis (e) in space, the

quaternion representation of this operation is


q = [q4] Norm(q) = 1 [4-7]


q =e sin( q4 =cos( [4-8]


e = [e e2 eT Norm(e) = 1 [4-9]


Notice that only one sine and one cosine function call is needed to calculate a

quaternion, where an euler rotation matrix would require three sine and three cosine

function calls, one each for roll, pitch, and yaw. Since trigonometric function calls are

computationally expensive, this is a great savings.

Spacecraft Equations of Motion

In this section, the kinematic and dynamic equations of motion for a spacecraft are

presented.









Attitude Kinematics

A spacecraft's orientation in space is represented by a quaternion. The quaternion

kinematic equations of motion are given by Wertz [38] as,

1 1
q = G2q = E(q)co [4-10]
2 2

where



S= -Ox and E(q)= 433x3+q [4-11]
I-T 00^ T
-q

The terms with the x by them signify a skew matrix.

Spacecraft Dynamics

First the equations of motion will be developed for a satellite and then expanded to

include other dynamics. Thus a lengthy derivation is needed to derive something as simple

as Euler's equation. This, however, is necessary so that when reaction wheels and fuel

slosh are included, understanding how these dynamics are incorporated will be easier.

Simple Spacecraft

To model a satellite, we starting off with a couple of definitions for the total mass of

the spacecraft, first moment of inertia, and the second moment of inertia respectively.

N N N
m = mn c mnrn J mn(rTrnl-rnr ) [4-12]
n= 1 n= 1 n= 1

where mn is a point mass in the spacecraft, Nis the total number of points. These defini-

tions are in reference to Figure 4-3 (Adapted from Hughes [31] p. 43). The second

moment of inertia (J) is typically just referred to as the moment of inertia. The moment of

inertia has two important properties: symmetry and positive definite.














*^, (m.)

Ln


Imn



.-0' (m,)


Figure 4-3. A rigid body satellite.

The momentum and velocity of each point in the spacecraft is

N
mnvn = fn + n vn = V+rn [4-13]
m= 1

where R = vo, f, is the external forces acting on point, and fmn is the force exerted by

point m on point n.

The linear momentum (p) of the total spacecraft is

N N N
P = pn= mnvn n= v m +rn) = mv+c [4-14]
n=1 n=l n=1

This equation could be further simplified by realizing that a spacecraft1 is a rigid body,

and thus c = 0. Thus

N
p =f= Ifn [4-15]
n=



1. Actually since spacecraft are composed of light weight, flexible materials, some authors model
them as flexible structures.


*Y,(m,)









where

N N

SI fimn =0 [4-16]
n=lm=l

which comes from Newton's third law, fnn +fnm = 0.

The angular momentum of the spacecraft is:

N
ho= rn xPn [4-17]
n=

Now differentiating this equation with respect to time, results in the following equation.

N N
ho= (rn xpn + rxpn)= [(vn -v)x mn'n+ r xfn] [4-18]
n= n=1

Realizing that vn x mnvn = 0 results in the following equation.


ho = -v xp +go [4-19]

where go is the total external torque about the point O.

Now if the equations are defined about the center of mass of a rigid body spacecraft,

then c = 0, vo = vn, and the cross product with momentum in [4-19] disappears. The

equations collapse to the more familiar ones for linear and angular momentum.


p = mvy ho = go [4-20]

Now using the previous equations involving rotating reference frames, we can write

the linear and angular momentum equations for a spacecraft in a rotating body fixed

frame.


p = mvo +c+ xc [4-21]


p = -oxp+f


[4-22]









-mxho- v xp +go


As before, these equations can be simplified using certain assumptions. However, deriving

the equations without the assumptions makes including internal dynamics (i.e., fuel slosh

and reaction/momentum wheels) easier.

These equations put into a matrix form are as follows:


p= p h V= [v]T F= f goT [4-24]


p = MV M = ml -cx [4-25]


where Mis the mass matrix of the system. The kinetic energy Tof the system is given by:

1 T
T = VMV [4-26]
2


T = mVyo + oWT x vo + -TJom [4-27]


The dynamics can be obtained from the kinetic energy by using the following:

aT aT
p = ho [4-28]
av am

daT aT = d aT +aT (a T [4-29]
+ x -- =f -t ) + mx a+v )x = go [4-29]
dtK8v) K9v) dt\m) m v 9J

These equations are called the quasi-Lagrangian equations and are shown in Hughes [31]

and Meirovitch [39]. Meirovitch derives the equations using quasi-coordinates, which is

somewhat difficult to understand and very long. However, this method has been used by

some authors such as Miller et al. [40] to model various flexible appendages on a satellite

and moving submasses inside of a satellite.


[4-23]









Angular Velocity

The spacecraft's orbit in this work is assumed to be circular. Since we are using an

LVLH reference frame for our spacecraft, the angular velocity of the system must account

for the orbital rate in addition to the maneuvers that the satellite is conducting.

COB/I = OaiB/ + apol [4-30]


col1 = -Cn n = T [4-31]


where B refers to the body frame, O refers to the orbital frame, I refers to the fixed inertial

frame, a/10 refers to angular velocity of B relative to O, C is a rotation matrix, [t is the

gravitational parameter of the Earth, and Rc is the distance the spacecraft is from the cen-

ter of the Earth.

Internal Disturbances

Various factors complicate the dynamics of a satellite's dynamics. Moving masses

such as reaction wheels and sloshing fuel create additional linear and angular momentum

that needs to be accounted for. This section will show how the above equations for a

simple spacecraft need to be modified in order to account for these effects.

Spacecraft with Reaction or Momentum Wheels

Starting off with reaction wheels, we will introduce a method to account for any

number of moving wheels in any configuration. However, since this work only deals with

reaction wheels, we will not concern ourselves with momentum wheels which change

their orientation and not their rate of spin.











aV












o2,



Figure 4-4. A rigid body spacecraft with a single reaction wheel.


Now we will add the contribution of the reaction or momentum wheels to the

spacecraft dynamics. Looking at Figure 4-4 (adapted from Hughes [31] p. 66), we can

define the following terms:

m = msat + m c = csat +mb J = Jsat +Ia + m(bTb bbT) [4-32]

where terms with the subscript sat refer to the satellite and terms with the subscript w refer

to the wheels. The linear momentum is given by the following equations:


Psat = msatv + cmx sat [4-33]

p, = mv + m, x b [4-34]


P = Psat +Pw = my + x c [4-35]

where b is the vector from the point O to the wheel. The angular equations of motion are

given by the following equations:

hsat = sato+ csat xv [4-36]









h, = IsCo. [4-37]

h =hsat + ha [4-38]

h = cxv + Jo+ ah, [4-39]

where a is a vector that describes the direction of the wheel axis when dealing with a sin-

gle wheel. When there are multiple wheels, a becomes a matrix where the columns are the

individual direction vectors for each wheel. For example, a system with n wheels would

look like:

a = [{al} {aw2} {awn [4-40]

where

awi = e2 e 3]T Norm(aw) = 1 [4-41]

The wheel speed co, is a vector where each element is the speed of a specific wheel. Thus

for a system with n wheel, the vector would be:


O = [Owl ,w2 ... wnr [4-42]

The equations of motion can be derived from kinetic energy for the system by using

the quasi-Lagrangian. The kinetic energy for the spacecraft with wheel dynamics is:

T= mvTv+ -CTJ w+ -I ~mwlC x m+IwmoTaTC [4-43]


Putting these equations into a matrix form results in the following equations.


S= [ph h, h V= [vo o o]T F= fg ga] [4-44]










ml -c x 0
S= MV M= c x J Ia [4-45]
OT IaT I _,aT

The equations of motion, like before, are as follows:


p = -coxp+f [4-46]

ho = -wx ho-voxp+go [4-47]


ha = ga [4-48]

Although these equations look the same as the previous equations without wheel dynam-

ics, remember the equations for linear and angular momentum are different. The linear and

angular momentum of the spacecraft and wheels are now coupled in these equations.

Equations for spacecraft with wheel dynamics are often simpler than what are shown

here. However, if the wheels are assumed to be located at the center of mass (which is

physically impossible) then the equations are greatly simplified,

= -w x ho + g [4-49]

since p = mvo and vo xp = vo x myv = 0. Further more, h = Jco+h since c = 0

and substituting this into [4-49] results in:

Ji +h, = -o x (Jw +h) + go [4-50]


Jm = x Jo + x h hw + g [4-51]

JC = -ox Jw + o x h, ga + go [4-52]

This final equation is the common equation found in Wie [22] and Sidi [9].







60


Fuel Slosh

Fuel slosh, as stated in a previous chapter, is the unwanted movement of fuel inside a

spacecraft. Each sloshing mode is modelled as a mass, with the forces exerted on the

spacecraft body represented by a spring, and the effects of baffling in the fuel tank






5.0


Geostationary Orbit
4.6



-Cosmic Dust
S CM-CP= Gl m
4.2 CM-C m Gravity Gradient
E 0 ei.o*

ILJ
O _, Magnetic (Equator)
38 1.0 Amp-Turn
I-
rSpecular Reflection-
Diffuse Reflection-
Solar Radiation
CM-CP=& I rm
0 3-4 \




3.0

Aerodynamic

unspot Min. ym
2.6 L' ight /
nspo Mox. M eight



2.2


10'7 10-5 10"3 10-'
TORQUE, [N-m]

Figure 4-5. A plot of common environmental disturbances a space-
craft is subjected too.









represented by a dampener. The dynamics of this are incorporated into our simple satellite

model by the following equations which will represent one sloshing mode.


p = my-c x w+mpn [4-53]


h = cxw+Jw+ mfb xn [4-54]


pf = m(v- b x o+ n) [4-55]

where n is the position of the sloshing mode, mf is the mass of the fuel, and b is the vector

from the center of mass of the spacecraft to the center of the fuel tank. These equations

reflect how the linear, angular, and fuel slosh momentum change with the addition of fuel

slosh's mass-spring-damper. The mass matrix now becomes:


ml -c x mf
M = x J mb x [4-56]
mf -mjbx mf

As with the wheel dynamics, the differential equations of motion for the linear and angular

momentum remain the same. The new differential equation for the fuel slosh is:


pf = -mW x (v rd x w) cy kin [4-57]

where the cf and kf are the dampening coefficient and spring constant for the sloshing

mode.

Environmental Disturbances in Space

There are many environmental disturbances that plague a spacecraft in orbit. These

disturbances constantly effect the performance of the spacecraft which necessitates the use

of a control system that counter acts them. The disturbance torques produced by some of









these are shown in Figure 4-5 (adapted from Hughes [31] p. 271). A further explanation of

the important disturbances follows.

Gravity Gradient

Euler and d'Alemert in 1749 first pointed out how the Earth's gravitational field was

not uniform. Later in 1780, Lagrange used these ideas to explain why the moon always

had the same side facing the Earth. Since the gravitational field over a rigid body is not

uniform, the center of mass is not the center of gravity. Thus there are torques about the

center of mass which depend on the orientation of the spacecraft.

The gravity gradient torque disturbance is found in [22], [31], and [38] as


ggg 3m 23 lo [4-58]

where ome is the orbital rate of the spacecraft, 03 is the third column in the orbital frame to

body frame rotation matrix. When a satellite is in low Earth orbit, gravity gradient torques

become the predominate environmental disturbance.

Solar Snap

In a previous chapter, modelling solar snap was discussed. The disturbance torque due

to solar snap is a combination of the quasi-static displacement of the solar panel and its

vibration.

gss = gqs + [4-59]

Conclusion

This chapter showed the basic concepts of spacecraft and disturbances dynamics.

These dynamics produce a system where the satellite body, wheels, and fuel are all

coupled. The disturbances presented in this chapter were internal (fuel slosh) and external






63


(gravity gradient). A complete system can now be produced which has all of the dynamics

of interest: satellite, wheel, fuel slosh, and thermally induced vibrations.















CHAPTER 5
ATTITUDE CONTROL OF SPACECRAFT

This chapter will introduce the reader to attitude control of a satellite and specifically

to a variable structure method called sliding mode control. First however we will cover

some of the typical hardware used in a satellite to control its attitude. Next the traditional

proportional derivative controller and sliding mode controllers will be introduced. Later in

this chapter we will compare and contrast the two controllers.

Satellite Control Hardware

This section will provide a brief overview of hardware used by a satellite to control

attitude.

Control Computational Hardware

Satellite's are typically equipped with only the bare minimum in computational ability

in an effort to reduce the cost of the vehicle. Unfortunately the complexity and

performance of the control system is constrained only to the simpler algorithms. In fact

some spacecraft, such as the lunar explorer, had no on board control systems. Instead it

was remotely controlled from Earth, which was a very brittle solution at best. This was

hailed as a great step forward in reducing the cost of space flight, but what happens if you

lose the connect between ground control and the satellite?

Any control algorithm used must be computationally inexpensive so it does not

overburden the processor, and is typically designed to operate at 1 Hz. One of the

simplest, easiest, and most used controller available is the venerable proportional









derivative (PD) controller. This is one reason that the PD has flourished in spacecraft

control.

Control Actuators

In order to overcome environmental disturbances, satellites are typically equipped

with various types of actuators that can provide an attitude control effort. Some of these

are active control systems which utilized electrical power or consumable resources to

effect the spacecraft's attitude. While others are passive, and rely on conservation of

energy to effect the satellite's motion. Passive methods of control are generally relegated

to dampening out disturbances during flight. These methods do not produce the same

levels of torque as their active counterparts, but are typically employed along with active

controls in an effort to reduce the power consumption of the spacecraft.

Some of the most common types of active control actuators will now be discussed. A

quick comparison of their control torques is given below.

Table 5-1: Comparison of common active attitude control hardware.
Control Hardware Control Force Range (Nm)
Magnetic Torque Rods 10E-6 (geostationary)
2.5E-3 (400 km)
Reaction Wheels 0.05-2
Thrusters 0.1-30

Reaction and momentum wheels

Reaction and momentum wheels are cylindrical masses that rotate about their center of

mass driven by a motor. They are produced in a variety of types and sizes, some of which

are shown in Figure 5-1. These devices are capable of changing or stabilizing the

orientation of a spacecraft with respect to its axis of rotation through conservation of

angular momentum. These are the best for delivering smooth continues control efforts.





















.i.
un *, .FfP
i. 5- .
*Iwr-v --- **:r i,
*ijwii

\ IRONLESS ARMATURE
BRUSHLESS DC MOTOR

Figure 5-1. Two types of reaction wheels produced by Ithaco.


Wheels that have a nominal spin rate of zero1 are called reaction wheels. Wheels that have


a constant momentum are called momentum wheels.


Moment wheels are used to create a bias in the system which resists small parasitic


environmental torques. Sometimes these wheels are also mounted in a gimbaled system


that has one or two degrees of freedom. This allows the momentum wheels (with their


constant momentum) to change the direction of the their momentum vector.


Magnetic torque rods


Magnetic torque rods are another method of attitude control and are shown below in


Figure 5-2. They are based on the idea of generating torque by magnetic fields. Here the


rods produce a magnetic field which interacts with the magnetic field of the Earth.


Obviously this method of attitude control produces small control forces compared to


1. This statement is not always true. Sometimes wheels are biased by a small amount to fight fric-
tion in the wheels and bearings. These zero-crossing problems, dead zone around zero rpm, can
effect the accuracy of a satellite when it is gathering scientific data.






















41In) (151in)

Figure 5-2. Torque rods produced by Ithaco.


reaction wheels. Also the maximum amount of control effort produced by this method is

dependant on the altitude of the spacecraft since the Earth's magnetic field is stronger the

closer the spacecraft is to the center of the Earth.

In addition to small torque capabilities, electromagnetic shielding to protect satellite

components must be considered. This has the additional disadvantage of adding weight

and increasing the complexity of the design process.

Gas jets

Gas jets or control jets (shown in Figure 5-3) are the most powerful hardware available

for attitude control, but these are not the most desirable to use for two reasons. First they

consume fuel which is a limited, non-renewable resource, and second they do not have

variable control but rather "on-off'. Thus bang-bang control schemes are created to

provide pulse-width modulation control efforts to mimic a variable type of device.

Proportional Derivative Control

First let's start off with one of the simplest to design and computationally efficient to

implement, the proportional derivative (PD) controller. Wie and Barba [41] developed

several computationally efficient control schemes for large angle maneuvers. Many of






















Figure 5-3. Spacecraft engine produced by Boeing.

these stabilizing schemes utilize quaternion and angular velocity feedback. The use of

quaternion representation allows for more realistic, large angle maneuver control schemes.

These schemes are formulated based on Lyapunov analysis, which produces a range of

positive stable gains for that control law. Thus in order to meet desired performance,

engineers must iterate through a significant number of gain combinations to obtain the

desired response in a clean simulation. However, even when a satisfactory response is

finally obtained, there are no guarantees how the satellite will behave in the presence of

disturbances, noise, or uncertainties.

Bong Wie [22] defines three different PD controllers utilizing quaternions in his book

where the origin of the quaternion system is (0,0,0,1).


u = -Kq-C C [5-1]

K = kl C = diag(c1, c2, c3) [5-2]

k
K = C = diag(cl, c2, c3) [5-3]
q4

K = ksgn(q4)1 C = diag(c1, c2, C3) [5-4]

where u is the control effort and all k's and c's are positive.









Sliding Mode Control Theory and Background

This section will present an overview of sliding mode to the reader. Further in depth

explanation can be obtained from Slotine [42] or DeCarlo et al. [43].

What is Sliding Mode?

Nonlinear model based control systems offer a level of dynamic capabilities which

linear techniques are incapable of providing when dealing with parameter uncertainties

and unmodeled dynamics. Sliding mode, which has been studied in the Soviet Union for

many years, is categorized as a variable structure control. This control structure has

excellent stability, robustness, and disturbance rejection characteristics. Sliding mode is

not a new control technique either, many researchers have utilized its robust properties to

control a variety of systems such as missiles [44], mechanical systems [45], robotic

manipulators [46][47][48],and submarines [49][50][51][52].

Theory

First let's look at a simple example of controlling a system. Given the system below,


y = u(t) [5-5]

u(t) = -k y(t) [5-6]

a solution to this differential equation can be found and is plotted in Figure 5-4. Notice

that neither of the two solutions drive the system to zero at the origin, but rather remains

oscillating with a stable behavior. The trick is to switch between the two solutions depend-

ing on which quadrant the system is in.


u(t) = -k y(t) if yy < 0[57]
-k2 y(t) otherwise












y


(a) (b)
Figure 5-4. Two possible solutions for the double integrator, where (a) is
k = k1 and (b) is k = k2. These plots are typically referred to as phase
portraits. They plot position and velocity in this example, but later (when
we are doing controls with full state feed back) the portraits will represent
error and error velocity.

By switching, we can drive the system to the origin in a stable manor, as shown in

Figure 5-5. Slotine took this idea and developed a systematic method for designing sliding

mode controllers. Looking at Figure 5-5, we can develop the following table:



Table 5-2: Summary of the switching law.


Quadrant y y u(t)

1 + + -k2

2 + k

3 + kl

4 k2


This table can be summarized by the switching law.


sgn(s) =-1 s(y,y)>0
l1 s(y, y) < 0


where


s(y,y) = y+ [y


[5-8]


[5-9]














yy





Figure 5-5. Switching between the two controllers results
in driving the system to zero in a stable manor.

is called your sliding surface and the control effort is given by:

u(t) = -K- sgn(s) [5-10]

where K is typically a scalar or diagonal matrix of positive definite gains.

In order to design a controller, we will now replace our states (y, y) with our errors

(x, x). This will now (as previously stated) drive our system errors to zero in a stable

manor. Second, we will modify the control effort equation to include a feed forward aspect

as follows.


u(t) = u-K- sgn(s) [5-11]

The term u is called the equivalent control effort. It is model based, thus if we have a

perfect model of the system, we can calculate the required control effort to produce the

desired results. However, that is never the case. There are always disturbances and

unmodelled dynamics that influence our system. Thus the job of the second term in [5-11]

is to reduce the error and error velocity produced by the equivalent control in a stable

manor. Also note that the X is a weighting factor between the error and error velocity. This

term can be adjusted depending on which is more important, and its effects are depicted














s ,,x>i
S Y s=O, k>1 y




s = O, 0 = 45' \
x = 1 s=0, x<1
Figure 5-6. Graphical representation of how x effects the sliding surface's
orientation in the phase plane. The dotted lines represent the sliding surface (s
= 0).


below in Figure 5-6. Typically this term is constant, but there is no reason it could not be

dynamic.

Simple Example

For a simple example, we will develop a sliding mode controller for a mass-spring-

dampener system. The dynamic equations for this system is

mx + cx + kx = u. [5-12]

Our sliding surface and equivalent control are


s = x + x, [5-13]


s = x +x = 0 [5-14]


m-(u -cx-kx)-xd + Xx = 0 [5-15]


u = m(x )+cx + kx [5-16]


Now looking at [5-16] we notice that if there is no error in the modelling of the system

or disturbances, then this control effort is sufficient to control the system. Also if we had a

model of any disturbances in the system, we could include them in our equivalent control









term. However in this example, there are no disturbances included in the system or we do

not have a model for any disturbances seen by the system during run-time. Unfortunately

there will always be error in modelling, noise, or disturbances. Thus we must include

another term to account for these problems.


u = u-K- sgn(s) [5-17]

Other Sliding Mode Designs for Spacecraft

Many authors have suggested many types of control architectures for satellite attitude

control [53-59]. However, the one of the most interesting and simpler to implement types

of control is sliding mode. The use of sliding mode control is not new to satellite attitude

control. Lo and Chen [60] designed a sliding mode controller scheme which avoids the

inverse of the inertia matrix, and smooths the control effort of the controller. Such a

strategy provides for a more efficient use of fuel. Their simulations involved small

uncertainties in the spacecraft inertia matrix and a sinusoidal disturbance. Simulations

showed favorable results. The control effort was


u = ox Jw +Jvd Jq- ksgn(s) [5-18]

where the vd a dd for all intensive purposes.

Boskovic et al. [61] and [62] developed a sliding mode controller, but specifically

designed the controller such that it did not saturate the control effort. They were able to

produce accurate results, and even tested the controller on a system with a larger inertia

matrix than what the controller was designed with. Unfortunately they drew an incorrect

conclusion that the controller is independent of the inertia matrix it is presented with.

Sliding mode will be able to control a system with a larger inertia matrix, but may not

control a system with a smaller inertia matrix. This is because the control effort is model








based, and thus it is possible to make a satellite go unstable by issuing too much control

effort. The simulations utilized a square wave disturbance and small uncertainties in the

inertia matrix. The controller developed was as follows:


ki
u = xJ XJq -diag sgn(sl) sgn(s2) sgn(3)] [5-19]



Crassidis et al. [63] developed an optimal variable-structure (sliding mode) controller

for large-angle maneuvers on satellites. A cost function was developed, which when

minimized, resulted in the optimal control effort. The the minimization lead to a two point

boundary value problem. The cost function and the resulting optimal control effort were

1 ~T~
n(co) = ft (pq q+ o)dt [5-20]


u = xJ( + J(X sgn(q4)[E(q)E(qd)(d- (qd)E(q)W(]) +d ksgn(s) [5-21]


Simulations of the controller on the MAP satellite showed favorable results.

Coleman and Godbole [64] conducted a performance trade study between fuzzy logic,

PID, and sliding mode control. The controllers were tuned for and tested on three similar

linear plants. In all three cases, the sliding mode controller out-performed the fuzzy logic

controller. This is a typical result, where sliding mode tends to provide a superior model

based performance compared to fuzzy logic, assuming the model is accurate enough.

However all of these authors utilized simple satellite dynamics with a simple external

disturbance in their simulations. Fuel slosh, which is coupled with the dynamics of the

satellite, is a more complex problem. It is difficult to accurately model the real disturbance

and requires a control system that can account for this uncertainty.









Satellite Controller Design For This Work

This section will cover the sliding mode controller developed for this research.

Controller Derivation

This formulation, which is a full-state feedback technique, utilizes the existing

dynamic model and compensates for uncertainty while formulating a control effort that

tracks the desired trajectories. The equations of motion for a satellite's attitude (in

quaternion space) are given by [22]

Jio = -QJm +f+ u [5-22]

1 1
qq + -q4m [5-23]


where


0 -03 (02 q
2 = 03 0 -( q = q [5-24]

-_(2 (0 0 q

Also f in [5-22] represents the modelling error and potential disturbances in the

system. This term will provide a feed forward aspect to the controller and should yield

better results. The q is a vector composed of the three imaginary elements of a quaternion.

However, in this formulation, we can take the state variable to be A or small euler angle

orientation, and equate this to q.



6y = 2 qy [5-25]



Thus the sliding surface (s) is given by:










s = + x = x + x = o+ kq [5-26]


where the error terms are defined as:


X = X -Xdesired [5-27]

In order to calculate the equivalent control effort (u), we need to take the derivative of

the sliding surface and set it equal to zero.


s = co+ q [5-28]

Now multiplying both sides by the inertia matrix (J) will allow us to substitute in the equa-

tions of motion ([5-22] and [5-23]).


Js = JC+ Jq = Jo+ Jq Jdesired [5-29]


Js = QJa+f+ + +Jq J desired [5-30]

The hat over particular variables denotes that they are estimates, since accurate values

for the inertia matrix and disturbances may not be known. Since the sliding surface does

not move, its time derivative is zero, and thus the equation becomes


Js = J+f+u + Jq-JOdesired = 0. [5-31]

Now solving for the equivalent control effort (u) from [5-31] yields:


u = -QJ- f- Jq + Jmdesired. [5-32]

The quaternion error rate in the above equation is defined as:

1~ ~ 1~ ~
q 2q + -q4co [5-33]
2 2










0 -(03 (02
03= 03 0 -C [5-34]

-(2 (1 0

The estimated disturbance in the above equation is equal to the estimated solar snap and

fuel slosh disturbance.


f = fss +ffs [5-35]

where f,, is the disturbance due to solar snap and ff, is the disturbance due to fuel slosh.

Finally, the sliding mode control effort is given by:

u = u-K-sgn(s) [5-36]

The K term in [5-36] will be discussed when stability of the controller is covered. This

switching term provides the appropriate control action to quickly drive the trajectories

back onto the sliding surface. Since real actuators do not have an infinite bandwidth, the

sign function tends to cause excessive chatter in the control effort. The sign function is

typically replaced by a saturation function. The saturation function not only smooths the

control response, it also reduces the amount of fuel used by the controller.

Fuel Slosh Disturbance

The effects of the fuel slosh are incorporated into the controller by the ffs term in the

controller. A discrete state space system is used to track the effects of the fuel slosh:

xfs = Axfs +Bu [5-37]


ffs = Cxfs [5-38]


where the state vector is:










Xf= x [5-39]
x

where x is the position of the fuel slosh and x is its velocity. The input to the system u is

the previous control effort used by the satellite. Also the state space matrices are discrete

time representations to reduce the computational load of the controller.

Stability of Controller

The following Lyapunov candidate function is proposed. This function is clearly

positive definite.


V sTJs [5-40]
2

Taking its time derivative and substituting in previous equations of motion and control

effort gives:


V = sTJs [5-41]

where the true dynamics of the system are inserted into the equations are


Js = -c x Jo + u +f+ hJq [5-42]

The modelled dynamics, which define the controller output are given by:


u = u-Ksgn(s) and u = cxJcw-kJq [5-43]

After some algebra, the dynamics cancel each other out, and the only terms left are the dif-

ference between the true disturbance forces and the estimated/modelled ones and the last

terms on the control effort equation. Thus [5-41] becomes:


V = sT(F-Ksgn(s)) [5-44]


where









F = f- [5-45]

Note that any differences between the true inertia matrix and the real inertia matrix, plus

any other modelling error or disturbance, is contained in f Where f is the true distur-

bances, modelling errors, unmodelled higher order terms, etc.

For a stable system, the lyapunov derivative needs to be negative definite. All

variables are known except for the scalar gain value K. Two of the three remaining terms

are positive definite, K and F. Thus we need to do further simplifications and algebra to

solve for the gain.

sTF sTKsgn(s) < 0 [5-46]

sTF-K(s) <0 [5-47]

where

(s) = sTsgn(s) [5-48]

This operation is a cross between an absolute value operation (if it were a scalar) and a

norm of a vector (if the sign function was not present). Thus this operation will be denoted

by the angled brackets and not confused with either absolute value or norm of a vector. A

quick example is shown below, and note that the result will always be scalar and positive.



1 -3 4] sgn -3 -3 4] = 8 [5-49]


Thus the scalar gain K is

K sT F
K > [5-50]
(S)








In order to guarantee the above inequality, we will add a little offset (r1) that will always

be positive.


K +s [5-51]
(s)

Now substituting this gain back into [5-47] results in:


sTF ( F + (S)

This equation has two results depending on whether s is positive or negative. If s is nega-

tive then the two sTF's sum together and the entire equation is negative definite. If s is

positive, then the two cancel each other out and the equation is still negative definite.

2sTF-q(s)<0 or -q(s)<0 [5-53]

Saturation

The sign function used in the sliding mode controller leads to an excessive amount of

chattering in the control effort. A common modification to the sliding mode controller is to

change the sign function to a saturation function with the following properties.


sat(a) = < a < a [5-54]
else sgn(a)

Thus our sliding mode control equation changes to


u = u-Ksat( [5-55]


where the ( is a scalar value that can be used to adjust when the saturation function will

saturate. All other sliding mode equations remain the same.









Comparing The Two Types of Controller

Slotine's sliding mode formulation and a PD controller have some very similar

attributes when the sliding mode is in the sliding region. Simplifying the system some,

assume that the equivalent control effort (u) is zero, since this is really a feed forward

term. Thus the control effort is now:


u = -Ksat( [5-56]


The sliding surface s is very similar to a PD control effort.


s x = x + kx [5-57]


where x is a rate error and x is a position error. When the errors are on the sliding surface

(meaning when the errors are small, they are in the saturation bounds), the saturation func-

tion returns a scaled version of the values passed to it (as opposed to a 1 or -1 when out-

side of the saturation bounds). Thus the control effort becomes scaled errors and scaled

rate error multiplied by a gain K.


u=-K [5-58]



u = --(x + x [5-59]


Then the PD parameters can be written as follows:

hK 1
k T, [5-60]


Conclusion

This chapter introduced PD and sliding mode theory for satellite attitude control. An in

depth discussion of sliding mode design, stability, and a comparison between sliding mode






82


and PD was provided. This laid the foundation for the similarities and differences between

the two controller which will soon manifest themselves when the results are presented.















CHAPTER 6
SIMULATION AND RESULTS

This chapter will discuss the simulation, the gains chosen for the PD controller, and

finally discuss the results. The first section will describe the specifics of the satellite and

the simulation. All of the variable from proceeding chapters will be filled in and

explained. Next the simulation, equations of motion, and the controllers that were

developed are covered. Finally the results of the simulation are discussed.

Simulation Dynamics and Parameters

The simulation used to evaluate the performance of the controllers was written in C++.

The first simulations were written in Matlab because of the ease of working with vectors

and matrices. The problem with Matlab, however, was the amount of time each simulation

took. Once the simulation were rewritten for C++, they were magnitudes faster. The reader

is referred to Appendix C, where some of the mathematical code developed is presented.

Satellite

The satellite shown in Figure 6-1 and Figure 6-2 is Earth Observing 1 (EO-1), and will

be the satellite used in all simulations. EO-1 was primarily constructed by Swales











Figure 6-1. EO-1 satellite with solar panel deployed.























Figure 6-2. Diagram of satellite's internal components.


Aerospace and is part of the New Millennium Program which is dedicated to validating

new technologies. The spacecraft was launched aboard a Boeing Delta II rocket from

Vandenberg Air Force Base on November 21, 2000. Some of the new technology being

tested is the Advanced Land Imager, hyperspectral imaging spectrometer called Hyperion,

X band phased array, pulsed plasma thrusters, and formation flying with Landsat 7

spacecraft which is already in orbit.



Table 6-1: Technical Specifications for the EO-1 spacecraft
Dry Mass 568 kg
Volume of Bus 1.5m x 1.5m x 2m
Inertia Tensor [443 179 429] kg m^2
ACS Zero Momentum, 3 axis stabilized
Command and Data Han- Mongoose V, Rad Hard at 12 Mhz,
dling Bus Architecture RISC Architecture
Solar Arrays 3 panel / Si w/GaAs / Articulating
Bus Structure Hexagonal with aluminum honeycomb
Propulsion 1 fuel tank with 4 thrusters
Propellent Capacity 23 kg
Mission Design Life 1.5 years