# Path Planning for Rapid Large-Angle Maneuver of Satellites Based on the Gauss Pseudospectral Method

Authors: Wan Zhang; Yao Zhang; Wenbo Li; Youyi Wang
Journal: Mathematical Problems in Engineering
 (2016)

---

## Abstract

A Gauss pseudospectral method is proposed in this study to solve the optimal trajectory-planning problem for satellite rapid large-angle maneuvers. In order to meet the requirement of rapid maneuver capability of agile small satellites, Single Gimbal Control Moment Gyros (SGCMGs) are adopted as the actuators for the attitude control systems (ACS). Because the singularity problem always exists for SGCMGs during the large-angle maneuvering of the satellites, a trajectory optimization method for the gimbal rates is developed based on the Gauss pseudospectral method. This method satisfies the control requirement of satellite rapid maneuvers and evades the singularity problem of SGCMGs. The method treats the large-angle maneuver problem as an optimization problem, which meets the boundary condition and a series of the physical constraints including the gimbal angle constraint, the gimbal rates constraint, the singularity index constraint, and some other performance criteria. This optimization problem is discretized as a nonlinear programming problem by the Gauss pseudospectral method. The optimal nonsingularity gimbal angle trajectory is obtained by the sequence of quadratic programming (SQP). This approach avoids the difficulties in solving the boundary value problem. The simulations reveal that the Gauss pseudospectral method effectively plans an optimal trajectory and satisfies all the constraints within a short time.


---

## Body

## 1. Introduction

With the rapid development of satellites, an increasing number of space missions require that these satellites own the capability of large-angle maneuvering in a short time, typically referred to as the agility of attitude [1]. New satellites, which are also called agile satellites, have wider range of application and greater efficiency. At the same time, agile satellites can provide longer data transmission times. The above-mentioned abilities can achieve these functions, which are hard to reach by traditional techniques, especially for acquiring remote information and approaching on-orbit servicing noncooperative targets [2].Agile satellite attitude maneuvering and control are topics that have recently attracted much attention. In [3], the authors mainly discussed the problem of agile satellite attitude tracking, and an adaptive control law was designed for restraining the disturbance and the uncertainty of the inertia matrix. Wie et al. [4] researched the problem of rapidly capturing the target of an agile satellite with Single Gimbal Control Moment Gyros (SGCMGs), and a nonlinear feedback control law was proposed based on the sign function. Reference [5] discussed how to transfer an attitude-tracking problem with a time-variable boundary condition into a fixed-boundary problem by means of line-of-sight coordinates, rotational quaternion, and relative error feedback signals.Considering that the satellites require rapid maneuvering capability, it is necessary to choose proper actuators to provide the control moment. At present, Control Moment Gyro (CMG) is the most globally studied actuator. According to the amount of framing, CMG can be divided into Single Gimbal Control Moment Gyro (SGCMG) and Double Gimbal Control Moment Gyro (DGCMG). SGCMG can output larger moment, quicker responses, less power consumption, greater operational time, easier vibration isolation detection, and so on [6]. However, the special singularity problem of system configuration limits its application. Therefore, the singularity problem must be considered during SGCMGs steering law design. The existing steering laws have the weakness of consuming too much energy-consumed medium, destroying SGCMGs, affecting the system, and so on [7].In 1978, Margulies and Aubrun [8] published the classic work about the singularity problem, laying a good foundation for related research. Bedrossian et al. [9] applied the singularity robust steering law to the control of SGCMGs. In 2000, Ford and Hall [10] proposed a steering law for avoiding singularity to reduce the error of the output torque. Wie et al. [11, 12] presented a generalized singularity robust steering law and a nondiagonal singularity robust steering law escaping the singular surface, but a remarkable torque error was produced. In [13–15], the SGCMGs steering laws were analyzed and compared for a spacecraft attitude control system based on the singular value decomposition (SVD) theory. The author also developed a new mixture steering law and an error control steering law to avoid the singularity problem.The singularity problem always exists for SGCMGs during the large-angle maneuvering of satellites, so a trajectory optimization method for the gimbal rates is developed in this study. Two methods exist for solving the optimal control problem. One method is the direct approach; the other is the indirect approach. The indirect approach solves a two-point boundary value problem of a set of differential equations. This method has high precision, but its convergence region is very limited, and the solving process is complicated and cumbersome. Furthermore, it is extremely difficult to solve the problem with inequality constrains. In the direct approach, a parameterization method is used to translate the optimization problem into a discrete nonlinear programming problem. The optimal solution can then be developed according to the performance index. This direct approach has a larger convergence region and has high demands for initial estimates. However, disadvantages still exist with this method, such as lower solution precision, slower convergent rate, and more estimators [16]. As one of the pseudospectral methods, the Gauss pseudospectral method is more frequently applied to trajectory optimization due to its high solution precision and fast convergent rate.In this paper, an attitude maneuver optimization method for the gimbal rates is developed based on the Gauss pseudospectral method. It satisfies the control requirement of satellite rapid maneuvering and evades the singularity problem of SGCMGs. This optimization problem is discretized as a nonlinear programming problem by the Gauss pseudospectral method. The optimal nonsingularity gimbal angle trajectory is then obtained by the sequence of quadratic programming (SQP). Finally, the mathematic simulations based on MATLAB/Simulink are conducted.


