1.Introduction
The asteroid explorer HAYABUSA was launched on
May,9,2003. Along the way, HAYABUSA has lost two
of three wheels, its attitude is currently stabilized by
only one-axis bias-momentum wheel. Based on an ideal
theory, such satellites are stable for 3 axes since one axis
is stabilized through exchange of angular momentum
with the wheel under a PD control, and the other two
axes are stabilized by the gyroscopic stiffness derived
from the momentum of the wheel.
HAYABUSA had been regarded as such an ideal satellite
as well. However, with a little products of inertia
caused by the fuel leakage, the wheel axes are different
from the principal axes, and thus the dynamics differs
from ideal ones. Because of this misaligment, a nutational
divergence happened to HAYABUSA in March
2007.
This paper analyzes the cause of this phenomenon
and discuss the influences of such misalignment upon
a bias-momentum one wheel satellite like HAYABUSA.
2.Theory
2-1 Equation of motion
To analyze the dynamics of HAYABUSA, consider a
typical bias-momentum satellite with only one momentum
wheel. We define an inertial reference frame [a0]
and a body-fixed reference frame [a1] and a principal
axes-fixed reference frame [a2]. A direction cosine vector
of the wheel axis which expresses a difference with
the principal axis is described by T = (1, 2, 3) in
[a2](Fig.1). We suppose the wheel is pointed at around
z-axis of [a2].
Fig.1: Tilt of the wheel
from the principal axes-fixed frame
Consider Aij,ij,ij as a direction cosine matrix, an angular
velocity vector and an angular vector, respectively.
They are described in [ai] relative to [aj]. The whole angular
momentum vector L in [a1] is written as
L = J110 + Hw (1)
where J1 is the inertia matrix in [a1] and Hw is the angular
momentum vector of the wheel. For simplicity,
neglecting disturbance torque, Euler equation is written
as
J1? 10 + ? 10(J110 + Hw) = ! ?H w (2)
where ?Hw is the control torque and ? x is a matrix describing
cross product x.
We use rotation vector to express the attitude. Using e
as a direction cosine of Eulerfs Eigenaxis and 0 E E
as rotation angle about the axis, rotation vector and
its time derivative is defined as1)
T = (x, y, z) = eT (3)
?
= + 1
2
?
+ ?2
=
2 sin 2
! cos 2
22 cos 2
1
12 E E
1
2 (0 E E )
(4)
Eqs.(2)(4) is written as
J1? 10 + ? 10(J110 + Hw) = ! ?Hw (5)
?
10 = 10 + 1
2
?
1010 + ?2
1010 (6)
in [a1].
In general, we can use below equations for transforms
among reference frames.
A21 J1A12 = J2
A2110 + 21 = 20
A2110 + 21 = 20
A21 ? xA12 = gA21x
(7)
Since the definitions about [a1], [a2] and , J2 is a principal
inertia matrix which consists of principal moments
of inertia Jx, Jy, Jz and we can write
21 = 0
?
21 = 0
(8)
A21Hw = h
A21 ?Hw = ?h
(9)
where h is the magnitude of Hw.
By combining A21 and Eqs.(5)(6), and using above
equations, we obtain
J2? 20 + ? 20(J220 + h) = !?h (10)
d
dt
(20 ! 21) =20 + 1
2
(20 ! 21) 20
+[(20 ! 21) f(20 ! 21) 20g]
(11)
2-2 Control law
Consider a control law of the wheel rotation. The
control law of HAYABUSA is PD control using angular
velocity and angle about z-axis of [a1].
That is to say, the feedback values for the PD control
are the z-component of the angular velocity and the
angle observed in [a1], or 10z, 10z. We can write the
magnitude of control torque ?h = j ?Hwj as
?h
= k110z + k210z (12)
where k1, k2 are the feedback gains.
We define a vector T = (1, 2, 3) such as
A12 =
(
? ?
)
(13)
This describes products of inertia in [a1]. Using and
, we can express the tilt of the wheel (or wheel array
error) and the products of inertia independently(Fig.2).
With above equations about transform of the references,
?h
is written as
?h
= k1T20 + k2T(20 ! 21) (14)
Fig.2: Description of the tilt and the products of inertia
in the principal axes-fixed frame
2-3 Characteristic equation
To obtain the characteristic equation, we linearize
Eqs.(10)(11) by restricting a range of motion properly.
Defining
J = J2
= 20
= 20 ! 21
(15)
we obtain the below system of differential equations.
J? + ? (J + h) = !k1ʃT ! k2ʃT (16)
?
= + 1
2
?
+ ?2 (17)
Consider restricting the range of angular velocity, angle
and the wheel rotation as near its zero, a certain
angle 0 and bias angular momentum h0, respectively.
Namely we suppose
=
= 0 +
h = h0 + h
(18)
Substituting to Eqs.(16)(17) and neglecting 2 terms, we
obtain the below system of linear differential equations
about ,.
J? + h0? = !k1ʃT ! k2ʃT (19)
?
= + 1
2
?
0 + 0 ?0
2
(20)
Its matrix expression is
d
dt
(
)
=
(
!J!1(k1ʃT ! h?) !J!1k2ʃT
E + 12
?0 + ?0
2 0
)
(21)
Finally, we obtain the characteristic equation
det
(
sE + J!1(k1ʃT ! h?) J!1k2ʃT
!E ! 12
?0 ! ?0
2 sE
)
= 0 (22)
where E is a unit matrix.
3.Stability analysis of the system
In this section, we analyze the stability of the system
(,) by eigenvalue analysis using Eq.(22).
3-1 Conditions about stability
Though Eq.(22) is a sixth-order equation about s, it
has double root at s = 0. Its eigenvectors are
xT
1 =
(
0 0 0 3 0 !1
)
xT
2 =
(
0 0 0 2 !1 0
) (23)
which means some deviation about angle remain and
the double root do not affect the stability. We can consider
the stability by analysis of the other 4 eigenvalues.
The characteristic equation is now fourth-order like
s4 + a1s3 + a2s2 + a3s + a4 = 0 (24)
Using Hurwitz criterion, we can say all real parts of the
eigenvalues are negative if and only if2)
ai > 0 (25)
and
3 = det
?
?
a1 a3 0
1 a2 a4
0 a1 a3
?
? > 0 (26)
3-2 Approximation by supposing magnitudes
It is quite hard to analyze the sign of ai, 3 because actual
expanded forms are quite complex and large with
some characteristic constants and the angle 0, where
characteristic constants mean J, , , k1, k2 and h0. So we
suppose the magnitudes of these constants and extract
only dominant terms in ai and 3.
When the tilt of the wheel and the products of inertia
is small, and is approximated as
T =
(
1 2 1
)
T =
(
1 2 1
) (27)
where 1, 2, 1 and 2 is also approximated small. Besides
we suppose j0j < 1[rad]. Based on at Eq.(4) as
the basis of the magnitudes, we suppose
i ? i
i ? i
h0
Ji ?
k1 ? Ji
k2 ? h0
(28)
and neglect terms which are less or equal to 3i
.
With these suppositions, ai is always positive. Meanwhile,
3 is not necessarily positive. It is written as
first-order form about T
0 = (0x, 0y, 0z) like
3 = A0x + B0y + C0z + D (29)
where
A = h0k1k2
2Jy Jz
[
xz21 + yz
1
Jx
2
]
B = !
h0k1k2
2Jx Jz
[
yz22 + xz
1
Jy
1
]
C = !
h0k1k2
2Jz
[ xz2
Jy
21
+ yz2
Jx
22
+ h0k1
Jx Jy Jz
xy12
]
D = k1h01
Jz
[ h30
k11
J2
x J2
y
xz +
(
!
h20
k2
1
J2
x Jy Jz
+
h20
k2
J2
x Jy !
k2
2
Jx Jz
)
yz2
]
+ k1h02
Jz
[ h30
k12
J2
x J2
y
yz +
( h20
k2
1
Jx J2
y Jz !
h20
k2
Jx J2
y
+
k2
2
Jy Jz
)
xz1
]
!
h20
2
1
Jx Jy
[ xz
Jy
21
+ yz
Jx
22
]
+
h0k1yx12
J2
z
[ Jz2
1
Jx Jy ! k22
]
(30)
where
ij = 1
Ji !
1
Jj
1 = h0k1
Jz
2 = k2
Jz !
h20
Jx Jy
(31)
We can not determine the sign of 3 unique, which
means a boundary about the stability exists depending
on the constants in Eq.(29). Some examples of
the boundaries are shown at Fig.3. We describe them
only on a plane (0x, 0y) approximately, because C is
smaller than A and B. The boundary varies its shape
depending on the signs of A and B.
Based on above analysis, a system (or an attitude) of
certain satellite with and may go into stable or unstable
depending on its 0 and it may converge or diverge.
We can say that the cause of these phenomena is the
composition of the feedback value. It includes x and yangle
(see Eq.(14)) which does not relate to control low
if the satellite does not have some misalignment i, i,
and this results changing the dynamics of the system.
The biased feedback value about angle takes place these
boundaries. So we can also say that the way of the definition
of 0 ,or the definition of inertial frame varies the
dynamics.
Fig.3: Examples of boundary about stability
4.Influence of the stability upon the motion
4-1 Analysis of the divergent motion of HAYABUSA
Based on the above stability analysis, we investigate
the divergent motion of HAYABUSA shown in Section
1. When this phenomenon happened, HAYABUSAfs
characteristic constants were
J =
?
?
352.4 0 0
0 268.2 0
0 0 428.3
?
?[kgm2]
h0 = !2.90[kgm2/s]
k1 = 114
k2 = 15.35
T =
(
0.0823 !0.0100 0.9966
)
T =
(
0.0823 !0.0100 0.9966
)
(32)
Here, = means the wheel axis is equal to the z-axis
of [a1]. And it happened after an attitude maneuver
about the x-axis, its rotation vector was
0 =
(
0.393 0.021 0.000
)
[rad] (33)
In this settings, 3 in Eq.(29) is calculated as
3 = A0x + B0y + C0z + D
= !3.34 10!80x + 6.37 10!90y@
+ 2.81 10!90z + 1.58 10!13
= !1.30 10!8
(34)
which means the boundary shape is like upper right of
Fig.3, and the angle (0x, 0y) = (0.393, 0.021) is placed
in the unstable region. We can say this is why the unstable
nutation divergence occurred.
To compare a numerical calculation using
Eqs.(16)(17)(32)(33) with HAYABUSAfs angle record,
we show the result of the calculation at Fig.4. The uppers
of Fig.4 is the records, the lower is the calculation.
These graph consist of time biased to zero when the
divergence occurred for the transverses and z-angle
described by rotation vector for verticals.
Compared with the record, the calculation corresponds
well.
-0.4
-0.35
-0.3
-0.25
-0.2
-0.15
-0.1
-0.05
0
0.05
0.1
0 50000 100000 150000 200000 250000
angle_z [deg]
time [sec]
-0.45
-0.4
-0.35
-0.3
-0.25
-0.2
-0.15
-0.1
-0.05
0
0.05
0.1
0 50000 100000 150000 200000 250000
angle_z [deg]
time [sec]
Fig.4: Comparison of HAYABUSA and simulation
about angle z
4-2 Nutational damping using the stability
We showed unstable motion was caused because its
angle was placed in the unstable region. In reverse, it
may be possible to damp some nutational motion when
we keep its angle in the stable region.
To show numerically, we place the angle as
0 =
(
!0.800 0.800 0.000
)
[rad] (35)
and run a numerical calculation with the constants
Eq.(32) again. We saw that Eq.(32) draws its stability
like upper right of Fig.3. And Eq.(35) instead of Eq.(33)
means the angle is in the stable region. For example, the
result is shown in Fig.5, where nutation damping can be
observed. It can be said that one momentum wheel can
act as a nutation damper if the attitude is kept in the
stable region.
0.755
0.76
0.765
0.77
0.775
0.78
0.785
0.79
0.795
0 1 2 3 4 5 6 7 8 9 10
angle_y [rad]
time [10^4 sec]
Fig.5: simulation in stable region
4-3 Stabilization strategy at usual operation of
HAYABUSA
The unstable motion of HAYABUSA occurred after
the large attitude maneuver. As is seen in section 4-2,
the targeted angle 0 was in the unstable region and
the attitude motion diverged. This is an unusual case,
and at usual operation, a target angle is set to zero after
every attitude maneuver, and the difference from zero is
recognized as error angle to be feedbacked. This section
shows that the unstable motion can be avoided around
this angle (0 = 0).
Consider the stability in this case. Eq.(29) is now written
as
3 = D (36)
The stability depends on only the sign of D. D = 1.58
10!13 on HAYABUSA means stable at 0 = 0.
We consider the magnitude of angle d which expresses
the distance between 0 = 0 and the boundary. It is
written as
d = D
pA2 + B2 + C2
(37)
Here, d = 2.65 10!4[deg] for HAYABUSA, and this
means HAYABUSA is always quite close to the unstable
region. This means it is possible that the angle goes
into the unstable region easily because the zero 0 =
0 is always near the boundary since generally jDj ?
jAj, jBj.
Fig.6: the distance between 0 = 0 and the boundary
However, the speed of divergence or convergence is
quite slow around the boundary. For example, some
numerical calculations show the amplitude of the nutational
motion grows only less than 1[deg] in one week.
At usual operation, HAYABUSA is not be kept left for
a week without attitude correction. Based on this assumption,
the amplitude of the nutational motion does
little increase or decrease around 0 = 0, and thus
HAYABUSA will be stabilized enough if its angle is
reset to 0 = 0 after every attitude maneuver. Here,
resetting to 0 = 0 means resetting the inertial frame
properly, we can do the maneuver easily. It is a guarantee
that no problems happened when the angle had
been around 0 = 0 at usual operation.
5.Use for nutational damper
In section4-2, we saw that taking an intentional angle
0 in stable region can damp nutational motion.
We also saw the speed of the amplitude of the nutational
motion decreasing is very slow around 0 =
0. Of course it is suitable as a nutational damper
that the speed is fast, and this becomes possible by
setting 0 at center of the stable region like Eq.(35).
However, the result of the numerical calculation at
section4-2 showed that the angle converges to around
T = (!0.800, 0.776, 0.074)[rad] (only y-angle is shown
at Fig.5) which has no concerns with the initial angle
Eq.(35). The reason is that when the system converges,
the control torque ?h also converges to zero, or satisfies
T = 0 (38)
(see Eq.(14)). This shapes a plane in the space
(x, y, z) and it is possible to converge wherever on
the plane. We can not control the convergent point,
which takes place some unforeseen errors. It is not
suitable for satellites when the onboard equipments require
high attitude accuracy. We should consider better
method of taking the intentional angle 0 to minimize
the error.
6.Conclusion
To analyze the nutational divergence on HAYABUSA,
we analytically derived the stability condition of a satellite
which is controlled by one bias-momentum wheel.
We showed that when the wheel is controlled under
PD control using the attitude-angle and angular velocity
about the wheel axis, the feedback value is biased
according to products of inertia, and this results in the
change of stability. Because feedbacked attitude-angle
varies the stability, its nutational motion diverge or
converge depending on its attitude-angle. We showed
angle-boundaries about stability, and based on this assumption,
we concluded the cause of the divergence on
HAYABUSA is that its attitude-angle was placed in the
unstable region.
To avoid this unstability, we considered two strategies.
One was to set the attitude-angle zero, the other
was to set the attitude-angle in the stable region. Setting
the attitude-angle is done by changing the definition of
the inertial frame. The former resulted eliminating the
biased feedback value and the nutational motion does
not diverge and converge. The latter, we could see the
motion converges, and we proposed one of the way to
use one momentum wheel as a nutatinal damper.
References
[1] B.Sacleux, : Rotation Vector-Based Attitude Control
Design,AIAA Guidance, Navigation,and Control Conference
and Exhibit, Portland, OR; UNITED STATES;9-
11 Aug. 1999. pp. 1941-1948. (1999)
[2] gPv, 䑺:㐧_, W, pp.94 (1994)