Spacecraft Attitude Stabilization with Piecewise-constant Magnetic Dipole Moment

Spacecraft attitude control can be obtained by adopting several actuation mechanisms. Among them electromagnetic actuators are widely used for generation of attitude control torques on satellites flying low Earth orbits. They consist of planar current-driven coils rigidly placed on the spacecraft typically along three orthogonal axes, and they operate on the basis of the interaction between the magnetic moment generated by those coils and the Earth’s magnetic field; in fact, the interaction with the Earth’s field generates a torque that attempts to align the total magnetic moment in the direction of the field. Magnetic actuators, also known as magnetorquers, have the important limitation that control torque is constrained to belong to the plane orthogonal to the Earth’s magnetic field. As a result, sometimes a reaction wheel accompanies magnetorquers to provide full three-axis control (see [1, Section 7.4]); moreover, magnetorquers are often used for angular momentum dumping when reaction or momentum-bias wheels are employed for accurate attitude control (see [1, Section 7.5]). Lately, attitude stabilization using only magnetorquers has been considered as a feasible option especially for low-cost micro-satellites and for satellites with a failure in the main attitude control system. In such scenario many control laws have been designed, and a survey of various approaches adopted can be found in [2]; in particular, Lyapunov based design has been adopted in [3, 4, 5, 6, 7]. In those works feedback control laws that require measures of both attitude and attitude-rate (i.e. state feedback control laws) are proposed; moreover, in [5] and [7] feedback control algorithms which need measures of attitude only (i.e. output feedback control algorithms) are presented, too.

All control laws presented in the cited works are continuos-time; thus, in principle they require a continuous-time measurement of the geomagnetic field and a continuos-time generation of coils’s currents. However, in practical implementations the time in which Earth’s magnetic field is measured (measurement time) has to be separated from the time in which coils’ currents are generated (actuation time). The latter separation is necessary because currents flowing in the spacecraft’s magnetic coils generate a local magnetic field that impairs the measurement of the geomagnetic field; as a result, when the Earth magnetic field is being measured no currents must flow in the coils, and consequently no magnetic actuation is possible. As a result, in practical applications a periodic switch between measurement time and actuation time is implemented. Usually measurement time is set much shorter than actuation time, so that the spacecraft is left unactuated for a very short fraction of the period; thus, the measurement process can be modeled as instantaneous, and during each period the magnetic dipole moment is kept constant since only the measurement at the beginning of the period is available. Thus the separation illustrated before translates into the constraint of being able to generate only constant magnetic dipole moments during each period. Then, magnetic control laws compatible with the latter constraint could be obtained by simply discretizing continuos-time control algorithms (see e.g. [8]); however, here a different approach is proposed; it consists in taking into account of the latter constraint directly during the control design phase. It will be shown that the proposed approach leads to a systematic selection of the sampling interval; the latter aspect represents an advantage with respect to discretization of continuous-time control laws where the sampling interval is usually determined by trial an error. The proposed method is presented for both a state and an output feedback .

In order to describe the attitude dynamics of an Earth-orbiting rigid spacecraft, and in order to represent the geomagnetic field, it is useful to introduce the following reference frames.

1. 1.

   Earth-centered inertial frame $\mathcal{F}_{i}$. A commonly used inertial frame for Earth orbits is the Geocentric Equatorial Frame, whose origin is in the Earth’s center, its $x_{i}$ axis is the vernal equinox direction, its $z_{i}$ axis coincides with the Earth’s axis of rotation and points northward, and its $y_{i}$ axis completes an orthogonal right-handed frame (see [1, Section 2.6.1] ).

2. 2.

   Spacecraft body frame $\mathcal{F}_{b}$. The origin of this right-handed orthogonal frame attached to the spacecraft, coincides with the satellite’s center of mass; its axes are chosen so that the inertial pointing objective is having $\mathcal{F}_{b}$ aligned with $\mathcal{F}_{i}$.

Since the inertial pointing objective consists in aligning $\mathcal{F}_{b}$ to $\mathcal{F}_{i}$, the focus will be on the relative kinematics and dynamics of the satellite with respect to the inertial frame. Let $q=[q_{1}\ q_{2}\ q_{3}\ q_{4}]^{\textrm{T}}=[q_{v}^{\textrm{T}}\ q_{4}]^{\textrm{T}}$ with $\|q\|=1$ be the unit quaternion representing rotation of $\mathcal{F}_{b}$ with respect to $\mathcal{F}_{i}$; then, the corresponding attitude matrix is given by