## 2. Satellite Description

This paper only considers a satellite that uses SGCMGs with a four-pyramid configuration as its actuator. The total angular momentum of the system can be written as(1)
H
=
I
b
w
b
+
h
c
=
I
t
w
b
+
A
g
I
c
g
δ
˙
+
A
s
I
w
s
Ω
,where vector w
b is the inertial angular rate of the satellite expressed in the body frame. Matrix I
b is the general inertia tensor of the satellite without SGCMGs. Matrix I
t is the inertia tensor of the satellite. Matrix I
c
g is the inertia tensor of the SGCMGs about its spin axis of the gimbal, and Matrix I
w
s is the inertia tensor of the SGCMGs about its spin axis of the rotor. Matrix A
s
=
s
1
s
2
⋯
s
n gives the orientation of the spin axis of rotor about the SGCMGs, and Matrix A
t
=
t
1
t
2
⋯
t
n gives the orientation of the output torque about the SGCMGs. Matrices A
s and A
t will change along with the frame angle. Vector Ω is the angular rate of the rotor relative to the actuator frame.Thus, by applying the momentum theorem to (1), the following relation is obtained:(2)
I
˙
t
w
b
+
I
t
w
˙
b
+
A
g
I
c
g
δ
¨
+
A
t
I
w
s
d
Ω
δ
˙
+
w
b
x
I
t
w
b
+
A
g
I
c
g
δ
˙
+
A
s
I
w
s
Ω
=
T
e
x
t
,where vector T
e
x
t is the resultant external momentum and matrix A
g gives the orientation of the gimbal rate.Equation (2) can also be expressed as(3)
I
t
w
˙
b
+
w
b
x
I
t
w
b
+
A
s
I
w
s
Ω
=
T
c
+
T
e
x
t
.In this equation, vectorT
c is the output control torque about the SGCMGs, (4)
T
c
=
-
B
δ
¨
+
C
1
δ
˙
,where(5)
B
=
A
g
I
c
g
,
(6)
C
1
=
C
+
R
c
w
b
+
w
b
x
A
g
I
c
g
,
(7)
C
=
A
t
I
w
s
d
Ω
,
(8)
R
c
w
b
=
I
˙
t
w
b
.The torque due toδ
¨ is less than the torque due to δ
˙, and w
b in C
1 is tiny compared with Ω; thus, the simplification of (4) is(9)
T
c
=
-
h
0
A
t
δ
˙
,where it is supposed that the nominal angular momentum about the rotor of each gyro is equal and can be written as h
0.The attitude kinematics are modeled using the quaternion differential equation(10)
q
˙
0
=
-
1
2
w
T
q
,
q
˙
=
-
1
2
w
x
q
+
1
2
q
0
w
.


## 3. Description of the SGCMGs’ Singularity Problem

