Stability of representative -pairs.
This chapter demonstrates analysis and control of the attitude motion of a gravity-gradient stabilized spacecraft in eccentric orbit. The attitude motion is modeled by nonlinear planar pitch dynamics with periodic coefficients and additionally subjected to external periodic excitation. Consequently, using system state augmentation, Lyapunov-Floquet (L-F) transformation, and normal form simplification, we convert the unwieldy attitude dynamics into relatively more amenable schemes for motion analysis and control law development. We analyze the dynamical system’s periodicity, stability, resonance, and chaos via numerous nonlinear dynamic theory techniques facilitated by intuitive system state augmentation and Lyapunov-Floquet transformation. Versal deformation of the normal forms is constructed to investigate the bifurcation behavior of the dynamical system. Outcome from the analysis indicates that the motion is quasi-periodic, chaotic, librational, and undergoing a Hopf bifurcation in the small neighborhood of the critical point-engendering locally stable limit cycles. Consequently, we demonstrate the implementation of linear and nonlinear control laws (i.e., bifurcation and sliding mode control laws) on the relatively acquiescent transformed attitude dynamics. By employing a two-pronged approach, the quasiperiodic planar motion is independently shown to be stabilizable via the nonlinear control approaches.
- gravity gradient
- nonlinear attitude control
- sliding mode
- Lyapunov-Floquet transformation
- normal forms
Ever since the launch of Sputnik—the first artificial satellite put into earth’s orbit in 1957—mankind has increasingly become dependent on space-based technology in many areas of our daily lives. For instance, space-based technology performs an indispensable role in telecommunications, navigation, personal entertainment, weather forecasting, farming, security, defense, scientific exploration, research, innovation, etc. Undoubtedly, the prominence of space technology in shaping humanity’s future is unequivocal. The success or failure of a given space mission is largely contingent upon the complex system analysis and design methodologies exerted in converting the initial idea into an elaborate functioning enterprise . It is for this reason that reliable and efficacious methodologies and tools are consistently utilized in space mission formulation and implementation. Thus, there is a need to continuously examine the effectiveness of prevailing space mission analysis and design methodologies. This is in order to improve prevailing tools and approaches that shall expedite relatively simpler, more reliable, and accurate mission modeling and analysis.
Space systems are required to function nominally in their designated orbital locations, maintain appropriate orientation, and conform to planned trajectories despite the ambient perturbing space environment. Strict mission pointing requirements normally constrain spacecraft in orbit around a large body to maintain a fixed stable orbital position and orientation during operation. However, perturbing space-environment torques act to dislodge positioned spacecraft and disorient stabilized ones [2, 3].
Modeling, analyzing, and controlling dynamics of space systems are therefore a crucial component of space mission design. The quest for relatively simpler, more accurate, and more reliable analytical methodologies and tools to represent, scrutinize, and manipulate the dynamics of space systems is therefore a worthwhile undertaking.
Inopportunely, dynamics of space systems tend to be commonly represented by coupled analytical models that possess complex structures encompassing nonlinearity, parameter-variant coefficients, and periodic external excitation terms [4, 5, 6, 7, 8, 9]. The requisite analysis essential to fathom such motion is not a trivial undertaking—except for few special cases, the general solution for such dynamical systems cannot be found. The complex structures of the motion’s analytical models characteristically point to nondeterministic and potentially chaotic systems over a range of initial conditions and system parameters. Therefore, to analyze dynamical space systems, we often have to be content with nonautonomous, nonlinear, and periodic differential equations [10, 11]. This presents an immense analysis challenge. For instance, time-varying eigenvalues of the periodic linear system matrix cannot determine the system stability. Consequently, methods such as linearization , averaging , and perturbation techniques [14, 15] have been consistently used to analyze such complex nonlinear, periodic motion. However, the two latter approaches tend to be limited to minimally excited systems (parameter multiplying the periodic terms is small), while linearization is restricted to small domains about the operating point. Further, such methods are inclined to be relatively cumbersome and normally augmented with numerical approaches to analyze dynamical systems and accomplish real-life applications [16, 17, 18].
The presence of perturbing torques in the ambient space environment tends to disorient an already stabilized spacecraft and further alters the orbital motion [2, 3]. This is contrary to the prevailing strict pointing requirements that constrain the satellite in orbit to maintain a stable orbital motion and fixed orientation during mission operation. A number of strategies are employed to stabilize spacecraft attitude motion and maintain a desired orientation despite the presence of perturbing torques in the space environment. The most common attitude control and stabilization approaches are three-axis stabilization, spin stabilization, and gravity-gradient stabilization. To provide the determined control input required to offset undesired attitude deviations, these methods employ either active control systems (e.g., thrusters, magnetic torquers, reaction wheels) or passive control systems (e.g., booms). Unlike passive controllers, active controllers utilize an external source of energy to drive the attitude control actuators [17, 19].
Among the stated attitude stabilization methods, gravity-gradient stabilization of spacecraft attitude is attractive due to its relatively intrinsic simplicity, reliability, and low cost . However, it is mostly feasible in low earth orbit due to its principle of operation as discussed in Section 2.
The motion about COM of a rigid gravity-gradient stabilized spacecraft is libratory about the pitch axis. This axis is normal to the orbital plane in an inverse-square gravity field. The satellite will oscillate about a position of stable relative equilibrium if the work done by external perturbing forces is greater than the rotational kinetic energy. The sufficient conditions for stability of relative equilibrium are explained in Section 2. The complete formulation of COM motion for a gravity-gradient stabilized satellite in eccentric orbit consists of six coupled, nonlinear second-order differential equations. This system of equations is considered analytically unsolvable in closed form [5, 7, 21, 22, 23].
This chapter aims to first investigate the periodicity, quasi-periodicity, and chaotic behavior of the gravity-gradient stabilized attitude motion. Moreover, the motion stability, resonances and bifurcation behavior will also be examined. Subsequently in Section 4, we synthesize suitable controllers to adequately offset the attitude perturbations experienced by the gravity-gradient stabilized spacecraft in an eccentric orbit. Requisite assumptions made to facilitate the attitude motion analyses will be explicitly stated and qualified.
To model, analyze, and control the nonlinear motion with parameter-variant coefficients and periodic forcing terms, we intend to use approaches based on:
System state augmentation
Lyapunov-Floquet (L-F) transformations
Normal form (NF) theory
The fitting use of the aforementioned transformations and techniques enables dynamical system analysis and control law development in transformed, parameter-invariant, and more tractable coordinates that preserve the original system Lyapunov stability properties [24, 25]. Consequently, we intend to exploit this propitious attribute in our investigation. Applications of L-F transformations in spacecraft dynamics have been previously investigated by authors such as [8, 26]. The former demonstrates how L-F theory enhances the representation of relative spacecraft dynamics in elliptical orbits, while the latter further proposes an orbit control law based on L-F theory. On the other hand, this chapter focuses on the dynamics of a rigid body about its COM while in elliptical orbit around a central large mass.
2. Gravity-gradient attitude stabilization in eccentric orbit
A gravity-gradient stabilized spacecraft attains a state of stable relative equilibrium when its points along the radius vector, points along the normal to the orbit plane, and is along the tangent to the orbit in the LVLH frame (Hill frame) as shown in Figure 1a. In addition, the condition must similarly be satisfied.
If the work done by external perturbing torques is greater than the rotational kinetic energy of the spacecraft about its COM, motion of the spacecraft in an elliptical orbit will be libratory as illustrated in Figure 1b. Equations representing the spacecraft orbital motion are identical with those of a point mass in an inverse-square law force field. To analyze the attitude dynamics, the spacecraft orbital motion (motion of COM) can be reasonably assumed to be independent of the spacecraft attitude motion (motion about COM). This assumption is justifiable because the satellite is small compared to the dimensions of the orbit. Under this assumption, the spacecraft’s orbital motion can hence transfer energy to the attitude motion, but the converse is assumed not to be possible. Thus, orbital parameters as determined functions of time are considered in analyzing attitude motion [7, 27].
When the spacecraft is considered as a rigid body in an inverse-square force field along an elliptical orbit, a complete formulation of equations of motion can be derived [5, 7, 21]. The resulting six second-order differential equations of motion are nonlinear and coupled. These equations of motion cannot be solved analytically in this exact form.
Ignoring other torques such as aerodynamic, magnetic, thermal bending, and solar radiation pressure, we can derive the equations of spacecraft attitude motion under the influence of inverse-square force field in an elliptical orbit. Additional assumptions are an ideal, perfect sphere earth without oblateness; largest spacecraft dimensions are extremely small compared to the orbit radius, and the spacecraft mass is negligible compared to the mass of the central body .
We further assume that the exact equations of motion can be linearized in small-angle motion characterization. Subsequently, the attitude dynamic models may be considered to consist of two equations with coupled roll-yaw angles and a third uncoupled equation describing the pitch angle dynamics. The pitch motion equation is hence independent of roll-yaw motion. The coupled roll-yaw -equations are homogeneous and can be solved for [5, 7, 21, 27].
Consequently, the exact problem is reduced to the equation of pitch motion with orbital parameters as functions of time and spacecraft mass parameters shown in Eq. (1):
where is a dimensionless ratio of the spacecraft’s principal moments of inertia given by
To analyze attitude motion in eccentric orbit, we substitute time with true anomaly, where f is the independent variable. Moreover, the COM will obey the following Keplerian orbit relations:
This is the well-known equation of plane pitch angle libratory motion in elliptical orbit [5, 7, 21, 27]. The primes indicate differentiation with respect to f. The planar pitch attitude motion equation is hence nonlinear with periodic coefficients in f. Analysis of this motion and subsequent synthesis of a fitting controller is not a trivial task. We hence intend to analyze this motion and synthesize suitable controllers to stabilize the system.
3. Attitude motion analysis
In general, L-F transformation techniques facilitate obtaining solutions of dynamical systems with periodic coefficients, evaluate periodically forced responses, and design feedback control laws. We shall augment these capabilities with normal form techniques (simplifies nonlinearity) and state augmentation (converts nonautonomous to autonomous system). Subsequently, the emanating synergies serve to accomplish the objectives of this chapter.
Floquet theory enables stability and response analysis of linear systems with periodic coefficients, i.e., . Floquet theory utilizes knowledge of characteristic exponents of the state transition matrix (STM) to infer that if the solution of a system is obtained over a full principle period, then the solution is known for all time [28, 29].
Further, the the Lyapunov-Floquet (L-F) transformation reduces the original nonautonomous linear differential system to an autonomous one with the form where is an constant matrix. In the LFT matrix, can be approximated via the methodology described by [25, 30] using shifted Chebyshev polynomials of the first kind. Chebyshev polynomials are chosen because they produce better approximation and convergence than other special orthogonal polynomials . The approximation of LFT matrix via Chebyshev polynomials contains elements with truncated Fourier representations as shown in Eq. (8). has similar Fourier representation:
Time-independent normal form (TINF) simplification facilitates construction of relatively lesser complex but qualitatively equivalent models of the original nonlinear dynamical systems. On the other hand, time-dependent normal form (TDNF) simplification considers original nonlinear dynamical systems with periodic coefficients by utilizing Lyapunov-Floquet transformation. This simplification is generally implemented on equations arising from Taylor series expansion via successive application of nonlinear near-identity transformations. Such a transformation entails preservation of the original system’s stability and bifurcation characteristics by the transformed models. The fundamental concept behind normal forms methodology is to simplify the system by eliminating as many nonlinear terms as possible. This is accomplished via application of successive series of near-identity transformations on the original system. The near-identity coordinate transformations are nonlinear and local. The reader is directed to the well-documented literature on normal form theory found in works by authors such as [32, 33, 34, 35].
To normalize nonlinear systems subjected to external periodic excitation, several authors such as [34, 35, 36, 37] either utilize approaches that introduce equation variables and/or detuning parameters or incorporate a bookkeeping parameter in their methodology. However, the augmenting parameters involved seemingly lack a uniform explicit connection to the terms in the dynamic equations under consideration. Consequently, we shall utilize a relatively more straightforward and intuitive approach that involves augmenting the system states by converting the periodic external excitation into a system state. The state augmentation approach here is intuitive because the augmented states directly emanate from the periodic forcing term(s); hence, they are neither ad hoc nor arbitrary. Moreover, neither a detuning parameter nor a bookkeeping parameter is required.
Strictly speaking the behavior of the attitude motion as shown in Eq. (7) is characterized in terms variation in true anomaly, f. However, the true anomaly similarly varies with time; hence, the reference as implicit time history is preferred. To demonstrate implicit time history behavior, we initially select and . Later on in Section 3.3, we shall analyze the impact of different -pair values on the motion.
Both the motion in original coordinates, Eq. (7) and the corresponding state augmented system, will be scrutinized. Figure 2 shows the implicit time history behavior of the motion in original coordinates.
Similarly, we shall scrutinize the history behavior of the system in Eq. (7) after augmenting the system states. In accordance with binomial expansion theorem , since and , then the magnitude of the terms in the binomial series progressively become smaller. Therefore, the binomial expansion of the term can be approximated as . Eq. (7) becomes
Therefore, the augmented system state-space representation with is shown in Eq. (12):
Figure 3 shows the augmented state system history behavior. Allowing for the expected minor discrepancies due to series and binomial expansion approximations, the system state response is comparable to that of the system in original coordinates shown in Figure 2.
From Figures 2 and 3aandb, the attitude motion is quasiperiodic as characterized by the absence of closed trajectory attractors in the phase space. The pitch angle librates roughly between , while the pitch angle rate of change varies between . The orbits generally appear to follow a “heart-shaped” path starting at the origin with two non-closing lobes on either side. Conversely, the augmented states are periodic as characterized by the closed circular limit cycle attractor centered at the origin in Figure 3d. The motion behavior discussed here indicates that and are susceptible to instability, unpredictability, and chaos. Consequently, these aspects of the motion flow will be investigated next. Neither the eigenvalues of the linear periodic term of Eq. (7) nor of the state augmented system Eq. (12) can be used to determine stability. Consequently, we’ll have to construct the Floquet transformation matrix (FTM) to analyze the dynamical system’s stability.
3.1 Stability and chaos
Stability analysis is preceded by computation of the dynamical system’s state transition matrix, . We first prepare the motion in Eq. (7) for expansion via shifted Chebyshev polynomials of the first kind by normalizing the principal period. This is because shifted Chebyshev polynomials are only valid for the period interval .
where , also series coefficient , as
We can ignore the terms of order greater than 7 in Eq. (13) without significant loss of accuracy because the follow-on terms have relatively small successive coefficients that rapidly approach zero. For instance, the 9th order term has the coefficient , while the 11th order term’s coefficient is Substituting the expanded trigonometric product in Eq. (7), we obtain
To normalize the principal period, let . Equivalently, represents duration within the principal period. Let represent the principal period, and then this implies that for a periodic term, . It follows
After substituting and from Eq. (15) into Eq. (14), the obtained state-space representation of the normalized attitude motion is given in Eq. (16), where and . Further, constitutes the initial conditions of the represented motion which correspond to pitch librations at position 1 of Figure 1b:
It is clear that Eq. (16) is of the form
3.1.1 Floquet multipliers and exponents
To facilitate computation of STM, FTM, and L-F transformation matrices using Chebyshev polynomials, we utilized the Chebfun software package on MATLAB™ . Summarily, Chebfun applies piecewise Chebyshev polynomial interpolation to construct smooth functions over the interval . Recall that and are the normalized principal period; hence, the computed is
The computed Floquet multipliers are critical since they lie on the unit circle with values of as shown in Figure 4. Consequently, this reveals a marginally stable system for the chosen -pair. It follows that the corresponding Floquet exponents are purely imaginary. This is consistent with the quasiperiodic system phase portraits that illustrate a librational motion.
The motion is hence stable in the sense of Lyapunov, but the inherent oscillations are disruptively significant to jeopardize nominal execution of the spacecraft mission.
3.1.2 Poincaré map
Figure 5 shows the constructed Poincaré section of the flow. There is a discernible main cluster of points in close proximity to the origin but restricted to the positive side of . Relatively scanty, isolated discrete points occupy the lower bottom half of the plot bound by
The Poincaré section composition suggests two possible flow behaviors. The groupings suggest a quasiperiodic trajectory. On the other hand, the scanty random points devoid of clustering could be due to transient behavior or chaos. A chaotic motion can briefly dwell on a near periodic trajectory before changing to a disparate trajectory with a period that is a multiple of the preceding motion. Consequently, it is needed to further investigate the presence of chaos in the attitude motion dynamics.
We define chaos as a bounded aperiodic steady-state motion behavior that is not in equilibrium and is sensitive to initial conditions. A minuscule divergence in the input rapidly grows to spawn an overwhelming difference in the system response. We begin investigating chaos in the attitude motion by plotting the system implicit time history with a minute divergence to the initial condition of the state, . We set the first initial condition to zero; the divergent second initial condition is obtained by adding to the zero initial condition. The obtained implicit time history of the two curves reveal the onset and progression of an overwhelming difference in the system response as illustrated in Figure 6. The slightly divergent initial condition results in an overwhelming difference in response that begins in the second half of the 4th orbit and then rapidly grows in the subsequent orbits. The attitude motion is hence chaotic.
3.1.4 Lyapunov exponents
To determine the average rate of divergence between the initially neighboring trajectories defined locally in the state space, we shall scrutinize the dynamic behavior of the motion’s Lyapunov exponents. Lyapunov exponent stability analysis affords a means of quantifiable expression for initial conditions sensitivity and dependence (i.e., chaos), by describing the exponential rate of growth or decay of a perturbation to a trajectory as time progresses . Lyapunov exponent is numerically expressed as
where is the tiny separating perturbation vector between the trajectories. The value of Lyapunov exponent will distinguish the nature of the trajectory according to the following criteria: (i) : Trajectory is stable and the motion is asymptotically stable. (ii) : Trajectory is neutral and the motion is characterized by some sort of steady-state. (iii) : Trajectory is unstable and chaotic. The Lyapunov exponent behavior for the motion given Eq. (16) is illustrated in Figure 7.
The computed Lyapunov exponents are equal in magnitude but opposite in sign with increasing periods because the flow in Eq. (16) is nonautonomous Hamiltonian. Hamiltonian systems are conservative. Therefore, the magnitude of which measures expansion in one direction is equal to the magnitude of which measures contraction in another direction. Since always, then as prescribed by the above distinction for , the attitude motion is chaotic. This outcome is consistent with the preceding chaos analyses that scrutinized motion sensitivity to initial conditions.
3.1.5 Stability charts
The orbit eccentricity and spacecraft’s ratios of principal moments of inertia are, respectively, defined as and . On the other hand, our stability analysis so far has been illustrated in the courtesy of arbitrarily set values of . Consequently is it essential to holistically scrutinize the motion behavior for all possible values of and . Constructing stability charts which partition the -plane into stable and unstable regions enables scrutiny of motion stability as and vary simultaneously. Transition curves in stability charts constitute frontiers that separate stable regions from unstable regions. We can derive transition curves in closed form via the FTM. Floquet theory requires a stable system to have a Floquet multiplier of magnitude, . It can hence be proven that the transition from stability to instability occurs when both Floquet multipliers are equal to 1 or both are equal to −1 (see Figure 4). Therefore, the transition curves in the -plane where the solution to the linear periodic term of Eq. (17) changes from stable to unstable (or vice versa) are determined by the conditions or .
Even though we can only construct transition curves associated with the linear term of the attitude motion equation, the outcome provides insightful perusal into the behavior of the whole equilibrium solution. The complete solution behavior can be arrived at by augmenting the evaluated linear periodic stability behavior with the combined effect of the nonlinear term and forcing term in Eq. (17). For instance, if a given -pair is initially located in the unstable region of the FTM-dependent stability chart, then the nonlinear and forcing terms will tend to exacerbate this instability, rendering the complete solution unstable for that particular combination. A similar argument can be made for an -pair located in a stable region. The constructed stability chart of the attitude motion is shown in Figure 8.
The darker regions constitute stable points while the lighter regions are unstable. The stability chart has a slightly larger stable region than the unstable region. Unstable solutions appear to be dominated by two regions approximately defined by (i) and increasing values of e plus (ii) and increasing values of . Essentially, a spacecraft with mass distribution such that is amenable to a wider range of eccentricity values above 0.2 to achieve an intended stable pitch angle motion.
For instance, we have previously considered the pair ; this pair is located in the stable region corresponding to critical Floquet multipliers confirmed in Figure 4. The motion at this location is stable in the sense of Lyapunov. The stability chart further accords a means of scrutinizing generalized stability behavior trends or commonalities between disparate -pairs. By picking representative -pairs from different regions, we tabulate the resultant illustrative Floquet multipliers as shown in Table 1.
|Stability region||e||σ||Floquet multipliers (|
From Table 1, we note that both selected marginal and stable regions are associated with critical Floquet multipliers. This implies that the pertinent -pair characterizes a motion stable in the sense of Lyapunov; i.e., the pitch angle librations are bound by . This is consistent with the dynamics presented in Figures 1 and 3. However, in the unstable regions, Floquet multipliers have magnitudes , implying that the pitch angle wanders beyond
The attitude dynamics are dominated by the linear and forcing terms delineated in Eq. (17). This is because if we consider the motion composed of only the linear and forcing terms, i.e., , the numerical solution is unbounded as shown in Figure 9. This is not the case if we consider motion composed of any of the following term combinations: , , or .
Moreover, if we consider the L-F transformed motion, Floquet exponents (eigenvalues of ) represent frequencies associated with the linear term. Periodic elements of the nonlinear matrix are a product of the truncated Fourier series matrices and . On the other hand, periodic elements of the forcing matrix are likewise multiplied by . Subsequently, resonance between the Floquet exponents and any of the periodic terms in the forcing matrix elements will trigger instability in the motion as well.
Bifurcation triggers the system’s equilibrium solutions to transition between the disparate regions of the stability chart. The orbit eccentricity, e, is the bifurcation parameter for the attitude motion (see Eq. (16)). This is because the general on-orbit spacecraft mass properties represented by typically tend to be constant for a gravity-gradient stabilized spacecraft. Therefore, it is essential to analyze the equilibrium solution dynamics as small increments are applied on the bifurcation parameter. We develop the normal form of our dynamics in the next section to facilitate bifurcation behavior analysis. Normal forms are not unique, consequently near identity transformation for the state-augmented system, and the L-F transformed system will be undertaken separately.
3.3 Versal deformation of the normal form and bifurcation analysis
Versal deformation refers to embedding the system in a parameterized family of systems containing all possible dynamics that can occur near the bifurcation point. Moreover, the family of systems should be transverse to the bifurcation surface with the number of parameters equal to the codimension of the bifurcation . The attitude motion undergoes a codimension-one bifurcation because only one parameter e is responsible for the loss of stability (for gravity-gradient stabilization to be maintained, has to remain fixed). Because our critical Floquet multipliers are complex and lie on the unit circle, this system will experience a Hopf bifurcation. Further, it is “well-known” that Hopf bifurcation is a codimension-one bifurcation. Firstly, we shall formulate the normal forms of nonlinearities up to the cubic order in Eqs. (16) and (12) to demonstrate the intended approach. Normalization of dynamics with higher order nonlinearities can be accomplished via the same techniques.
3.3.1 State-augmented system
This system now possesses 4th order nonlinearity. is in the Jordan canonical form. The normal form is evaluated by successive application of near identity transformation of the form
where is an homogeneous vector of monomials in of degree r. The state-augmented system is independent of periodic coefficients; hence, we solely obtain the TINF given in Eq. (22):
When the external forcing term is augmented as a system state, the magnitude of the external forcing frequency appears as solitary, linear, imaginary conjugate coefficients in the normal form, i.e., 1 (see Eq. (9)). Eigenvalues of the linear matrix in Eq. (12) constitute the conjugate coefficients in the linear terms of the reduced nonlinearity normal forms (i.e., ).
Moreover, after obtaining the straightforward solutions of and and then substituting these values in the equations for and , the last two equations of the normal form can be expressed as
where C1 and C2 are the integration constants originating from the analytical solutions of and , respectively.
We shall investigate the bifurcation of Eq. (12) via its normal form given above. Because the periodic system maintains the same general structure, we may treat the respective limit cycles as equilibria and study their bifurcations. We utilize the versal deformation of the normal form to investigate the change in the stability structure of the dynamics in the neighborhood of the critical point of the bifurcation parameter. Essentially, construction of versal deformation of the normal form facilitates characterization of system dynamics at the critical point and its small neighborhood. Therefore, we handily gain complete understanding of the qualitative phase space dynamics of the dynamical system in the neighborhood of the critical point.
We define the normal form versal deformation parameter as . The parameter represents a small change in the eigenvalues of the normal form corresponding to a small change in the bifurcation parameter in the original system coordinates. It is a prerequisite condition to obtain a relationship between the versal deformation parameter and the original system bifurcation parameter . Incorporating in Eq. (22), we obtain Eq. (23):
By defining small increments on the bifurcation parameter as , we can write to represent the disparate sets of bifurcation parameter in the neighborhood of the critical parameter . We employ the least squares, curve fitting technique proposed by  to obtain the relationship between and , as . NB values of and in Eq. (23) were evaluated by forward action transformations of initial conditions in the original coordinates, i.e.,
The closed-form analytical solutions for and in the versal deformation in Eq. (23) are straightforward. To obtain and , we introduce the complex changes of variable and followed by the polar coordinates and . The last two equations in Eq. (23) become
Similarly, and are the integration constants originating from the analytical solutions of and , respectively. The values of these integration constants are evaluated from the initial conditions specified in the original coordinates. After back transformation of the normal form closed-form analytical solutions above, we obtain the motion in the original coordinates.
The back transformed and constitute the augmented states given in Figure 3c. Moreover, but is generally small in the order of magnitude . The integration constants on evaluation are imaginary whose magnitudes are close to identity. Consequently, from Eq. (25), back transformation of the sinusoidal and will result in the trigonometric augmented states whose amplitude is determined by the magnitude of the integration constants. From Eq. (24), we can express the transient solution of as . Since , the motion is characterized by a locally stable limit cycle in the neighborhood of the bifurcation point. The limit cycle is stable in the sense of Lyapunov but not asymptotically stable. Post-bifurcation attractors that transform into quasiperiodic attractors portraying a limit cycle in the original coordinates are obtained via back transformation and subsequently shown in Figure 10 ().
3.3.2 L-F transformed system
Here, we also demonstrate analysis of bifurcation behavior subject to different values of -pair. Consequently, by utilizing the developed stability chart (Figure 8), we select . This -pair lies on a transition curve; hence, the system is guaranteed to be bifurcating. Again, by considering up to the cubic nonlinearity, the history behavior from Eq. (16) is shown in Figure 11.
The system likewise possesses critical Floquet multipliers that lie on the unit circle of values () and purely imaginary Floquet exponents, . The computed FTM and matrices are given in Eq. (26):
The computed periodic L-F transformation matrix and are plotted in Figure 12.
We consider the L-F transformed dynamics in Eq. (27) up to the cubic nonlinearity. After applying the L-F transformation given in Eq. (19) on the attitude motion in Eq. (16), we obtain the system in Eq. (27):
To normalize this externally excited motion, the system states are augmented to convert the system from nonautonomous to autonomous. We define additional states shown in Eq. (28) and plotted in Figure 13:
After substituting the above augmented states into Eq. (27), we obtain the system shown in Eq. (29)—whose order of nonlinearity increases to four. The transformed denominator term has been approximated by the binomial expansion equivalent, i.e.,
Apart from raising the order of nonlinearity, state augmentation further introduces periodic linear terms, and in Eq. (29). Consequently, the augmented dynamics with linear parameter-variant coefficients necessitate a second L-F transformation to convert the linear term to parameter invariant:
The computed parameters corresponding to the second L-F transformation are as follows: critical Floquet multipliers, , and purely imaginary Floquet exponents, . The computed second FTM and constant matrices are given in Eqs. (30) and (31), respectively:
The computed second periodic L-F transformation matrix, , and its inverse, , are similarly presented as plots in Figure 14.
We designate the second L-F transformation as . Here, . After applying the second L-F transformation to the state augmented periodic system, we obtain
We first evaluate the TDNF and then average out the periodic terms to obtain the simplified TINF. The closed-form analytical solutions for and are constants. Variables and in the and differential equations are substituted by their respective computed constants. This computation is carried out through the forward action transform of the L-F, modal, and near-identity transformations of the initial conditions declared in the original coordinates:
The Floquet exponents are conjugate coefficients in the linear terms of the normal forms before being multiplied by the substituted constant values equal to and .
Similarly, we define the normal form versal deformation parameter as and incorporate it into Eq. (34):
After defining small increments on the bifurcation parameter again as , we express the disparate sets of bifurcation parameter in the neighborhood of the critical parameter as, . The relationship between and is evaluated via the procedure previously stated in Section 3.3.1 to yield .
The closed-form analytical solutions for and of the versal deformation normal form are straightforward. To obtain and , we introduce the complex changes of variable, and followed by the polar coordinates and . The last two equations in (35) become
() are the respective constants of integration whose value is evaluated from the initial conditions specified in the original coordinates. and are real, whereas and are complex. is similarly small in the order of magnitude . and back transformation via inverse near-identity, modal, and single L-F transformations forms the augmented states given in Eq. (28) and is plotted in Figure 13. Eq. (36) yields a steady-state solution of the limit cycle amplitude as . When , the solution of and results in locally stable limit cycle with amplitude corresponding to . Consequently, the quasiperiodic attractors in the original coordinates delineating a limit cycle are obtained after back transformation as shown in Figure 15 ().
Solutions of the versal deformation equations enable investigation of the post-bifurcation steady-state behavior in the small neighborhood of the bifurcation point. However, as observed by , this method is only useful for local analysis. This is because minor errors introduced by back transformation close to the bifurcation points significantly grow as you move further away.
4. Attitude motion feedback control
After setting and , the motion in Eqs. (9) and (16) is then numerically integrated to obtain the uncontrolled responses shown in Figure 16a and b, respectively (NB Figure 16a is the same as Figure 2a). The slight difference in the long-term motion behavior in Figure 16b may be attributed to the approximation of the trigonometric product term by a truncated series in Eq. (13) and possibly fidelity of the numerical integrator used in the computation. As established in Section 3.1, the attitude motion for the considered -pair is quasiperiodic, marginally stable, and chaotic. Despite the system being stable in the sense of Lyapunov, the inherent oscillations are disruptively significant and require stabilization if the spacecraft is expected to successfully conduct its mission.
Motion controller design shall be undertaken on augmented state system, L-F transformed, and near-identity transformed coordinates. We shall hence first transform the system dynamics into these more amenable but topologically equivalent dynamical structures that retain the Lyapunov stability and bifurcation properties of the original system. Augmentation of the attitude dynamic states has been conducted in Section 3. Control law development will be first considered in L-F transformed coordinates and then followed by near-identity transformed coordinates of Eq. (16).
4.1 Lyapunov-Floquet transformation
Prior to computing the L-F transformation matrix and its inverse, matrices (for the case), we computed the FTM and matrices over the interval via shifted Chebyshev polynomials of the first kind for the system in Eq. (16). The evaluated aforementioned matrices are shown below. Here we alternatively present elements of and as truncated Fourier series described in Eq. (8). Previously in Figures 12 and 14, we have presented the periodic plots of these series. Further, recall that and are the normalized principal period; hence, :
After applying the L-F transformation, , to the attitude motion in Eq. (16), it becomes
The Lyapunov stability properties are preserved in the new coordinates after the system is transformed by the L-F transformation matrix. The L-F transformation theory guarantees that a suitable controller realized in the L-F transformed coordinates will be correspondingly efficacious after back transformation into the original system coordinates. Consequently, we shall endeavor to systematically synthesize suitable controllers to stabilize the motion in the transformed coordinates. Our control synthesis strategy will first consider linear control laws before exploring nonlinear control strategies.
In order to formulate appropriate control laws that would stabilize the quasiperiodic motion analyzed in Section 3, we introduce a control input in Eq. (1) as shown below:
The control action will generally represent torque per unit moment of inertia as a function of true anomaly. Eq. (42) is first used to synthesize linear control laws followed by nonlinear control law development.
4.2 Linear control
Though linear control law principles are conventionally intended for controlling linear parameter invariant systems , we initially consider them to control our nonlinear dynamics as an initial analysis step. Since most linear control methods tend to be relatively simpler to analyze and implement compared to nonlinear control methods, it is prudent to ascertain the suitability of linear control prior to embarking on relatively more complicated techniques. To implement linear control, we shall consider pole-placement approach to determine the negative feedback gain required to stabilize the system.
4.2.1 State-augmented system
where constitutes the nonlinear terms, is the linear matrix, and is the state vector. To synthesize the parameter-invariant linear state feedback controller, Eq. (42) becomes
The linear state feedback controller is of the form , and the control input scaling vector is . Though is full rank, the linear pair is not controllable. This is because the system controllability matrix, , shown in Eq. (45) has a rank of 2, instead of rank 4:
Consequently, a linear state feedback controller cannot stabilize the system dynamics associated with the -pair considered. It is to be noted that the two states, in Eq. (44), are virtual states serving to simplify the system but are not accessible in the actual system dynamics. This is illustrated by the fact that controllability matrix does not have full rank.
4.2.2 L-F transformed system
In this case, the parameter-invariant linear state feedback controller is similarly of the form . The control input is scaled by the matrix in L-F transformed coordinates. Back transformation of the product via inverse L-F transformation will guarantee a single control input in the system original coordinates as will be demonstrated in Eqs. (54) and (55). is full rank and the linear pair is controllable. The L-F transformed Eq. (42) will be shown in Eq. (46):
Therefore, the system closed-loop dynamics subjected to a linear control law will be of the form
We initially place poles at . Then we evaluate the corresponding matrix to realize this pole placement. Several stable double pole locations with a decreasing factor of 10 (i.e., −1, −10, −100, −1000…) were considered. None of these pole-placement locations demonstrated notable success in stabilizing the system. For instance, poles placed at produce a response for a duration of slightly beyond 1.5 cycles before the states abruptly become indeterminate at about as shown in Figure 17.
In this analysis, the system response in the original coordinates is realized via the back transformation . Therefore, similar to the state augmentation case, the L-F transformed nonlinear system demonstrates inability to be stabilized by a linear control law. The presence of periodic coefficients (elements of ) associated with the nonlinear and forcing terms renders the system untenable to be controlled via LTI system control approaches. In general, the “region of application” of linear control for nonlinear systems is dependent on magnitude of nonlinearity and initial conditions. Many times, linear control may stabilize nonlinear systems locally, but this is not guaranteed.
4.3 Nonlinear control
Nonlinear control appears more suitable than linear control to stabilize the attitude motion. However, conventional nonlinear techniques are not readily amenable to dynamics with periodic coefficients and periodic external excitation. Hence, requisite system state augmentation, L-F, or near-identity transformations will be undertaken prior to controller design. We shall first consider sliding mode control (SMC), and then bifurcation control will be implemented on the marginally stable system to stabilize post-bifurcation response.
4.3.1 Sliding mode control
Sliding mode control is a robust nonlinear feedback control methodology that is suitable for achieving accurate tracking for a class of nonlinear systems. SMC methodology is based on variable structure control law that results in the state trajectory “sliding” along a discontinuity surface in the state space [44, 45]. Though SMC is deterministic, nonlinear, and robust, its implementation is prone to undesirable “chattering” along the sliding surface . Design of SMC involves (i) selection of the switching function (stable hyperplane in the state space on which the dynamics will be restricted) and (ii) control law synthesis.
22.214.171.124 State-augmented system
Here, we implement a SMC that tracks a desired null pitch angle via a negative rate of growth. Dynamics in the original coordinates possess periodic coefficients rendering the dynamics unwieldy and unfavorable to synthesize a SMC. Therefore, we develop the SMC law based on the augmented state dynamics—which are liberated from periodic coefficients. To design a sliding mode controller for the state augmented systems, we designate the switching function as given in Eq. (48). The switching function represents the actual system state (i.e., attitude pitch angle) reference error (difference between desired and actual pitch angle) that the controller desires to maintain at zero. Therefore, when , as :
where represents the control input. Derivative of the switching function with and substituted from Eq. (49) is
Setting , we derive the following control input:
A sigmoid function, is preferred instead of the signum function to reduce chattering around the sliding surface typical of sliding mode controllers. is generally small. Employing a direct Lyapunov approach, stability of the sliding mode controller applied here is ascertained by setting to be the Lyapunov function. Hence, . The guaranteed negative definiteness of the Lyapunov function derivative demonstrated in Eq. (53) points to a stable controller. This equation is obtained after subsequent substitution for and in the equation of :
Figure 18 shows the sliding mode controlled system response.
It is observable that the sliding mode controller in the state-augmented system achieves stabilization of the motion throughout any number of orbits. Both and are adequately confined to zero. The augmented states and remain unaffected.
126.96.36.199 L-F transformed system
where . We define the sliding function according to Eq. (56) to ensure that when , as . The sliding surface represents the reference pitch angle error. The controller attempts to maintain a zero error throughout, i.e., :
From Eq. (57), we set and derive the following sliding mode control law:
To reduce chattering around the sliding surface typical of sliding mode controllers due to fast switching of the signum function, a sigmoid function is similarly preferred. We again apply the direct Lyapunov approach to analyze the sliding mode controller stability by selecting as the Lyapunov function. Asymptotic stability will be guaranteed if the sliding function derivative is negative definite. Hence, the switching function derivative is . Substituting for with the control input likewise substituted, we obtain the stability-criteria satisfying relationship below:
Figure 19 shows the sliding mode controlled system in L-F transformed coordinates. The response in Figure 19 is back transformed via the inverse L-F transformation resulting in controlled states in the original coordinates shown in Figure 20.
We observe that, similar to the state augmented case, the sliding mode controller stabilizes the L-F transformed motion as well by invariably confining the states to zero as desired. Though specific values for the -were used to demonstrate this technique, stabilization of the planar pitch motion by SMC approach is independent of the assigned -values. However, the possibility of a synthesized sliding mode controller being impractical to implement exists if the required control effort is colossally prohibitive.
4.3.2 Bifurcation control
The critical Floquet multipliers corresponding to purely imaginary Floquet exponents indicate that the system is in the stability boundary. Consequently, it is essential to stabilize the system post bifurcation apart from modifying other motion characteristics such as rate of growth. To achieve these objectives, we engage nonlinear bifurcation control with full state feedback. Synthesis of such a controller is facilitated by the normalized dynamics which are relatively more tractable compared to the dynamics as represented in the original coordinates. Periodic coefficients and complexity in structure of the dynamic equations in the original coordinates drastically convolute synthesis of bifurcation control law. Dynamics of the states in the original coordinates will eventually be obtained via back transformation of the normal form, modal, and L-F transformations. Location of the complex Floquet multipliers on the unit circle (Figure 4) indicates that the pitch attitude motion is undergoing a Hopf bifurcation with a limit cycle attractor of controllable radius. Therefore, the structure of the normal form will also verify a Hopf bifurcation occurring in the neighborhood of the critical point of the bifurcation parameter (i.e., orbit eccentricity).
To illustrate the intended approach, we shall formulate the normal forms of nonlinearities up to the cubic order in Eqs. (12) and (16). Normalization of dynamics with higher order nonlinearities can be accomplished through the same techniques. Similar to the preceding cases, we’ll consider the augmented states and L-F transformed systems separately.
188.8.131.52 State-augmented system
In Section 3.1.1, we demonstrated how to obtain the TINF of the state-augmented system—shown in Eq. (23). Obtaining the closed-form analytical solutions for and in Eq. (23) is straightforward. On the other hand to evaluate and , we introduce the complex changes of variable, and followed by the polar coordinates and . The last two equations in Eq. (23) become
where C1 and C2 are the integration constants obtained when solving for and , respectively.
We solve Eq. (61) and then utilize the results to complete the closed-form analytical solution of Eq. (23). The closed-form solutions of and are then back transformed to the original coordinates producing the uncontrolled motion behavior shown in Figure 20. The system response in Figure 21 is a cognate approximation of the originally obtained numerical solution in Figures 3a and 16. Again, a quasiperiodic motion is characterized by non-closed curves and is observed in the corresponding phase portrait. Moreover, a codimension-one Hopf Bifurcation is verified by the normal form structure.
To synthesize a bifurcation control law of the normal form, we first add a control input in Eq. (62):
Let the scaling matrix and control input, respectively, be of the form
Back transformation of the product via inverse normal form and modal transformations will guarantee a single control input in the system original coordinates as demonstrated in Eqs. (54) and (55). The proportional gains are custom tuned to and . is a scalable parameter meant to suppress strange trajectory behavior according to Poincaré-Bendixson theorem in the system phase space. The resulting response of the bifurcation-controlled augmented state system is shown in Figure 22. The augmented states remain unaffected as previously shown in Figure 18b.
The libratory amplitude of the quasiperiodic pitch angle motion is tremendously stabilized and confined to a significantly diminished limit cycle attractor as illustrated in Figure 21.
184.108.40.206 L-F transformed system
As already indicated in Section 3.3.2, in addition to synthesizing bifurcation control law via L-F transformed dynamics, we shall also demonstrate analysis of the spacecraft attitude dynamics due to different values of e and . Therefore, and is once again considered in this section. L-F transformation analysis of the attitude dynamics associated with these values of e and has been comprehensively conducted in Section 3.3.2. Subsequently, the normalized TINF system was obtained in Eq. (34).
In Eq. (34), the closed-form analytical solutions for and are constants. Variables and in the and differential equations are substituted by their respective computed constants. This computation is carried out through the forward action transform of the L-F, modal, and near-identity transformations of the initial conditions declared in the original coordinates.
The Floquet exponents are conjugate coefficients in the linear terms of the normal forms before being multiplied by the substituted constant values equal to and .
To obtain and , we introduce the complex changes of variable and followed by the polar coordinates and . The last two equations in Eq. (34) become
Results from Eq. (64) which is easier to solve are then used to obtain the closed-form analytical solutions to Eq. (34). Then, are then back transformed to the original coordinates, producing the uncontrolled motion shown in Figure 23. The system response in Figure 23 (with nonzero initial conditions) is a cognate approximation of the originally obtained numerical integration solution in Figure 11. The back transformed augmented states are similarly shown in Figure 24 corresponding to Eq. (28) where the amplitude of is times that of p.
Let the scaling matrix and control input be of the form
The proportional gains are custom-tuned to and . Figure 25 shows dynamic behavior of the implemented bifurcation control in original coordinates with nonzero initial conditions.
The oscillating amplitude of the quasiperiodic pitch angle motion is tremendously stabilized relative to the initial behavior illustrated in Figure 23. This hence demonstrates successful control of the post-bifurcation attitude dynamics about the spacecraft center of mass.
Bifurcation control is a nonlinear control technique that affects the behavior of the closed-loop system by modifying nonlinearity and post-bifurcation behavior. Therefore, the location of Floquet multipliers (exponents) is generally preserved post-bifurcation control. Figure 26 shows the preserved locations of the Floquet multipliers before and after bifurcation control (and ). This location of Floquet multipliers is consistent with the limit cycle shown in Figure 25b corresponding to a simply stable system with relatively subdued librations.
In this chapter, we illustrated techniques for analyzing and stabilizing the attitude motion of a gravity-gradient stabilized spacecraft. The motion dynamics are shown to be nonlinear with periodic coefficients and subjected to external periodic excitation. Methodologies employed here utilize state augmentation, Lyapunov-Floquet transformation theory, and normal forms to realize relatively more tractable dynamical systems that are amenable to conventional controller synthesis techniques. Floquet theory was used to investigate system stability. State augmentation facilitated analysis via normal forms by transforming the dynamical system from nonautonomous to an autonomous one.
Outcome from the analysis showed that the attitude motion is quasiperiodic, chaotic, and stable in the sense of Lyapunov for the particular -pairs considered. Subsequently, the motion stability chart that was constructed facilitated prediction of -combination leading to stable or unstable dynamics. The stable regions of the stability curves were found to predict marginal and not asymptotically stable dynamics. However, the emanating librations need to be stabilized for nominal mission operations to be realized. Conversely, the -combinations located in the unstable regions resulted in aperiodic unstable dynamics. The computed Lyapunov exponents indicate that the chaotic dynamics also depend on initial values of pair, not just on the magnitudes of -pairs.
Both outcomes of the twofold versal deformation analyses (disparate values of -pairs considered) indicate establishment of locally stable limit cycles by the quasiperiodic flow post bifurcation. Since the eccentricity varies as , relatively small deviations from the critical point of the order trigger a significant topological change in the structure of the motion flow.
The quasiperiodic, nonlinear, and periodically forced pitch attitude motion is challenging to control. The synthesized linear controller served as starting point for developing more adept control laws. Not surprisingly, the linear controller failed to stabilize the complexly structured nonlinear dynamical system. As stated, in general the “region of application” of linear control for nonlinear systems is dependent on magnitude of nonlinearity and initial conditions. Many times, linear control may stabilize nonlinear systems locally, but this is not guaranteed.
On-orbit perturbations cause disturbing torques that bifurcate the attitude motion; it is hence imperative to stabilize the system attitude dynamics in the small neighborhood of the bifurcation parameter’s critical point. Local nonlinear bifurcation control law implemented on the attitude motion undergoing a Hopf bifurcation was shown to stabilize the attitude motion. The bifurcation controller which modifies the nonlinearity and post-bifurcation behavior further prevents the attitude motion from becoming chaotic because bifurcation is the path to chaos. Implemented in the TINF, the bifurcation control law would subsequently stabilize the secular and periodic attitude perturbations experienced by a spacecraft in elliptical orbit about its nominal operating point.
Sliding mode control law was based on driving both system states to zero on the sliding surface when the sliding surface reference error is equal to zero. The SMC law was similarly shown to be successful by invariably restricting the pitch angle to zero.
Future work will consider torques generated by sources such as magnetism and oblateness of the earth, atmospheric drag, solar radiation pressure, thermal bending, etc. Further, nonlinearities beyond the cubic term in the L-F transformation and TDNF case of near identity transformation would also be analyzed. In addition, physical implementation of the controllers and derivation of TDNF-based control laws requires future scrutiny. As demonstrated, all the control effort inputs are single torques per unit moment of inertia which for instance can be implemented via thrusters. Consequently, sizing and implementation of the control effort are an essential task.
Financial support for this research was partially provided by the Interplanetary Initiative of Arizona State University.
Conflict of interest
The authors declare no conflict of interest.
eccentricity true anomaly, radians orbit radius, m center of mass identity matrix principal moment of inertia about the roll, pitch, and yaw axes, respectively, kgm2 local-vertical/local-horizontal semilatus rectum Lyapunov-Floquet transformation matrix (LFT) time-dependent normal forms time-independent normal forms earth gravitational parameter, spacecraft angular velocity, rad/s spacecraft pitch angle, spacecraft roll angle, spacecraft yaw angle, state transition matrix (STM) Floquet transition matrix (FTM)
true anomaly, radians
orbit radius, m
center of mass
principal moment of inertia about the roll, pitch, and yaw axes, respectively, kgm2
Lyapunov-Floquet transformation matrix (LFT)
time-dependent normal forms
time-independent normal forms
earth gravitational parameter,
spacecraft angular velocity, rad/s
spacecraft pitch angle,
spacecraft roll angle,
spacecraft yaw angle,
state transition matrix (STM)
Floquet transition matrix (FTM)
- These authors contributed equally.