(see [9, Section 5.4]).

Let

Then the relative attitude kinematics is given by

where $\omega\in\mathbb{R}^{3}$ is the angular rate of $\mathcal{F}_{b}$ with respect to $\mathcal{F}_{i}$ resolved in $\mathcal{F}_{b}$ (see [9, Section 5.5.3]).

The attitude dynamics in body frame can be expressed by

where $J\in\mathbb{R}^{3\times 3}$ is the spacecraft inertia matrix, and $Q\in\mathbb{R}^{3}$ is the control torque expressed in $\mathcal{F}_{b}$ (see [9, Section 6.4]).

The spacecraft is equipped with three magnetic coils aligned with the $\mathcal{F}_{b}$ axes which generate the magnetic control torque

where $m_{coils}\in\mathbb{R}^{3}$ is the vector of magnetic dipole moments for the three coils, and $B^{b}$ is the geomagnetic field at spacecraft expressed in body frame $\mathcal{F}_{b}$.

Let $B^{i}$ be the geomagnetic field at spacecraft expressed in inertial frame $\mathcal{F}_{i}$. Note that $B^{i}$ varies with time at least because of the spacecraft’s motion along the orbit. Then

which shows explicitly the dependence of $B^{b}$ on both $q$ and $t$.

Grouping together equations (4) (5) (6) the following nonlinear time-varying system is obtained

in which $m_{coils}$ is the control input.

In order to design control algorithms, it is important to characterize the time-dependence of $B^{b}(q,t)$ which is the same as characterizing the time-dependence of $B^{i}(t)$. Assume that the orbit is circular of radius $R$; then, adopting the so called dipole model of the geomagnetic field (see [10, Appendix H]) we obtain

In equation (9), $\mu_{m}=7.746\ 10^{15}$ Wb m is the total dipole strength (see [11]), $r^{i}(t)$ is the spacecraft’s position vector resolved in $\mathcal{F}_{i}$, and $\hat{r}^{i}(t)$ is the vector of the direction cosines of $r^{i}(t)$; finally $\hat{m}^{i}$ is the vector of direction cosines of the Earth’s magnetic dipole moment expressed in $\mathcal{F}_{i}$. Here we set $\hat{m}^{i}=[0\ 0\ -1]^{\textrm{T}}$ which corresponds to setting the dipole’s colelevation equal to $180^{\circ}$. Even if a more accurate value for that angle would be $170^{\circ}$ (see [11]), here we approximate it to $180^{\circ}$ since this will substantially simplify future symbolic computations.

Equation (9) shows that in order to characterize the time dependence of $B^{i}(t)$ it is necessary to determine an expression for $r^{i}(t)$ which is the spacecraft’s position vector resolved in $\mathcal{F}_{i}$. Define a coordinate system $x_{p}$, $y_{p}$ in the orbital’s plane whose origin is at Earth’s center; then, the position of satellite’s center of mass is clearly given by

where $n$ is the orbital rate, and $\phi_{0}$ an initial phase. Then, coordinates of the satellite in inertial frame $\mathcal{F}_{i}$ can be easily obtained from (10) using an appropriate rotation matrix which depends on the orbit’s inclination $incl$ and on the right ascension of the ascending node $\Omega$ (see [1, Section 2.6.2]). Plugging into (9) the expression of the latter coordinates, an explicit expression for $B^{i}(t)$ can be obtained; it can be easily checked that $B^{i}(t)$ is periodic with period $2\pi/n$. Consequently system (8) is a periodic nonlinear system with period $2\pi/n$.

As stated before, the control objective is driving the spacecraft so that $\mathcal{F}_{b}$ is aligned with $\mathcal{F}_{i}$. From (2) it follows that $C(q)=I$ for $q=[q_{v}^{\textrm{T}}\ q_{4}]^{\textrm{T}}=\pm\bar{q}$ where $\bar{q}=[0\ 0\ 0\ 1]^{\textrm{T}}$. Thus, the objective is designing control strategies for $m_{coils}$ so that $q_{v}\rightarrow 0$ and $\omega\rightarrow 0$.

In [7] both a static state feedback and a dynamic output feedback are presented; both feedbacks were obtained as modifications of those in [5, 6].

The static state feedback control proposed in [7] is given by

and it is shown that picking $k_{1}>0$, $k_{2}>0$, and choosing $\epsilon>0$ and small enough, local exponential stability of $(q,\omega)=(\bar{q},0)$ is achieved.

The dynamic output feedback proposed in [7] is given by