This paper considers a representative SGCMG system, which is in a pyramid array of four single gimbal CMGs, as shown in Figure1. This arrangement has minimum redundancy. Each gimbal axis is perpendicular to the side of the pyramid. In this paper, we considered beta as 54.73 degrees.Figure 1
Pyramid mounting arrangement of four single gimbal CMGs.The output moment of the SGCMGs is shown as (4). In fact, simplified equation (9) is the first choice for use in the steering law because it makes the singularity analysis and steering law design easier and more convenient. According to the geometric meaning of A
t, the singularity problem corresponds to the situation when the rank of A
t is less than 3; that is, d
e
t
A
t
A
t
T
=
0. Therefore, the singular value can be defined as follows:(11)
D
=
det
⁡
A
t
A
t
T
.


## 4. Gauss Pseudospectral Method

Pseudospectral methods, also known as “orthogonal collocation methods” in optimal control, arose from spectral methods that were traditionally used to solve fluid dynamics problems [17, 18]. The Gauss pseudospectral method (GPM) is a direct transcription method for discretizing a continuous optimal problem into a nonlinear programming problem (NLP). The Gauss pseudospectral method is an orthogonal collocation discretization method, and it has a faster convergent rate than other methods [19, 20].


### 4.1. Optimal Control Problem in Bolza Form

Without generality loss, the following optimal control problem is described in Bolza form, thus finding the optimal control variableu
t
∈
R
m to minimize the cost functional. Consider(12)
m
i
n
J
=
Φ
x
t
0
,
t
0
,
x
t
f
,
t
f
+
∫
t
0
t
f
g
x
t
,
u
t
,
t
d
t
.Subject to the dynamic constraints, the inequality path constraints and the boundary conditions(13)
s.t.
x
˙
t
=
f
x
t
,
u
t
,
t
,
C
x
t
,
u
t
,
t
≤
0
,
ϕ
x
t
0
,
t
0
,
x
t
f
,
t
f
=
0
,where x
t
∈
R
n is the state and t
0 and t
f are the initial time and final time, respectively. The variable t
∈
t
0
t
f. f is a vector function of n dimensions, C is a vector function of c dimensions, and ϕ is a vector function of q dimensions.For proper handling, the time domain requires transformation fromt
∈
t
0
t
f into τ
∈
-
1
1. They are related as(14)
τ
=
2
t
t
f
-
t
0
-
t
f
+
t
0
t
f
-
t
0
.After transformation, the optimal control problem in Bolza form can be written in the following form:(15)
J
=
Φ
x
-
1
,
x
1
,
t
0
,
t
f
+
t
0
-
t
f
2
∫
-
1
1
g
x
τ
,
u
τ
,
τ
,
t
0
,
t
f
d
τ
s.t.
x
˙
τ
=
t
0
-
t
f
2
f
x
τ
,
u
τ
,
τ
,
t
0
,
t
f
ϕ
x
-
1
,
x
1
,
t
0
,
t
f
=
0
C
x
τ
,
u
τ
,
τ
,
t
0
,
t
f
.


### 4.2. State Approximation Using the Gauss Pseudospectral Method

The Gauss pseudospectral method takes theK order Legendre-Gauss points and the end points as nodes. Then, the state is approximated using a basis of K
+
1 Lagrange interpolating polynomials,(16)
x
τ
≈
X
τ
=
∑
i
=
0
K
X
i
l
i
τ
,
X
i
=
X
τ
i
.The control variable is approximated as(17)
u
τ
≈
U
τ
=
∑
i
=
1
K
U
i
l
~
i
τ
,
U
i
=
U
τ
i
,where l
i and l
~
i are the Lagrange polynomial.Next, the expression in (1) is differentiated with respect to τ that gives(18)
x
˙
τ
k
≈
X
˙
τ
k
=
∑
i
=
0
K
x
τ
i
l
˙
i
τ
k
k
=
1
,
…
,
K
=
∑
i
=
0
K
x
τ
i
D
k
i
τ
k
,where matrix D
∈
R
k
×
(
k
+
1
) is a differentiation matrix.Then, the dynamic constraint is transcribed into algebraic constraints as the following form:(19)
∑
i
=
0
K
X
i
D
k
i
τ
k
-
t
0
-
t
f
2
f
X
k
,
U
k
,
τ
k
,
t
0
,
t
f
=
0
.Next, because the final stateX
f satisfied the dynamic constraint, (18) can be calculated by Gaussian integration; that is,(20)
X
f
-
X
0
-
t
f
-
t
0
2
∑
i
=
1
K
w
i
f
X
τ
i
,
U
τ
i
,
τ
i
,
t
0
,
t
f
=
0
,where τ
i is the LG point and w
i is the LG weight.Finally, the optimal control problem can be expressed as(21)
J
=
Φ
x
-
1
,
x
1
,
t
0
,
t
f
+
t
0
-
t
f
2
∑
i
=
1
K
w
i
g
X
i
,
U
i
,
τ
i
,
t
0
,
t
f
s.t.
∑
i
=
0
K
X
i
D
k
i
τ
k
-
t
0
-
t
f
2
f
X
k
,
U
k
,
τ
k
,
t
0
,
t
f
=
0
,
X
f
-
X
0
-
t
f
-
t
0
2
∑
i
=
1
K
w
i
f
X
τ
i
,
U
τ
i
,
τ
i
,
t
0
,
t
f
=
0
,
ϕ
X
t
0
,
t
0
,
X
t
f
,
t
f
=
0
,
C
X
t
,
U
t
,
t
≤
0
.


## 4.1. Optimal Control Problem in Bolza Form

Without generality loss, the following optimal control problem is described in Bolza form, thus finding the optimal control variableu
t
∈
R
m to minimize the cost functional. Consider(12)
m
i
n
J
=
Φ
x
t
0
,
t
0
,
x
t
f
,
t
f
+
∫
t
0
t
f
g
x
t
,
u
t
,
t
d
t
.Subject to the dynamic constraints, the inequality path constraints and the boundary conditions(13)
s.t.
x
˙
t
=
f
x
t
,
u
t
,
t
,
C
x
t
,
u
t
,
t
≤
0
,
ϕ
x
t
0
,
t
0
,
x
t
f
,
t
f
=
0
,where x
t
∈
R
n is the state and t
0 and t
f are the initial time and final time, respectively. The variable t
∈
t
0
t
f. f is a vector function of n dimensions, C is a vector function of c dimensions, and ϕ is a vector function of q dimensions.For proper handling, the time domain requires transformation fromt
∈
t
0
t
f into τ
∈
-
1
1. They are related as(14)
τ
=
2
t
t
f
-
t
0
-
t
f
+
t
0
t
f
-
t
0
.After transformation, the optimal control problem in Bolza form can be written in the following form:(15)
J
=
Φ
x
-
1
,
x
1
,
t
0
,
t
f
+
t
0
-
t
f
2
∫
-
1
1
g
x
τ
,
u
τ
,
τ
,
t
0
,
t
f
d
τ
s.t.
x
˙
τ
=
t
0
-
t
f
2
f
x
τ
,
u
τ
,
τ
,
t
0
,
t
f
ϕ
x
-
1
,
x
1
,
t
0
,
t
f
=
0
C
x
τ
,
u
τ
,
τ
,
t
0
,
t
f
.


## 4.2. State Approximation Using the Gauss Pseudospectral Method