with $\delta\in\mathbb{R}^{4}$. It is shown that selecting $k_{1}>0$, $k_{2}>0$, $\alpha>0$, $\lambda>0$, and choosing $\epsilon>0$ small enough, local exponentially stability of equilibrium $(q,\omega,\delta)=(\bar{q},0,\frac{1}{\epsilon\lambda}\bar{q})$ is achieved.

Both control laws (11) and (12) seem implementable in practice since $B^{b}(q(t),t)$ can be measured using magnetometers, $q_{v}(t)$ can be measured using attitude sensors, and $\omega(t)$ can be measured by attitude rate sensors. However, as explained in the introduction, the time in which $B^{b}$ is measured should be separated from the time in which control action is applied; in fact, magnetic control torque is obtained by generating currents flowing through magnetic coils, and those currents create a local magnetic field which impairs the measurement of $B^{b}$. As a result, in order to take into account of the previous constraint, we should consider only feedback laws in which $B^{b}$ is instantaneously measured at the beginning of some fixed length interval, $m_{coils}$ is held constant over that interval, and the latter operations are repeated periodically. In the sequel $T$ will denote the sampling interval, that is the length of each of those intervals.

A static state feedback law that fulfills the previous constraint can be designed by simply discretizing (11) through a sample-and-hold operation, thus obtaining the following control algorithm

Similarly, a dynamic output feedback law compatible with the constraint of generating a piecewise-constant signal for $m_{coils}$ can be obtained by discretizing (12) using the forward differencing approximation (see [12]), thus obtaining the following control law

with $\delta\in\mathbb{R}^{4}$.

In this section we will focus on state-feedback law (13) and show how to choose the parameters $k_{1}$, $k_{2}$, $\epsilon$, and the sampling interval $T$ so that equilibrium $(q,\omega)=(\bar{q},0)$ is locally exponentially stable for closed-loop system (8) (13). In order to derive those design criteria, it suffices considering the restriction of closed-loop system (8) (13) to the open set $\mathbb{S}^{3+}\times\mathbb{R}^{3}$ where

Since on the latter set $q_{4}=(1-q_{v}^{\textrm{T}}q_{v})^{1/2}$ then the considered restriction is given by the following nonlinear time-varying hybrid system

where

and

and where (7) has been used. It is simple to show that if $(q_{v},\omega)=(0,0)$ is a locally exponentially stable equilibrium for (16), then $(q,\omega)=(\bar{q},0)$ is a locally exponentially stable equilibrium for (8) (13).

Next, consider the linear approximation of (16) about the equilibrium $(q_{v},\omega)=(0,0)$ as defined in [13]. The latter approximation is given by

Note that since $B^{i}$ is bounded (see (9)), then Theorem II.1 in [13] applies to the nonlinear time-varying hybrid system (16); consequently, $(q_{v},\omega)=(0,0)$ is a locally exponentially stable equilibrium for (16) if and only if it is an exponentially stable equilibrium for (19) (20).

Next consider the continuous-time system (19) and sample its state $[q_{v}(t)^{\textrm{T}}\omega(t)^{\textrm{T}}]^{\textrm{T}}$ at $t=(k+1)T$. The following discrete-time system is thus obtained

where

and

By [14, Proposition 7] it follows that if the linear time-varying discrete-time system (20) (21) is exponentially stable then the linear time-varying hybrid system (19) (20) is exponentially stable.

Based on the previous discussion, our objective has become studying stability of linear time-varying discrete-time system (20) (21). For that purpose it is useful to perform the following change of variables $z_{1}(k)=q_{v}(kT)$, $z_{2}(k)=\omega(kT)/\epsilon$ obtaining

Note that since $G_{i}(k,0)=0\ \ i=1,2$, then $G_{i}(k,T)$ can be expressed as

where

Let $L(k,T)=H_{2}(k,T)B^{i}(kT)^{\times T}$. Note that since $B^{i}(t)$ is periodic, then the following average

is well defined and indipendent of integer $s$. The expression of $L_{av}$ has been computed symbolically using Mathematica^TM, but here it is omitted to save space. It turns out that $L_{av}^{0}\triangleq\lim_{T\rightarrow 0}L_{av}(T)$ is symmetric, and that $L_{av}^{0}$ is positive definite if and only if the orbit is not equatorial, that is if its inclination satisfies $incl\neq 0$. Thus, we make the following assumption.

At this point, consider the average system of (22) as defined in [15] which is given by

It can be verified that all assumptions of [15, Theorem 2.2.2] are fulfilled by system (22). Thus, the following proposition holds true.