The Gauss pseudospectral method takes theK order Legendre-Gauss points and the end points as nodes. Then, the state is approximated using a basis of K
+
1 Lagrange interpolating polynomials,(16)
x
τ
≈
X
τ
=
∑
i
=
0
K
X
i
l
i
τ
,
X
i
=
X
τ
i
.The control variable is approximated as(17)
u
τ
≈
U
τ
=
∑
i
=
1
K
U
i
l
~
i
τ
,
U
i
=
U
τ
i
,where l
i and l
~
i are the Lagrange polynomial.Next, the expression in (1) is differentiated with respect to τ that gives(18)
x
˙
τ
k
≈
X
˙
τ
k
=
∑
i
=
0
K
x
τ
i
l
˙
i
τ
k
k
=
1
,
…
,
K
=
∑
i
=
0
K
x
τ
i
D
k
i
τ
k
,where matrix D
∈
R
k
×
(
k
+
1
) is a differentiation matrix.Then, the dynamic constraint is transcribed into algebraic constraints as the following form:(19)
∑
i
=
0
K
X
i
D
k
i
τ
k
-
t
0
-
t
f
2
f
X
k
,
U
k
,
τ
k
,
t
0
,
t
f
=
0
.Next, because the final stateX
f satisfied the dynamic constraint, (18) can be calculated by Gaussian integration; that is,(20)
X
f
-
X
0
-
t
f
-
t
0
2
∑
i
=
1
K
w
i
f
X
τ
i
,
U
τ
i
,
τ
i
,
t
0
,
t
f
=
0
,where τ
i is the LG point and w
i is the LG weight.Finally, the optimal control problem can be expressed as(21)
J
=
Φ
x
-
1
,
x
1
,
t
0
,
t
f
+
t
0
-
t
f
2
∑
i
=
1
K
w
i
g
X
i
,
U
i
,
τ
i
,
t
0
,
t
f
s.t.
∑
i
=
0
K
X
i
D
k
i
τ
k
-
t
0
-
t
f
2
f
X
k
,
U
k
,
τ
k
,
t
0
,
t
f
=
0
,
X
f
-
X
0
-
t
f
-
t
0
2
∑
i
=
1
K
w
i
f
X
τ
i
,
U
τ
i
,
τ
i
,
t
0
,
t
f
=
0
,
ϕ
X
t
0
,
t
0
,
X
t
f
,
t
f
=
0
,
C
X
t
,
U
t
,
t
≤
0
.


## 5. Simulation Results and Analysis

In this section, the mathematic simulations based on MATLAB/Simulink are conducted to prove the effectiveness of the Gauss pseudospectral method. The parameters are given as follows.The inertia tensor of the satellite is chosen asI
=
diag
⁡
228
188
176 kg·m2.The state isx
=
q
0
,
q
1
,
q
2
,
q
3
,
w
x
,
w
y
,
w
z
,
δ
1
,
δ
2
,
δ
3
,
δ
4
T. The control variable is u
=
δ
˙
1
,
δ
˙
2
,
δ
˙
3
,
δ
˙
4
T.The satellite rotates 45° about the axis of rolling. This paper considers the time and energy consumption optimal problem. Thus, the cost functional is(22)
J
=
∫
t
0
t
f
1
+
u
T
R
u
d
t
.The first state isx
0
=
1
0
0
0
0
0
0
0
0
0
0
T and the final state is x
f
=
[
0.9239
0.3824
0
0
0
0
0
0
0
0
0
]
T. The path constraints are as follows:(23)
C
1
=
q
0
2
+
q
1
2
+
q
2
2
+
q
3
2
=
1
,
C
2
=
w
≤
1
8
r
a
d
/
s
,
C
3
=
δ
˙
i
≤
1
r
a
d
/
s
,
i
=
1,2
,
3,4
,
C
4
=
A
t
A
t
T
≥
0.2
,
C
5
=
T
≤
25
N
m
.The tolerance error isε
≤
10
-
3 and the time of calculation is 4.1 min. The simulation results are shown in Figures 2–7. The satellite rotates 45° around the axis of rolling, and it only costs 7.836 s.Figure 2
Time response of quaternion.Figure 3
Time response of angular velocity.Figure 4
Time response of gimbal angles.Figure 5
Time response of gimbal angular speeds.Figure 6
Time response of singular value rate.Figure 7
Time response of control torque.


## 6. Controller Design