For the average system (24) the following stability result holds.

Based on the previous theorem, fix $0<T<T^{*}$, and consider the corresponding value of $\epsilon_{0}$ given by (28); from Proposition 2 it follows that there exists $0<\epsilon_{1}\leq\epsilon_{0}$ such that system (22) is exponentially stable for all $0<\epsilon\leq\epsilon_{1}$. In conclusion, from Theorem 1, Proposition 1, and the preceding discussion, the following main proposition is obtained.

In this section we will focus on output-feedback law (14) and show how to choose parameters $k_{1}$, $k_{2}$, $\alpha$, $\lambda$, $\epsilon$, and sampling interval $T$ so that equilibrium $(q,\omega,\delta)=(\bar{q},0,\frac{1}{\epsilon\lambda}\bar{q})$ is locally exponentially stable for closed-loop system (8) (14).

Following the same approach as in state-feedback design, we determine the linear time-varying discrete-time system given by the interconnection of (21) with the following system

where $\delta_{v}=[\delta_{1}\ \delta_{2}\ \delta_{3}]^{\textrm{T}}$ and $\tilde{\delta}_{4}=\delta_{4}-\frac{1}{\epsilon\delta}$. Similarly to the state feedback case, if system (29) is exponentially stable, then $(q,\omega,\delta)=(\bar{q},0,\frac{1}{\epsilon\lambda}\bar{q})$ is an exponentially stable equilibrium for (8) (14).

For system (29), perform the change of variables $z_{1}(k)=q_{v}(kT)$, $z_{2}(k)=\omega(kT)/\epsilon$, $z_{3}(k)=q_{v}(kT)-\epsilon\lambda\delta_{v}(kT)$, $z_{4}(k)=\tilde{\delta}_{4}(kT)$ obtaining

The average system of (22) as defined in [15] is given by

Then, a proposition parallel to Proposition 1 holds true, and for the average system (31) the following stability result holds.

Thus, the following main proposition holds.

We consider the same case study presented in [7] in which the spacecraft’s inertia matrix is given by $J=\text{diag}[27\ 17\ 25]\ \text{kg m}^{2}$; the inclination of the orbit is given by $incl=87^{\circ}$, and the orbit’s altitude is 450 km; the corresponding orbital period is about 5600 s. Without loss of generality the right ascension of the ascending node $\Omega$ is set equal to 0, whereas the initial phase $\phi_{0}$ (see (10)) has been randomly selected and set equal to $\phi_{0}=0.94$ rad.

Consider an initial state characterized by attitude equal to to the target attitude $q(0)=\bar{q}$, and by the following high initial angular rate

due for example to an impact with an object.

In [7] the continuous-time feedback (11) has been designed setting $k_{1}=2\ 10^{11}$, $k_{2}=3\ 10^{11}$, and $\epsilon=10^{-3}$. Here, we keep $k_{1}=2\ 10^{11}$, $k_{2}=3\ 10^{11}$ and parameters $T$ and $\epsilon$ are chosen as follows. First, studying numerically the behavior of the eigenvalues of $A_{s}(T)$ with $T$ (see (26)), we determine the value $T^{*}=1490$ s for which it occurs that $A_{s}(T)$ is Hurwitz for all $0<T<T^{*}$. Then, the sampling time has been set equal to $T=20$ to which there corresponds $\epsilon_{0}=1.3\ 10^{-3}$ (see (28)). By Proposition 2, it occurs that setting $\epsilon\leq 1.3\ 10^{-3}$ and small enough, equilibrium $(q,\omega)=(\bar{q},0)$ is locally exponentially stable for (8) (13). Here, proceeding by trial and error, we fix $\epsilon=10^{-3}$. Simulations results reported in Fig. 1 show that actually acquisition of the desired attitude is achieved.

A similar approach has been followed to select parameters and sampling interval for output-feedback (14); the obtained results validate Proposition 3, but they are omitted for lack of space.

In this work two magnetic control laws for spacecraft attitude stabilizations have been presented. Both magnetic control laws possess the feature of generating piecewise-constant magnetic dipole moments; the latter aspect makes them easier to implement than continuos-time control laws; in fact, implementation of continuos-time laws requires a discretization process, in which the discretization interval is usually selected by trial and error; the latter trial and error process in not necessary using the approach here proposed since the sampling interval can be determined in a systematic way.

The author acknowledges Prof. J. L. Mancilla-Aguilar for his help and Prof. A. Nascetti for fruitful discussions.