Next, the optimal data were connected with the model in MATLAB/Simulink, and the PD controller was chosen to compensate for the errors between the target and actual trajectory.The equation governing a satellite is expressed as follows:(24)
I
t
w
˙
b
+
w
b
x
I
t
w
b
+
A
s
I
w
s
Ω
=
T
c
+
T
e
x
t
.We define the quaternionQ
d
=
[
q
d
0
q
d
T
]
T, where q
d
=
[
q
d
1
q
d
2
q
d
3
]
T is the desired reference attitude. The quaternion for the attitude error is Q
e
=
[
q
e
0
q
e
T
]
T with q
e
=
q
e
1
q
e
2
q
e
3
T. Using the multiplication law for quaternions, we then obtain(25)
q
e
0
=
q
0
q
d
0
-
q
T
q
d
,
q
e
=
q
d
0
q
-
q
0
q
d
+
q
x
q
d
.Subject to the constraint(26)
Q
e
T
Q
e
=
1
.We denote byw
d
=
w
d
1
w
d
2
w
d
3
T the desired angular velocity and by w
e
=
w
-
w
d the angular speed error.The kinematic equation for the attitude error can be expressed as(27)
q
˙
e
0
=
-
1
2
q
e
T
w
e
,
q
˙
e
=
1
2
q
e
x
+
q
e
0
E
3
×
3
w
e
.In order to design this controller we first define the following Lyapunov function:(28)
V
=
k
1
q
e
T
q
e
+
1
2
w
e
T
I
w
e
,where k
1
>
0.The first time derivative ofV can be obtained as(29)
V
˙
=
2
k
1
q
e
T
q
˙
e
+
w
e
T
I
w
˙
e
=
2
k
1
q
e
T
q
˙
e
+
w
e
T
I
w
˙
-
w
˙
d
.Then substituting (1) and (6) into (8), we obtain(30)
V
˙
=
k
1
q
e
T
q
e
x
+
q
e
0
E
3
×
3
w
e
+
w
e
T
T
c
+
T
e
x
t
-
w
b
x
I
t
w
b
+
A
s
I
w
s
Ω
-
I
t
w
˙
d
=
w
e
T
k
1
q
e
0
E
3
×
3
-
k
1
q
e
x
q
e
+
w
e
T
T
c
+
T
e
x
t
-
w
b
x
I
t
w
b
+
A
s
I
w
s
Ω
-
I
t
w
˙
d
.Sinceq
e
x
q
e
=
0, (9) becomes(31)
V
˙
=
w
e
T
k
1
q
e
0
q
e
+
T
c
+
T
e
x
t
-
w
b
x
I
t
w
b
+
A
s
I
w
s
Ω
-
I
t
w
˙
d
.Let control torqueT
c be(32)
T
c
=
-
K
p
q
e
-
K
d
w
e
-
T
e
x
t
+
w
b
x
I
t
w
b
+
A
s
I
w
s
Ω
+
I
t
w
˙
d
,where K
p
=
k
1
q
e
0 and K
d is a symmetric positive matrix.Using the control law obtains(33)
V
˙
=
-
w
e
T
K
d
w
e
.So the stability of the satellite system is ensured. This completes the proof.The quaternion errors and angular speed errors are shown in Figures8 and 9.Figure 8
Time response of quaternion errors.Figure 9
Time response of angular speed errors.The results show that the errors between target and actual maneuvering quaternion are less than4
×
1
0
-
4. The errors between the target and the actual maneuvering angular speeds are less than 1
0
-
3
r
a
d
/
s in all directions. Therefore, the path planned by the Gauss pseudospectral method has high solution precision and is extremely quick.


## 7. Conclusions

In this paper, an attitude maneuver optimization method for gimbal rates is developed based on the Gauss pseudospectral method. It satisfies the control requirement of rapid maneuver for satellite and evades the singularity problem of SGCMGs. This method has higher solution precision and a fast convergence rate, and it avoids the difficulties in solving a boundary value problem. According to the definition of rapid maneuvering where the attitude rate of a satellite is 1~10°/s and the planning angular speed is 5.7°/s. Then, through mathematic simulations based on MATLAB/Simulink, the errors of quaternion and angular speed are obtained within4
×
1
0
-
4 and 1
0
-
3
r
a
d
/
s. In conclusion, the Gauss pseudospectral method can effectively plan an optimal trajectory, satisfying all constraints within a short time.



---

*Source: 1081267-2016-02-24.xml*