- Email: [email protected]
Line-of-sight nonlinear model predictive control for autonomous rendezvous in elliptical orbit
JID:AESCTE AID:4086 /FLA
[m5G; v1.218; Prn:30/06/2017; 16:23] P.1 (1-8)
Aerospace Science and Technology ••• (••••) •••–•••
1
Contents lists available at ScienceDirect
2
67 68
3
Aerospace Science and Technology
4 5
69 70 71
6
72
www.elsevier.com/locate/aescte
7
73
8
74
9
75
10
76
11 12 13 14 15
Line-of-sight nonlinear model predictive control for autonomous rendezvous in elliptical orbit
82
Department of Mechanical Engineering, York University, 4700 Keele Street, Toronto, Ontario M3J 1P3, Canada
22 23 24 25 26 27 28 29 30
83 84
19 21
79 81
18 20
78 80
Peng Li, Zheng H. Zhu ∗
16 17
77
85
a r t i c l e
i n f o
Article history: Received 25 August 2016 Received in revised form 19 June 2017 Accepted 23 June 2017 Available online xxxx Keywords: Spacecraft autonomous rendezvous Line-of-Sight Nonlinear model predictive control Quadratic programming Measurement uncertainties
31 32
a b s t r a c t
86 87
The paper investigates the trajectory planning and control of autonomous spacecraft rendezvous in the orbital plane with line-of-sight dynamics. The control problem, based on nonlinear model predictive control, is formulated in terms of line-of-sight range and azimuth angle. The state feedback with measurement uncertainties is introduced to form a closed-loop optimal control problem by integration of receding horizon strategy. Furthermore, the control input increment instead of total control input is considered in the cost function to generate a smooth transient response. The formulated nonlinear optimal control problem is then transformed into convex quadratic programming problems over the predictive horizon, leading to a computationally efficient algorithm implementable for spacecraft. The numerical results show that the newly proposed line-of-sight nonlinear model predictive control scheme is able to effectively generate optimized approach trajectories with satisfactory control accuracy and the proposed method is insensitive to the measurement uncertainties. © 2017 Elsevier Masson SAS. All rights reserved.
88 89 90 91 92 93 94 95 96 97 98
33
99
34
100
35
101
36 37
1. Introduction
38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59
The increasing population of space debris, such as nonfunctional satellites or spent upper stages of rocket, poses serious threats to space flight missions [1,2]. To address the threats, active space debris removal has emerged as an appealing strategy for the sustainable use of outer space [3,4]. In such missions, the chaser spacecraft (chaser) has to track the motion of space debris, or a non-cooperative target in general, and then approaches and rendezvouses with the target. Challenges arise in the trajectory planning and tracking control of the chaser, where the chaser is usually required to perform complicated and skewed maneuvers for rapid tracking and station keeping, especially in the proximity operation, to deal with the non-cooperative target. Many control methodologies and/or strategies have been devoted to generate optimal approaching trajectories to autonomously transfer from one elliptical orbit to another with various objectives, such as, efficient fuel consumption [5], shortest approaching time [6], high control accuracy or robustness [7], subject to operational constraints. Notably, Nonlinear Optimal Control (NOC) has been recognized as one of the most attractive methods to deal with the constrained optimization problems since it optimizes a specific cost function while satisfying the nonlinear equality and/or
60 61 62 63 64 65 66
*
Corresponding author. E-mail address: [email protected] (Z.H. Zhu).
http://dx.doi.org/10.1016/j.ast.2017.06.030 1270-9638/© 2017 Elsevier Masson SAS. All rights reserved.
inequality constraints [8,9]. However, to obtain a feasible solution to the closed-loop NOC in a fast manner is challenging, even for an unconstrained case [10]. Alternatively, Nonlinear Model Predictive Control (NMPC), which is based on receding horizon strategy (RHS) and re-planning of optimal trajectory in real time by solving the NOC at each sampling instant [11,12], has been proved as an effective method [13–16]. The resulting NOC problem can be further reduced to a quadratic programming (QP) problem that is computationally affordable for computers on-board spacecraft. Autonomous rendezvous with a target in near-field generally employs laser imaging detection and ranging system and advance video guidance system for relative navigation in terms of the relative distance and the Line-Of-Sight (LOS) angles with respect to the target. The aforementioned approaches were based on the LocalVertical Local-Horizontal (LVLH) formulation [5–7], including the LVLH based NMPC approach [13–16]. As a result, the relative navigation information has to be transformed from the LOS frame to the LVLH frame. The extra transformation between the LOS and LVLH frames complicates the derivation of guidance control and adds extra computational efforts for on-board computers. To reduce the computational requirement for the on-board computers, LOS based autonomous rendezvous were developed to employ the navigation directly [17,18]. The control law was developed based on the phase plane analysis. No attempt has been made to the LOS based NMPC in autonomous spacecraft rendezvous, to the best knowledge of authors.
102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132
JID:AESCTE
AID:4086 /FLA
[m5G; v1.218; Prn:30/06/2017; 16:23] P.2 (1-8)
P. Li, Z.H. Zhu / Aerospace Science and Technology ••• (••••) •••–•••
2
1 2
67
Nomenclature
68
3 4 5 6 7 8 9 10 11
69
a A( X ) Ak B Bk e f Fρ , Fθ
12 13 14 15
mc Nc , N P, Q
16 17
RT
semi-major axis of the elliptical orbit state-dependent system matrix discretized system matrix at time step k control input matrix discretized control matrix at time step k eccentricity of the elliptical orbit true anomaly of the elliptical orbit control force components along the LOS range and azimuth angle directions mass of the chaser numbers of control horizon and predictive horizon weight matrices on the control correction and control error distance from the center of the Earth to the target
Ts U max U x, y Xd U k
δ ρ , δθ θ, θd ρ , ρd
ρS σρ , σθ ω
sampling time interval maximum control force control input vector coordinates of the relative distance vector under orbital frame desired state vector incremental control variable measurement errors of LOS range and azimuth angle azimuth angle and desired azimuth angle LOS range and desired LOS range safety distance from the target standard deviations of Gaussian distributions in terms of LOS range and azimuth angle first order time derivative of the true anomaly
18
where R T = a(1 − e 2 )/(1 + e cos f ) is the distance from the center of the Earth to the target, f denotes the true anomaly, e is the eccentricity, a is the semi-major axis of the target, ω = δ˙ = μa(1 − e2 )/ R 2T and ω˙ = −2μe sin f / R 3T are the first and second order time derivatives of the true anomaly respectively, μ is the gravitational constant of the Earth, mc is the mass of the chaser, and F ρ and F θ are the force components in the ρ and θ directions of the LOS frame. The detailed derivation of Eq. (1) is given in Appendix. Introduce the new state vector X = {x1 , x2 , x3 , x4 } T with x1 = ρ , x2 = θ , x3 = ρ˙ and x4 = ρ θ˙ . Then, Eq. (1) is reduced to a set of first-order differential equations, such that,
20 21 22 23 24 25 26 27 28 29 30 31 32
⎡
34 35
Fig. 1. Schematic of the LOS frame.
37 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66
72 73 74 75 76 77 78 79 80 81 82 83
In the current work, the autonomous rendezvous is formulated in the LOS frame to directly use the relative navigation information. This intuitive description of the trajectory planning will result in concise formulation of the system dynamics and the controller design, especially for autonomous rendezvous with a passive noncooperative target. 2. Line-of-sight formulation of spacecraft dynamics Consider a chaser approaching a target in the orbital plane of an elliptical orbit, as shown in Fig. 1. To focus on the fundamentals of the LOS based NMPC, the current work is limited to the in-plane rendezvous with the target. The orbital frame, shown as O o − xo y o zo in Fig. 1, is defined with its origin at the Center of Mass (CM) of the target, where the y o -axis is along the orbital radius of the target, the xo -axis lies in orbital plane and is perpendicular to the y o -axis, and the zo -axis is normal to the orbital plane to complete a right-hand system. The LOS frame is formed by the range ρ and azimuth angle θ with its origin at the CM of the target, where the azimuth angle θ is measured from the xo -axis in the orbital plane as shown in Fig. 1. Accordingly, the equations of motion of the chaser in the LOS frame is expressed as per [19],
ρ¨ − ρ θ˙ 2 + 2ωρ θ˙ +
1 − 3 sin2 θ
μ R 3T
− ω2 ρ =
Fρ mc
μ Fθ ρ θ¨ + 2ρ˙ θ˙ − 2ωρ˙ − ω˙ + 3 3 sin θ cos θ ρ = RT
mc
(1a) (1b)
0
0
1
⎢ 0 0 0 ⎢ ⎢ μ A( X ) = ⎢ 2 2 0 ⎢ ω − R 3 (1 − 3 sin x2 ) 0 T ⎣ μ ω˙ + 3 R 3 sin x2 cos x2 0 2(ω − T ⎤ ⎡ and
85 86 87 88 89 90 91 92 93 94 95 96 97 98
˙ (t ) = A ( X ) X (t ) + BU (t ) X
33
38
71
84
19
36
70
0 0 ⎢0 0⎥ ⎥ B =⎢ ⎣1 0⎦ 0 1
0
x4 x1 x4 ) x1
⎤
(2)
99 100
⎥ 1 ⎥ x1 ⎥ − 2ω ⎥ ⎥ ⎦
101 102 103 104
0
105 106 107
(3)
where A ( X ) ∈ 4×4 and B ∈ 4×2 are the state-dependent system matrix and control input matrix, respectively, and U = { F x1 , F x2 } T /mc is the control input. It should be noted that the state x1 = ρ will never approach to zero in practice because the target has finite dimensions. The control objective of the trajectory planning for an autonomous spacecraft rendezvous is to approach the target with a desired state, X d = {x1d , x2d , 0, 0} T , smoothly to prevent sudden accelerations or decelerations of the chaser.
108 109 110 111 112 113 114 115 116 117 118 119 120 121
3. Nonlinear model predictive control
122 123
The control objective inevitably leads to an optimal control problem subject to constraint of thrust magnitude. The NMPC, also known as the Nonlinear Receding Horizon Control (NRHC), has been widely used to solve the constrained optimal control problem, due to its advantage of online generation of a set of feedback control commands by iteratively solving an open-loop NOC problem at each sampling instant. The receding horizon process is repeated by shifting the time one-step forward each time [20]. Accordingly, the optimal control problem of the trajectory planning
124 125 126 127 128 129 130 131 132
JID:AESCTE AID:4086 /FLA
[m5G; v1.218; Prn:30/06/2017; 16:23] P.3 (1-8)
P. Li, Z.H. Zhu / Aerospace Science and Technology ••• (••••) •••–•••
1 2 3 4 5 6 7 8 9 10 11 12 13 14
and tracking of spacecraft rendezvous can be formulated as a series of continuous-time nonlinear optimal control (CNOC) problems by NMPC at each sampling instant. Denote the starting and ending times of the rendezvous as t 1 and t f . Divide the total time interval into n subintervals evenly, such that, [t 1 , t 2 , · · · , tn+1 ] where tn+1 = t f . Furthermore, assume predictive (T p ) and control (T c ) horizons are constant with T p ≥ T c . Then, the receding predictive horizon [tk , tk + T p ] of the kth CNOC problem starts from tk ∈ [t 1 , tn+1 ] with the following definition. Problem Ψ C N O C : At any time instant tk ∈ [t 1 , tn+1 ], find the control input rate U˙ (t ) ∈ 2 and states X (t ) ∈ 4 that minimize the following quadratic cost function over a given predictive horizon T p ,
15
tk + T p
16 17 18
min J tk =
Q (t ) X (t ) − X d dt
tk+ T c
20 21
+
22
67
⎧ U k,k = U k + U k−1 ⎪ ⎪ ⎪ U k+1,k = (U k+1 + U k ) + U k−1 ⎪ ⎪ ⎪ ⎪ .. ⎪ ⎪ ⎪ . ⎪ ⎨ N c −1 U = U k + i + U k −1 ⎪ k + N − 1 , k c ⎪ ⎪ ⎪ ⎪ i = 0 ⎪ ⎪ .. ⎪ ⎪ ⎪ ⎪ ⎩. U k+i ,k = U k+ N c −1,k (i = N c , . . . , N )
73
U˙ (t ) T P (t )U˙ (t )dt
(4)
69 70 71 72
75 76 77
(10)
78 79 80 81 82 83 84
(11)
where U k+i ,k is the predicted control input at the (k + i)th time instant with respect to the kth time instant, U k−1 is the known N −1 control input at the kth time instant, U˜ k,k = i =c 0 U k+i ,k |i and
U˜ k =
68
74
˜ U˜ k + F˜ U k−1 U˜ k,k = M
tk
23
The resulting DNOC problem is computationally heavy and is transformed into a QP problem whose solution can be obtained as follows. Define the recursive relationship of the control input at any time instant with respect to the kth time instant in the control horizon as
or in a compact notation
tk
19
24
X (t ) − X d
T
3
N c −1
85 86 87 88 89 90
28
and the box constraint
˜ is U k+i |i (|· is the Ket vector notation), M ˜ (i , j ) = I 4 , (i ≥ j ), and F˜ = a lower triangular matrix with M Nc i =1 I 4 |i where I 4 is the identity matrix. Correspondingly, the set of states X k+ j ,k , ( j = 1, . . . , N ), in the
29
U j (t ) ≤ U max ∈
predictive horizon with respect to the kth time instant can be evaluated recursively, such that,
95
⎧ X k+1,k = A k X k + B k U k,k ⎪ ⎪ ⎪ X k+2,k = ( A k+1 A k ) X k + ( A k+1 B k )U k,k + B k+1 U k+1,k ⎪ ⎪ ⎪ ⎨. .. N −1 ⎪ N −1 N −1 ⎪ ⎪ ⎪ ⎪ A k+ j X k + A k+i ) B k+ j U k+ j ,k ( ⎪ ⎩ X k+ N ,k =
97
subject to the system dynamics
25 26 27
30 31 32 33 34 35 36 37 38 39 40 41 42 43
˙ (t ) = A ( X ) X (t ) + BU (t ) X
46
(6)
where U j (t ) is the jth element of the control input U (t ), U max is the maximum control force available, and ( Q (t ) ∈ 4×4 , P (t ) ∈ 2×2 ) are the time-varying and positive definite symmetric weight matrices, respectively. The solution of the kth CNOC problem is used as the initial condition for the (k + 1)th CNOC problem. To solve the CNOC problem, it is discretized into the discrete nonlinear optimal control (DNOC) problem over the given receding predictive horizon T p as follows. Problem Ψ DNOC : At any time instant k ∈ [1, n + 1], find the sequences of incremental control impulse U k = U k − U k−1 and state vectors X k that minimize the quadratic cost function over the given predictive horizon,
44 45
(5)
min J k =
k +N
{ X i − X d }T Q k { X i − X d } +
i =k+1
47
k+ N c −1 i =k
(7)
48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66
U iT P k U i
subject to the system dynamics
X i +1 = A i X i + B i U i ,
( i = k + 1, . . . , k + N )
(8)
and the box constraint
|U j |k ≤ U
max
(9)
where the corresponding predictive and control horizons T p = N T s and T c = N c T s in Eq. (4) are placed by N and N c with the assumption of N c ≤ N. Note that T s is the sampling time interval. Moreover, Eq. (8) is the discrete representation of Eq. (2) by zero-order holder, where the coefficient matrices A i and B i are assumed constant within each sampling time interval T s , but are updated at the beginning of each time interval with the states at the end of previous time interval. The value of weight matrices Q k and P k will be determined based on mission requirements or rendezvous strategy.
i =0
j =0
j =0
(12)
min J k =
− X˜ d
Q˜
˜ kp X
− X˜ d
T + U˜ k
(14)
⊕iN=c 1 P i
T
T
and P i , (i = 1, 2, ..., N c ), and thus Q˜ = Q˜ and P˜ = P˜ . ˜ d − φ X˜ kp − Defining and substituting an auxiliary vector E = X Γ U k−1 into Eq. (14) yield, [21]
=
1 2
T U˜ k
T
H U˜ k +
+ E T Q˜ E
105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 121 122 123 124 125 126 127 128 129
T Q˜ G U˜ k − E + U˜ k P˜ U˜ k T C T U˜ k
101
120
P˜ U˜ k
˜ d = i =1 X d |i , Q˜ = ⊕ N Q i and P˜ = where X (⊕ is the i =1 direct sum notation). It is worth noting that the weight matrices Q˜ and P˜ inherit the symmetric property from Q i , (i = 1, 2, ..., N)
100
104
N
min J k = G U˜ k − E
96
103
i =1 + j
vector is composed of the predicted states and the detailed expression of the coefficient matrices φ , Γ and G are given in Ref. [21]. Notably, the coefficient matrices φ , Γ and G are constant within each time step but are updated at the beginning of each step with the states at the end of previous step. For the sake of derivation convenience, the subscript k is omitted and the same is applied to the matrices E , H and C to be discussed later. Furthermore, the cost function in Eq. (7) is rewritten in a compact form,
94
102
p
T
93
99
˜ k = φ X k + Γ U k−1 + G U˜ k X (13) N p ˜ k = i =1 X k+i ,k |i and the superscript p indicates the where X
˜ kp X
92
98
Substituting Eq. (11) into Eq. (12) yields the set of predicted states X k+ j ,k , ( j = 1, . . . , N), in terms of U˜ k,k and in the compact notation as,
91
130 131
(15)
132
JID:AESCTE
AID:4086 /FLA
[m5G; v1.218; Prn:30/06/2017; 16:23] P.4 (1-8)
P. Li, Z.H. Zhu / Aerospace Science and Technology ••• (••••) •••–•••
4
1
67
2
68
3
69
4
70
5
71
6
72
7
73
8
74
9
75
10
76
11
77
12
78
13
79
14
80
15
81
16
82
17
83
18
84
Fig. 2. Block diagram of LOS based NMPC.
19 20 21 22 23 24 25 26 27 28
where H = 2[ G T Q˜ G + P˜ ] is the Hessian matrix that is positive definite and symmetric, C = −2G T Q˜ E is a column vector with the same dimension of U˜ k . Note that the vector E is constant within each time interval. Thus, the term E T Q˜ E in Eq. (15) is constant and does not affect the minimization of the cost function. It can be safely ignored from the cost function. Accordingly, Eq. (15) is equivalent to the following standard QP problem. Problem Ψ Q P , find U˜ k that minimizes the following cost function,
29 30 31 32 33 34 35
min J k =
1 2
T
U˜ k H U˜ k + C T U˜ k
(16)
The positive definite matrix H makes it a convex optimization problem where U˜ k is a global minimum. Moreover, the control magnitude constraint is expressed in terms of control increment U˜ k [21],
36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66
max ˜ U˜ k + F˜ U k−1 ≤ U˜ max −U˜ ≤M max
N
(17)
where U˜ = i =c 1 U max |i . Thus, Eq. (13) and Eq. (17) are the system dynamics and the control magnitude constraint, respectively. Compared our approach in Eqs. (1)–(17) with the approach based on the LVLH frame, the only differences are (i) the transformation of relative navigation information from LOS to LVLH frame is eliminated, (ii) the formation of state dependent system matrix. It is obvious the first difference results in saving of computational effort within each time step. For the second difference, the state dependent system matrix in the current approach is sparser than the matrix used in the common method, see the Ref. [22] in the revised manuscript. Thus, the associated computational effort will be lighter compared with that of Ref. [22]. Finally, measurement errors (δ ρ , δθ ) from navigation sensing systems are inevitable, such as the optical artifacts (e.g., glints, glares, saturations, hot pixels) in the sensors’ fields of view, which may result in sudden changes in control input. Theoretically, the sudden changes in control input could be indirectly rejected by smoothing the control input increment, which is the objective of the proposed optimal control cost function. Thus, the effectiveness of disturbance rejection will be demonstrated by adding the errors δ X = [δ ρ , δθ, 0, 0] T to the control input U k and state X k+1 at the beginning of each time instant. The navigation errors are assumed as Gaussian white noise with zero mean. The problem is then solved by the interior-point method due to its high convergence rate and the ease of implementation and the computational diagram of the close-loop control strategy is given as shown in Fig. 2.
85 86
Table 1 Parameters of weight matrices.
87 88
Parameter
Value
89
State weight matrix Control weight matrix
Q = diag(50, 5000, 1000, 1/ρ ) P = diag(100, 100)
90 91 92
4. Rendezvous strategy
93
The rendezvous of the chaser with the target is achieved by a judicious two-phased strategy to ensure the final azimuth direction at rendezvous to be satisfied precisely. In the first phase, the control task is to maneuver the chaser from its initial position to an intermediate point in the vicinity of the target for collision avoidance, such that, ρ > ρ S and ρ S is the safety distance from the target determined by the mission requirement. At the same, the chaser must adjust its azimuth angle from the initial condition θ0 to the desired value, θ = θd . In the second phase, the control task is to approach the target in a straight trajectory along the desired azimuth angle until ρ = ρd to accomplish the rendezvous mission. In the current work, it is assumed that ρ S = 20 m, θ0 and θd are determined by the particular rendezvous scenarios as specified in the following numerical examples. The above two-phased rendezvous strategy is achieved by the properly selected weight matrices, see in Table 1. Assume the chaser’s ability to maneuver in LOS range and angle directions is equally important in the rendezvous mission. This leads to the same penalties on the two control input (F ρ , F θ ) in the weight matrix P . The determination of weight matrix Q is heavily dependent on the rendezvous strategy. The control task in the first phase requires the chaser to aim the target in the required azimuth angle with a safe distance away. In this phase, the control accuracy for the azimuth angle is less important when the chaser is far away from the target. Accordingly, the penalty on the azimuth angle θ is set much larger than that on the LOS range ρ in this phase. In the second phase, it is critical to maintain a smooth trajectory in the proximity operation. Thus, the penalty for the third state (the approaching speed) is designed to slow down the relative velocity of the chaser for a safe approach. The last diagonal element in the weight matrix Q is the penalty for the fourth state – azimuth angular rate and is designed as state dependent, 1/ρ , instead of constant. Different from other three penalties, this penalty weight increases dramatically as the chaser approaches the target. Thus, the azimuth angular rate is reduced quickly to avoid the chattering of approaching trajectory because a stable azimuth angle at the final approach is mission critical. Again, it should be noted that the LOS range ρ does not decrease to zero in real situation even when
95
94 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132
JID:AESCTE AID:4086 /FLA
[m5G; v1.218; Prn:30/06/2017; 16:23] P.5 (1-8)
P. Li, Z.H. Zhu / Aerospace Science and Technology ••• (••••) •••–•••
1 2 3
Parameters
5 7 8 9
67
Table 2 Rendezvous conditions.
4 6
5
LOS range (ρ ), m Azimuth angle (θ ), ◦ Relative velocity (ρ˙ ), m/s Azimuth angle rate (θ˙ ), ◦ /s
68
Case 1
69
Case 2
Initial state
Desired state
Initial state
Desired state
70
80 0 0 0
1 0, 90, 180, 270 0 0
80 0, −90, −180, −270 0 0
1
72
−90
73
0 0
74
10 11 12 13 14 15
77
Table 3 Rendezvous navigation sensor uncertainties [23]. Measurands
Standard deviation (σ )
LOS range (ρ ), cm Azimuth angle (θ ), ◦
1.518 0.002787
78 79 80 81 82
17 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46
83
the rendezvous is completed. This is because the LOS frame is centered at the CM of the target and there is a minimum value of ρ measured from the CM to external surface of the target.
49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66
84 85 86 87
5. Case studies Fig. 3. In-plane approaching trajectories for Case 1.
The advantages of the proposed LOS NMPC in the spacecraft rendezvous are demonstrated by dynamic simulation of two rendezvous scenarios. The rendezvous conditions are listed in Table 2. The chaser is assumed to navigate by optical sensing system with direct measurement of the range and azimuth angle to the target. The sensor uncertainties, listed in Table 3 as per [23], are taken as a disturbance to the state, where σρ and σθ are the standard deviations of the Gaussian distributions in terms of LOS range ρ and the azimuth angle θ , respectively. The target is assumed moving in an elliptic orbit with a semimajor axis of 10,000 km and a true anomaly of 0◦ . The initial position of the target is assumed at perigee with an initial relative distance between the chaser and the target spacecraft of 80 m. To explore the application limit of the proposed approach, the eccentricity of 0.3 is used in the case studies. The predictive horizon and control horizons are defined as N = 10 and N c = 1. The mass of the chaser is assumed as 500 kg and the maximum available control thrust is limited to U max = 1 N. The simulation duration is set to 600 s and the sampling time is T s = 0.1 s. Note that the simulation results are shown in the LOS frame and the relative states are transformed to the orbital frame to show the approaching trajectory. 5.1. Approach target in different directions
47 48
75 76
16 18
71
In this case, the coplanar rendezvous with the target from various directions is studied. At t = 0, the chaser is assumed to running in the same orbit as the target, but 80 m ahead with 0◦ azimuth angle. The control objective is to approach the target along x0 -axis and y 0 -axis directions by controlling the azimuth angle θ equal to 0◦ , 90◦ , 180◦ and 270◦ , respectively. For the purpose of collision avoidance and subsequent proximity operations, the chaser is required to stop at 1 m away from the target in all directions at the end of approach. The detailed rendezvous conditions are listed as the Case 1 in Table 2. The simulation results of Case 1 are shown in Figs. 3–5. Fig. 3 shows the coplanar approaching trajectories in the orbital frame along x0 -axis and y 0 -axis directions, respectively. The chaser maneuvers around the target first to align the azimuth angle to the desired ones and then approaches to the desired position nearly in a straight line in the final stage as required, leading to a successful rendezvous. Next, Fig. 4 displays the time histories of the LOS range, relative velocity, azimuth angle and azimuth angle rate, respectively. Since
the penalty term 1/ρ is relatively small at the beginning, the control on the azimuth angle rate is weak and the angle approaches to the desired value at a high rate initially. As the LOS range reduces and the penalty term 1/ρ increases, the azimuth angle rate is reduced and the change of azimuth is flattened, see Fig. 4(c). As a result, the chaser adjusts its trajectory with respect to the target first and then approaches in a straight line. This approaching mode is widely adopted in space rendezvous missions to avoid collision with the target [6,?]. The azimuth angle reaches its desired value at 100 s while the chaser reaches its final position at 600 s. The proposed control scheme yields smooth trajectories in all cases, which is evident in Fig. 4. It shows there is nearly no overshoot in the time histories of the LOS range and the azimuth angle. To further demonstrate the effectiveness of the proposed control law, the control force profiles in the LOS range and azimuth angle channels are plotted in Fig. 5. It is interesting to note that the control force F ρ in the LOS range channel reaches its maximum output magnitude, −1 N, in the first tens of seconds, see Fig. 5(a), in order to quickly decrease the LOS range to reduce the rendezvous time. The negative symbol means the corresponding thrust is against the motion. The control force F θ in the azimuth angle channel converges to zero at around 100 s after the azimuth angle θ reaches its desired value, see Fig. 5(b), reflecting the effectiveness of the proposed two-phased rendezvous strategy. Furthermore, the chattering in the control force profiles is very small, which indicates that the proposed control law is effective in rejecting the measurement uncertainties to reach the desired values of the LOS range and azimuth angle smoothly, see Figs. 3–4. 5.2. Approach target from different initial azimuth angles
88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120
In this case, the coplanar rendezvous with the target is studied by assuming the chaser starting from different initial azimuth angles but with the same distance, 80 m, from the target. The control objective is to drive the chaser approaching the target in the − y 0 -axis direction from 4 different initial azimuth angles θ0 : 0◦ , −90◦ , −180◦ and −270◦ , respectively. Similarly to Case 1, the chaser is required to stop at 1 m away from the target in the − y 0 -axis direction at the end of rendezvous for subsequent proximity operations. The detailed rendezvous conditions are listed as the Case 2 in Table 2. The simulation results are shown in Figs. 6–8. Fig. 6 plots the rendezvous trajectories in the orbital frame. In all cases, the chaser
121 122 123 124 125 126 127 128 129 130 131 132
JID:AESCTE
6
AID:4086 /FLA
[m5G; v1.218; Prn:30/06/2017; 16:23] P.6 (1-8)
P. Li, Z.H. Zhu / Aerospace Science and Technology ••• (••••) •••–•••
1
67
2
68
3
69
4
70
5
71
6
72
7
73
8
74
9
75
10
76
11
77
12
78
13
79
14
80
15
81
16
82
17
83
18
84
19
85
20
86
21
87
22
88
23
89
Fig. 5. Control input: (a) LOS range and (b) azimuth angle.
24
90
25
91
26
92
27
93
28
94
29
95
30
96
31
97
32
98
33
99
34
100
35
101
36
102
Fig. 6. In-plane approaching trajectories for Case 2.
37 38
104
39 40 41 42 43 44 45 46 47
Fig. 4. Case 1: time histories of (a) LOS range (b) relative velocity (c) azimuth angle and (d) azimuth angle rate.
48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66
103
successfully maneuvers to the vicinity of the target in the − y 0 -axis direction first and then decreases the LOS range to the designated position in a straight line as required. Fig. 7 shows the responses of the LOS range, relative velocity, azimuth angle and azimuth angle rate, respectively. Similarly, the azimuth angle rate changes quickly in the first 70 s since the chaser is far away from the target and the penalty term 1/ρ is relatively small. As the chaser getting closer to the target, the penalty 1/ρ increases hyperbolically. Accordingly, the azimuth angle rate approaches zero quickly which implies the azimuth angle reaches its desired value. The chaser adjusts its orientation with respect to the target prior to approach a straight trajectory in the final approach as per the rendezvous strategy. The trajectories are smooth as required. Finally, Fig. 8 plots the control forces in the LOS range and azimuth angle channels, which show the same trends as those in Fig. 5. The simulations were conducted based on Eq. (2) and Eq. (A.1) under the simulation scenarios 1 and 2, respectively. The computational times were 134.26 s and 136.11 s under the environment
of i7-4510 2.6 GHz central processing unit and 16.0 GB RAM. The computational time of LOS model is slightly less as expected. It should be noted that the direct outputs of relative navigation system are the LOS range and the azimuth angle. Our comparison did not include the time delay due to the transformation of navigation information between the LOS and orbital frames, which depends on individual converter and circuitry. Therefore, the current approach is advantageous in saving computational time if one considers the possibility of eliminating navigation information transformation process under LOS frame. It is also worthy pointing out that the control consumption and state evolution are almost identical in the above cases for difference orbital eccentricities ranging from 0 to 0.3 although only the results of e = 0.3 are reported. This is because the trajectory of the target is almost straight during the short rendezvous period no matter what the eccentricity is, for the given initial condition of rendezvous, the rendezvous location and elliptical orbital parameters ( f = 0◦ and a = 10,000 km) in the current study. This is an interesting finding that indicates the eccentricity of target orbit has little impact on the rendezvous due to the short rendezvous period and the short distance between the target and the chaser compared with the distance from the target to the center of Earth.
105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127
6. Conclusions
128 129
This paper develops a LOS based NMPC for the coplanar autonomous rendezvous with non-cooperative targets. The LOS formulation simplifies the equations of relative dynamics of space-
130 131 132
JID:AESCTE AID:4086 /FLA
[m5G; v1.218; Prn:30/06/2017; 16:23] P.7 (1-8)
P. Li, Z.H. Zhu / Aerospace Science and Technology ••• (••••) •••–•••
7
1
67
2
68
3
69
4
70
5
71
6
72
7
73
8
74
9
75
10
76
11
77
12
78
13
79
14
80
15
81
16
82
17
83
18
84
19
85
20
86
21
87
22
88
23
89
24
90
25
91
Fig. 8. Control input: (a) LOS range and (b) azimuth angle.
26
92
27
93
Acknowledgements
28 29
94
This work is supported by the Discovery Grant and Discovery Accelerator Supplement Grant of the Natural Sciences and Engineering Research Council of Canada (NSERC).
30 31 32 33
37 38 39 40 41
42 43 44
x¨ y¨
=
45 47
Fig. 7. Case 2: time histories of (a) LOS range (b) relative velocity (c) azimuth angle and (d) azimuth angle rate.
50 51 53 54 55 56 57 58 59 60 61
craft using the relative navigation information from the measurement directly. The NMPC model is derived into a state dependent representation and the control is optimized for the control error and control smoothness simultaneously to achieve a smooth transient response. Furthermore, the proposed methodology is robust for various eccentricities of the target’s orbit because the trajectory of the target is almost straight during the short rendezvous period. Numerical simulations demonstrate the newly proposed method is effective and capable of achieving the control objective and it is insensitive to navigation measurement errors.
62
Conflict of interest statement
65 66
None declared.
1 mc
−2ω
0 Fx Fy
x˙ y˙
+
ω2 −ω˙ ω˙ ω2
x y
+
μ R 3T
−x
2y
102 103 104 105 106 107 108 109 110 111
(A.1)
112 113
where {x, y } T , {˙x, y˙ } T are the relative position and velocity between the chaser and the target spacecraft, ω is first order time derivative of the true anomaly, R T is the distance from the center of the Earth to the target, and { F x , F y } T is the applied force. The relationship between the orbital frame and the LOS frame can be expressed as,
114
120
x y
=
cos θ 0
0 sin θ
ρ
(A.2)
115 116 117 118 119 121 122 123
Accordingly, the first and second order derivatives of the relative states {x, y } T can be expressed in terms of LOS frame, such that,
63 64
0 2ω
+
46
52
98
101
The derivation of LOS model is based on Ref. [19]. Assume the distance between two spacecraft is sufficiently small compared to the orbit radius but sufficiently large compared to the largest dimension of either the target or the chaser. Then, the relative motion of spacecraft can be described by the linearized Tschauner– Hempel equation in the orbital frame, such that,
36
49
97
100
35
48
96
99
Appendix A. Derivation of line-of-sight model
34
95
x˙ = ρ˙ cos θ − ρ θ˙ sin θ y˙ = ρ˙ sin θ + ρ θ˙ cos θ
x¨ = ρ¨ − ρ θ˙ 2 cos θ − ρ θ¨ + 2ρ˙ θ˙ sin θ y¨ = ρ¨ − ρ θ˙ 2 sin θ + ρ θ¨ + 2ρ˙ θ˙ cos θ
Substituting Eq. (A.3) and Eq. (A.4) into Eq. (A.1) yields,
124 125 126
(A.3)
127 128 129
(A.4)
130 131 132
JID:AESCTE
AID:4086 /FLA
[m5G; v1.218; Prn:30/06/2017; 16:23] P.8 (1-8)
P. Li, Z.H. Zhu / Aerospace Science and Technology ••• (••••) •••–•••
8
1 2 3 4 5 6 7 8 9
⎧ Fρ μ μ ⎪ ⎪ ρ¨ − ρ (θ˙ − ω)2 + 3 ρ − 3 1 − 3 sin2 θ sin θ = ⎪ ⎨ mc RC RT 3 ⎪ R μ F θ ⎪ ⎪ ρ θ¨ + 2ρ˙ (θ˙ − ω) − ωρ ˙ − 2 1 − 3T cos θ = ⎩ RT
(A.5)
mc
RC
where { F ρ , F θ } T is the applied force expressed in the LOS frame and can be transformed from the orbital frame as,
10 11
Fρ Fθ
= A −1
Fx Fy
,
A −1 =
and
cos θ − sin θ
sin θ cos θ
12
Equation (A.5) can be further linearized by approximating the term
13
R 3T
14
R 3C
with first order Taylor expansion, such as,
15 16 17 18
R 3T R 3C
19
22 23 24 25 26 27
≈ 1+ ≈1−
20 21
30 31 32 33 34 35 36 37 38 39 40 41 42 43
RT
− 32
≈1−
3 2
·
2y RT
+
15 8
·
2y RT
2 + ···
3y
(A.6)
RT
and,
⎧ μ Fρ ⎪ 2 ⎪ ˙ 2 + 2ωρ θ˙ + 1 − 3 sin2 θ ¨ ρ − ρ − ω ρ= θ ⎪ ⎨ 3 mc RT ⎪ μ F θ ⎪ ⎪ ρ θ¨ + 2ρ˙ θ˙ − 2ωρ˙ − ω˙ + 3 sin θ cos θ ρ = ⎩ 3 RT
28 29
2y
(A.7)
mc
References [1] J.C. Liou, N.L. Johnson, A sensitivity study of the effectiveness of active debris removal in LEO, Acta Astronaut. 64 (2009) 236–243, http://dx.doi.org/10.1016/ j.actaastro.2008.07.009. [2] J.C. Liou, N.L. Johnson, N.M. Hill, Suboptimal power-limited rendezvous with fixed docking direction and collision avoidance, J. Guid. Control Dyn. 66 (2010) 648–653, http://dx.doi.org/10.1016/j.actaastro.2009.08.005. [3] W. Hao, Z.H. Zhu, D. Jin, H. Hu, Model predictive control with output feedback for a deorbiting electrodynamic tether system, J. Guid. Control Dyn. (2016), http://dx.doi.org/10.2514/1.G000535, in press. [4] N.-H. Chau, I. Sharf, Adaptive reactionless motion and parameter identification in postcapture of space debris, J. Guid. Control Dyn. 36 (2013) 404–414, http://dx.doi.org/10.2514/1.57856. [5] P. Lu, X. Liu, Autonomous trajectory planning for rendezvous and proximity operations by conic optimization, J. Guid. Control Dyn. 36 (2013) 375–389, http://dx.doi.org/10.2514/1.58436.
[6] G. Boyarko, O. Yakimenko, M. Romano, Optimal rendezvous trajectories of a controlled spacecraft and a tumbling object, J. Guid. Control Dyn. 34 (2011) 1239–1252, http://dx.doi.org/10.2514/1.47645. [7] F. Gavilan, R. Vazquez, E.F. Camacho, Chance-constrained model predictive control for spacecraft rendezvous with disturbance estimation, Control Eng. Pract. 20 (2012) 111–122, http://dx.doi.org/10.1016/j.conengprac.2011.09.006. [8] H.J. Peng, J. Xin, B.S. Chen, Optimal nonlinear feedback control of spacecraft rendezvous with finite low thrust between libration orbits, Nonlinear Dyn. 76 (2014) 1611–1632, http://dx.doi.org/10.1007/s11071-013-1233-9. [9] Z. Ma, O. Ma, N. Shashikanth, Optimal approach to and alignment with a rotating rigid body for capture, J. Astronaut. Sci. 55 (2007) 407–409, http://dx.doi. org/10.1007/BF03256532. [10] H.G. Zhang, H.L. Yan, D.R. Liu, Neural-network-based near-optimal control for a class of discrete-time affine nonlinear systems with control constraints, IEEE Trans. Neural Netw. 20 (2009) 1490–1503, http://dx.doi.org/10.1109/TNN.2009. 2027233. [11] E.F. Camacho, C. Bordons, Model Predictive Control, Springer-Verlag, London, 2004. [12] J.M. Maciejowski, Predictive Control with Constraints, Pearson Education, Harlow, UK, 2002. [13] E.N. Hartley, P.A. Trodden, A.G. Richards, J.M. Maciejowski, Model predictive control system design and implementation for spacecraft rendezvous, Control Eng. Pract. 20 (2012) 695–713, http://dx.doi.org/10.1016/j.conengprac.2012. 03.009. [14] L. Breger, J.P. How, Gauss’s variational equation-based dynamics and control for formation flying spacecraft, J. Guid. Control Dyn. 30 (2007) 437–448, http://dx.doi.org/10.2514/1.22649. [15] S.D. Cairano, H. Park, I. Kolmanovsky, Model predictive control approach for guidance of spacecraft rendezvous and proximity maneuvering, Int. J. Robust Nonlinear Control 22 (2012) 1398–1427, http://dx.doi.org/10.1002/rnc.2827.
67
[16] M. Leomanni, E. Rogers, S.B. Gabriel, Explicit model predictive control approach for low-thrust spacecraft proximity operations, J. Guid. Control Dyn. 37 (2014) 1780–1790, http://dx.doi.org/10.2514/1.G000477. [17] S.H. Yu, Range-rate control algorithms and space rendezvous schemes, J. Guid. Control Dyn. 20 (1997) 206–208, http://dx.doi.org/10.2514/2.4024. [18] S.H. Yu, Autonomous rendezvous in elliptical orbits, Acta Astronaut. 41 (1997) 95–101, http://dx.doi.org/10.1016/S0094-5765(97)00204-X. [19] T. Liu, J. Zhao, Spacecraft Dynamics, Harbin Institute of Technology Press, Harbin, China, 2003, pp. 86–87 (in Chinese). [20] A. Weiss, M. Baldwin, S.R. Erwin, I. Kolmanovsky, Model predictive control for spacecraft rendezvous and docking: strategies for handling constraints and case studies, IEEE Trans. Control Syst. Technol. 23 (2015) 1638–1646, http://dx.doi. org/10.1109/TCST.2014.2379639. [21] P. Li, Z.H. Zhu, S.A. Meguid, State dependent model predictive control for orbital rendezvous using pulse-width pulse-frequency modulated thrusters, Adv. Space Res. 58 (2016) 64–73, http://dx.doi.org/10.1016/j.asr.2016.04.022. [22] C.D. Karlgaard, Robust rendezvous navigation in elliptical orbit, J. Guid. Control Dyn. 29 (2006) 459–499, http://dx.doi.org/10.2514/1.19148. [23] L.C. Ma, X.Y. Meng, Z.Z. Liu, L.F. Du, Suboptimal power-limited rendezvous with fixed docking direction and collision avoidance, J. Guid. Control Dyn. 36 (2012) 229–239, http://dx.doi.org/10.2514/1.56449.
92
68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109
44
110
45
111
46
112
47
113
48
114
49
115
50
116
51
117
52
118
53
119
54
120
55
121
56
122
57
123
58
124
59
125
60
126
61
127
62
128
63
129
64
130
65
131
66
132
[m5G; v1.218; Prn:30/06/2017; 16:23] P.1 (1-8)
Aerospace Science and Technology ••• (••••) •••–•••
1
Contents lists available at ScienceDirect
2
67 68
3
Aerospace Science and Technology
4 5
69 70 71
6
72
www.elsevier.com/locate/aescte
7
73
8
74
9
75
10
76
11 12 13 14 15
Line-of-sight nonlinear model predictive control for autonomous rendezvous in elliptical orbit
82
Department of Mechanical Engineering, York University, 4700 Keele Street, Toronto, Ontario M3J 1P3, Canada
22 23 24 25 26 27 28 29 30
83 84
19 21
79 81
18 20
78 80
Peng Li, Zheng H. Zhu ∗
16 17
77
85
a r t i c l e
i n f o
Article history: Received 25 August 2016 Received in revised form 19 June 2017 Accepted 23 June 2017 Available online xxxx Keywords: Spacecraft autonomous rendezvous Line-of-Sight Nonlinear model predictive control Quadratic programming Measurement uncertainties
31 32
a b s t r a c t
86 87
The paper investigates the trajectory planning and control of autonomous spacecraft rendezvous in the orbital plane with line-of-sight dynamics. The control problem, based on nonlinear model predictive control, is formulated in terms of line-of-sight range and azimuth angle. The state feedback with measurement uncertainties is introduced to form a closed-loop optimal control problem by integration of receding horizon strategy. Furthermore, the control input increment instead of total control input is considered in the cost function to generate a smooth transient response. The formulated nonlinear optimal control problem is then transformed into convex quadratic programming problems over the predictive horizon, leading to a computationally efficient algorithm implementable for spacecraft. The numerical results show that the newly proposed line-of-sight nonlinear model predictive control scheme is able to effectively generate optimized approach trajectories with satisfactory control accuracy and the proposed method is insensitive to the measurement uncertainties. © 2017 Elsevier Masson SAS. All rights reserved.
88 89 90 91 92 93 94 95 96 97 98
33
99
34
100
35
101
36 37
1. Introduction
38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59
The increasing population of space debris, such as nonfunctional satellites or spent upper stages of rocket, poses serious threats to space flight missions [1,2]. To address the threats, active space debris removal has emerged as an appealing strategy for the sustainable use of outer space [3,4]. In such missions, the chaser spacecraft (chaser) has to track the motion of space debris, or a non-cooperative target in general, and then approaches and rendezvouses with the target. Challenges arise in the trajectory planning and tracking control of the chaser, where the chaser is usually required to perform complicated and skewed maneuvers for rapid tracking and station keeping, especially in the proximity operation, to deal with the non-cooperative target. Many control methodologies and/or strategies have been devoted to generate optimal approaching trajectories to autonomously transfer from one elliptical orbit to another with various objectives, such as, efficient fuel consumption [5], shortest approaching time [6], high control accuracy or robustness [7], subject to operational constraints. Notably, Nonlinear Optimal Control (NOC) has been recognized as one of the most attractive methods to deal with the constrained optimization problems since it optimizes a specific cost function while satisfying the nonlinear equality and/or
60 61 62 63 64 65 66
*
Corresponding author. E-mail address: [email protected] (Z.H. Zhu).
http://dx.doi.org/10.1016/j.ast.2017.06.030 1270-9638/© 2017 Elsevier Masson SAS. All rights reserved.
inequality constraints [8,9]. However, to obtain a feasible solution to the closed-loop NOC in a fast manner is challenging, even for an unconstrained case [10]. Alternatively, Nonlinear Model Predictive Control (NMPC), which is based on receding horizon strategy (RHS) and re-planning of optimal trajectory in real time by solving the NOC at each sampling instant [11,12], has been proved as an effective method [13–16]. The resulting NOC problem can be further reduced to a quadratic programming (QP) problem that is computationally affordable for computers on-board spacecraft. Autonomous rendezvous with a target in near-field generally employs laser imaging detection and ranging system and advance video guidance system for relative navigation in terms of the relative distance and the Line-Of-Sight (LOS) angles with respect to the target. The aforementioned approaches were based on the LocalVertical Local-Horizontal (LVLH) formulation [5–7], including the LVLH based NMPC approach [13–16]. As a result, the relative navigation information has to be transformed from the LOS frame to the LVLH frame. The extra transformation between the LOS and LVLH frames complicates the derivation of guidance control and adds extra computational efforts for on-board computers. To reduce the computational requirement for the on-board computers, LOS based autonomous rendezvous were developed to employ the navigation directly [17,18]. The control law was developed based on the phase plane analysis. No attempt has been made to the LOS based NMPC in autonomous spacecraft rendezvous, to the best knowledge of authors.
102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132
JID:AESCTE
AID:4086 /FLA
[m5G; v1.218; Prn:30/06/2017; 16:23] P.2 (1-8)
P. Li, Z.H. Zhu / Aerospace Science and Technology ••• (••••) •••–•••
2
1 2
67
Nomenclature
68
3 4 5 6 7 8 9 10 11
69
a A( X ) Ak B Bk e f Fρ , Fθ
12 13 14 15
mc Nc , N P, Q
16 17
RT
semi-major axis of the elliptical orbit state-dependent system matrix discretized system matrix at time step k control input matrix discretized control matrix at time step k eccentricity of the elliptical orbit true anomaly of the elliptical orbit control force components along the LOS range and azimuth angle directions mass of the chaser numbers of control horizon and predictive horizon weight matrices on the control correction and control error distance from the center of the Earth to the target
Ts U max U x, y Xd U k
δ ρ , δθ θ, θd ρ , ρd
ρS σρ , σθ ω
sampling time interval maximum control force control input vector coordinates of the relative distance vector under orbital frame desired state vector incremental control variable measurement errors of LOS range and azimuth angle azimuth angle and desired azimuth angle LOS range and desired LOS range safety distance from the target standard deviations of Gaussian distributions in terms of LOS range and azimuth angle first order time derivative of the true anomaly
18
where R T = a(1 − e 2 )/(1 + e cos f ) is the distance from the center of the Earth to the target, f denotes the true anomaly, e is the eccentricity, a is the semi-major axis of the target, ω = δ˙ = μa(1 − e2 )/ R 2T and ω˙ = −2μe sin f / R 3T are the first and second order time derivatives of the true anomaly respectively, μ is the gravitational constant of the Earth, mc is the mass of the chaser, and F ρ and F θ are the force components in the ρ and θ directions of the LOS frame. The detailed derivation of Eq. (1) is given in Appendix. Introduce the new state vector X = {x1 , x2 , x3 , x4 } T with x1 = ρ , x2 = θ , x3 = ρ˙ and x4 = ρ θ˙ . Then, Eq. (1) is reduced to a set of first-order differential equations, such that,
20 21 22 23 24 25 26 27 28 29 30 31 32
⎡
34 35
Fig. 1. Schematic of the LOS frame.
37 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66
72 73 74 75 76 77 78 79 80 81 82 83
In the current work, the autonomous rendezvous is formulated in the LOS frame to directly use the relative navigation information. This intuitive description of the trajectory planning will result in concise formulation of the system dynamics and the controller design, especially for autonomous rendezvous with a passive noncooperative target. 2. Line-of-sight formulation of spacecraft dynamics Consider a chaser approaching a target in the orbital plane of an elliptical orbit, as shown in Fig. 1. To focus on the fundamentals of the LOS based NMPC, the current work is limited to the in-plane rendezvous with the target. The orbital frame, shown as O o − xo y o zo in Fig. 1, is defined with its origin at the Center of Mass (CM) of the target, where the y o -axis is along the orbital radius of the target, the xo -axis lies in orbital plane and is perpendicular to the y o -axis, and the zo -axis is normal to the orbital plane to complete a right-hand system. The LOS frame is formed by the range ρ and azimuth angle θ with its origin at the CM of the target, where the azimuth angle θ is measured from the xo -axis in the orbital plane as shown in Fig. 1. Accordingly, the equations of motion of the chaser in the LOS frame is expressed as per [19],
ρ¨ − ρ θ˙ 2 + 2ωρ θ˙ +
1 − 3 sin2 θ
μ R 3T
− ω2 ρ =
Fρ mc
μ Fθ ρ θ¨ + 2ρ˙ θ˙ − 2ωρ˙ − ω˙ + 3 3 sin θ cos θ ρ = RT
mc
(1a) (1b)
0
0
1
⎢ 0 0 0 ⎢ ⎢ μ A( X ) = ⎢ 2 2 0 ⎢ ω − R 3 (1 − 3 sin x2 ) 0 T ⎣ μ ω˙ + 3 R 3 sin x2 cos x2 0 2(ω − T ⎤ ⎡ and
85 86 87 88 89 90 91 92 93 94 95 96 97 98
˙ (t ) = A ( X ) X (t ) + BU (t ) X
33
38
71
84
19
36
70
0 0 ⎢0 0⎥ ⎥ B =⎢ ⎣1 0⎦ 0 1
0
x4 x1 x4 ) x1
⎤
(2)
99 100
⎥ 1 ⎥ x1 ⎥ − 2ω ⎥ ⎥ ⎦
101 102 103 104
0
105 106 107
(3)
where A ( X ) ∈ 4×4 and B ∈ 4×2 are the state-dependent system matrix and control input matrix, respectively, and U = { F x1 , F x2 } T /mc is the control input. It should be noted that the state x1 = ρ will never approach to zero in practice because the target has finite dimensions. The control objective of the trajectory planning for an autonomous spacecraft rendezvous is to approach the target with a desired state, X d = {x1d , x2d , 0, 0} T , smoothly to prevent sudden accelerations or decelerations of the chaser.
108 109 110 111 112 113 114 115 116 117 118 119 120 121
3. Nonlinear model predictive control
122 123
The control objective inevitably leads to an optimal control problem subject to constraint of thrust magnitude. The NMPC, also known as the Nonlinear Receding Horizon Control (NRHC), has been widely used to solve the constrained optimal control problem, due to its advantage of online generation of a set of feedback control commands by iteratively solving an open-loop NOC problem at each sampling instant. The receding horizon process is repeated by shifting the time one-step forward each time [20]. Accordingly, the optimal control problem of the trajectory planning
124 125 126 127 128 129 130 131 132
JID:AESCTE AID:4086 /FLA
[m5G; v1.218; Prn:30/06/2017; 16:23] P.3 (1-8)
P. Li, Z.H. Zhu / Aerospace Science and Technology ••• (••••) •••–•••
1 2 3 4 5 6 7 8 9 10 11 12 13 14
and tracking of spacecraft rendezvous can be formulated as a series of continuous-time nonlinear optimal control (CNOC) problems by NMPC at each sampling instant. Denote the starting and ending times of the rendezvous as t 1 and t f . Divide the total time interval into n subintervals evenly, such that, [t 1 , t 2 , · · · , tn+1 ] where tn+1 = t f . Furthermore, assume predictive (T p ) and control (T c ) horizons are constant with T p ≥ T c . Then, the receding predictive horizon [tk , tk + T p ] of the kth CNOC problem starts from tk ∈ [t 1 , tn+1 ] with the following definition. Problem Ψ C N O C : At any time instant tk ∈ [t 1 , tn+1 ], find the control input rate U˙ (t ) ∈ 2 and states X (t ) ∈ 4 that minimize the following quadratic cost function over a given predictive horizon T p ,
15
tk + T p
16 17 18
min J tk =
Q (t ) X (t ) − X d dt
tk+ T c
20 21
+
22
67
⎧ U k,k = U k + U k−1 ⎪ ⎪ ⎪ U k+1,k = (U k+1 + U k ) + U k−1 ⎪ ⎪ ⎪ ⎪ .. ⎪ ⎪ ⎪ . ⎪ ⎨ N c −1 U = U k + i + U k −1 ⎪ k + N − 1 , k c ⎪ ⎪ ⎪ ⎪ i = 0 ⎪ ⎪ .. ⎪ ⎪ ⎪ ⎪ ⎩. U k+i ,k = U k+ N c −1,k (i = N c , . . . , N )
73
U˙ (t ) T P (t )U˙ (t )dt
(4)
69 70 71 72
75 76 77
(10)
78 79 80 81 82 83 84
(11)
where U k+i ,k is the predicted control input at the (k + i)th time instant with respect to the kth time instant, U k−1 is the known N −1 control input at the kth time instant, U˜ k,k = i =c 0 U k+i ,k |i and
U˜ k =
68
74
˜ U˜ k + F˜ U k−1 U˜ k,k = M
tk
23
The resulting DNOC problem is computationally heavy and is transformed into a QP problem whose solution can be obtained as follows. Define the recursive relationship of the control input at any time instant with respect to the kth time instant in the control horizon as
or in a compact notation
tk
19
24
X (t ) − X d
T
3
N c −1
85 86 87 88 89 90
28
and the box constraint
˜ is U k+i |i (|· is the Ket vector notation), M ˜ (i , j ) = I 4 , (i ≥ j ), and F˜ = a lower triangular matrix with M Nc i =1 I 4 |i where I 4 is the identity matrix. Correspondingly, the set of states X k+ j ,k , ( j = 1, . . . , N ), in the
29
U j (t ) ≤ U max ∈
predictive horizon with respect to the kth time instant can be evaluated recursively, such that,
95
⎧ X k+1,k = A k X k + B k U k,k ⎪ ⎪ ⎪ X k+2,k = ( A k+1 A k ) X k + ( A k+1 B k )U k,k + B k+1 U k+1,k ⎪ ⎪ ⎪ ⎨. .. N −1 ⎪ N −1 N −1 ⎪ ⎪ ⎪ ⎪ A k+ j X k + A k+i ) B k+ j U k+ j ,k ( ⎪ ⎩ X k+ N ,k =
97
subject to the system dynamics
25 26 27
30 31 32 33 34 35 36 37 38 39 40 41 42 43
˙ (t ) = A ( X ) X (t ) + BU (t ) X
46
(6)
where U j (t ) is the jth element of the control input U (t ), U max is the maximum control force available, and ( Q (t ) ∈ 4×4 , P (t ) ∈ 2×2 ) are the time-varying and positive definite symmetric weight matrices, respectively. The solution of the kth CNOC problem is used as the initial condition for the (k + 1)th CNOC problem. To solve the CNOC problem, it is discretized into the discrete nonlinear optimal control (DNOC) problem over the given receding predictive horizon T p as follows. Problem Ψ DNOC : At any time instant k ∈ [1, n + 1], find the sequences of incremental control impulse U k = U k − U k−1 and state vectors X k that minimize the quadratic cost function over the given predictive horizon,
44 45
(5)
min J k =
k +N
{ X i − X d }T Q k { X i − X d } +
i =k+1
47
k+ N c −1 i =k
(7)
48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66
U iT P k U i
subject to the system dynamics
X i +1 = A i X i + B i U i ,
( i = k + 1, . . . , k + N )
(8)
and the box constraint
|U j |k ≤ U
max
(9)
where the corresponding predictive and control horizons T p = N T s and T c = N c T s in Eq. (4) are placed by N and N c with the assumption of N c ≤ N. Note that T s is the sampling time interval. Moreover, Eq. (8) is the discrete representation of Eq. (2) by zero-order holder, where the coefficient matrices A i and B i are assumed constant within each sampling time interval T s , but are updated at the beginning of each time interval with the states at the end of previous time interval. The value of weight matrices Q k and P k will be determined based on mission requirements or rendezvous strategy.
i =0
j =0
j =0
(12)
min J k =
− X˜ d
Q˜
˜ kp X
− X˜ d
T + U˜ k
(14)
⊕iN=c 1 P i
T
T
and P i , (i = 1, 2, ..., N c ), and thus Q˜ = Q˜ and P˜ = P˜ . ˜ d − φ X˜ kp − Defining and substituting an auxiliary vector E = X Γ U k−1 into Eq. (14) yield, [21]
=
1 2
T U˜ k
T
H U˜ k +
+ E T Q˜ E
105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 121 122 123 124 125 126 127 128 129
T Q˜ G U˜ k − E + U˜ k P˜ U˜ k T C T U˜ k
101
120
P˜ U˜ k
˜ d = i =1 X d |i , Q˜ = ⊕ N Q i and P˜ = where X (⊕ is the i =1 direct sum notation). It is worth noting that the weight matrices Q˜ and P˜ inherit the symmetric property from Q i , (i = 1, 2, ..., N)
100
104
N
min J k = G U˜ k − E
96
103
i =1 + j
vector is composed of the predicted states and the detailed expression of the coefficient matrices φ , Γ and G are given in Ref. [21]. Notably, the coefficient matrices φ , Γ and G are constant within each time step but are updated at the beginning of each step with the states at the end of previous step. For the sake of derivation convenience, the subscript k is omitted and the same is applied to the matrices E , H and C to be discussed later. Furthermore, the cost function in Eq. (7) is rewritten in a compact form,
94
102
p
T
93
99
˜ k = φ X k + Γ U k−1 + G U˜ k X (13) N p ˜ k = i =1 X k+i ,k |i and the superscript p indicates the where X
˜ kp X
92
98
Substituting Eq. (11) into Eq. (12) yields the set of predicted states X k+ j ,k , ( j = 1, . . . , N), in terms of U˜ k,k and in the compact notation as,
91
130 131
(15)
132
JID:AESCTE
AID:4086 /FLA
[m5G; v1.218; Prn:30/06/2017; 16:23] P.4 (1-8)
P. Li, Z.H. Zhu / Aerospace Science and Technology ••• (••••) •••–•••
4
1
67
2
68
3
69
4
70
5
71
6
72
7
73
8
74
9
75
10
76
11
77
12
78
13
79
14
80
15
81
16
82
17
83
18
84
Fig. 2. Block diagram of LOS based NMPC.
19 20 21 22 23 24 25 26 27 28
where H = 2[ G T Q˜ G + P˜ ] is the Hessian matrix that is positive definite and symmetric, C = −2G T Q˜ E is a column vector with the same dimension of U˜ k . Note that the vector E is constant within each time interval. Thus, the term E T Q˜ E in Eq. (15) is constant and does not affect the minimization of the cost function. It can be safely ignored from the cost function. Accordingly, Eq. (15) is equivalent to the following standard QP problem. Problem Ψ Q P , find U˜ k that minimizes the following cost function,
29 30 31 32 33 34 35
min J k =
1 2
T
U˜ k H U˜ k + C T U˜ k
(16)
The positive definite matrix H makes it a convex optimization problem where U˜ k is a global minimum. Moreover, the control magnitude constraint is expressed in terms of control increment U˜ k [21],
36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66
max ˜ U˜ k + F˜ U k−1 ≤ U˜ max −U˜ ≤M max
N
(17)
where U˜ = i =c 1 U max |i . Thus, Eq. (13) and Eq. (17) are the system dynamics and the control magnitude constraint, respectively. Compared our approach in Eqs. (1)–(17) with the approach based on the LVLH frame, the only differences are (i) the transformation of relative navigation information from LOS to LVLH frame is eliminated, (ii) the formation of state dependent system matrix. It is obvious the first difference results in saving of computational effort within each time step. For the second difference, the state dependent system matrix in the current approach is sparser than the matrix used in the common method, see the Ref. [22] in the revised manuscript. Thus, the associated computational effort will be lighter compared with that of Ref. [22]. Finally, measurement errors (δ ρ , δθ ) from navigation sensing systems are inevitable, such as the optical artifacts (e.g., glints, glares, saturations, hot pixels) in the sensors’ fields of view, which may result in sudden changes in control input. Theoretically, the sudden changes in control input could be indirectly rejected by smoothing the control input increment, which is the objective of the proposed optimal control cost function. Thus, the effectiveness of disturbance rejection will be demonstrated by adding the errors δ X = [δ ρ , δθ, 0, 0] T to the control input U k and state X k+1 at the beginning of each time instant. The navigation errors are assumed as Gaussian white noise with zero mean. The problem is then solved by the interior-point method due to its high convergence rate and the ease of implementation and the computational diagram of the close-loop control strategy is given as shown in Fig. 2.
85 86
Table 1 Parameters of weight matrices.
87 88
Parameter
Value
89
State weight matrix Control weight matrix
Q = diag(50, 5000, 1000, 1/ρ ) P = diag(100, 100)
90 91 92
4. Rendezvous strategy
93
The rendezvous of the chaser with the target is achieved by a judicious two-phased strategy to ensure the final azimuth direction at rendezvous to be satisfied precisely. In the first phase, the control task is to maneuver the chaser from its initial position to an intermediate point in the vicinity of the target for collision avoidance, such that, ρ > ρ S and ρ S is the safety distance from the target determined by the mission requirement. At the same, the chaser must adjust its azimuth angle from the initial condition θ0 to the desired value, θ = θd . In the second phase, the control task is to approach the target in a straight trajectory along the desired azimuth angle until ρ = ρd to accomplish the rendezvous mission. In the current work, it is assumed that ρ S = 20 m, θ0 and θd are determined by the particular rendezvous scenarios as specified in the following numerical examples. The above two-phased rendezvous strategy is achieved by the properly selected weight matrices, see in Table 1. Assume the chaser’s ability to maneuver in LOS range and angle directions is equally important in the rendezvous mission. This leads to the same penalties on the two control input (F ρ , F θ ) in the weight matrix P . The determination of weight matrix Q is heavily dependent on the rendezvous strategy. The control task in the first phase requires the chaser to aim the target in the required azimuth angle with a safe distance away. In this phase, the control accuracy for the azimuth angle is less important when the chaser is far away from the target. Accordingly, the penalty on the azimuth angle θ is set much larger than that on the LOS range ρ in this phase. In the second phase, it is critical to maintain a smooth trajectory in the proximity operation. Thus, the penalty for the third state (the approaching speed) is designed to slow down the relative velocity of the chaser for a safe approach. The last diagonal element in the weight matrix Q is the penalty for the fourth state – azimuth angular rate and is designed as state dependent, 1/ρ , instead of constant. Different from other three penalties, this penalty weight increases dramatically as the chaser approaches the target. Thus, the azimuth angular rate is reduced quickly to avoid the chattering of approaching trajectory because a stable azimuth angle at the final approach is mission critical. Again, it should be noted that the LOS range ρ does not decrease to zero in real situation even when
95
94 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132
JID:AESCTE AID:4086 /FLA
[m5G; v1.218; Prn:30/06/2017; 16:23] P.5 (1-8)
P. Li, Z.H. Zhu / Aerospace Science and Technology ••• (••••) •••–•••
1 2 3
Parameters
5 7 8 9
67
Table 2 Rendezvous conditions.
4 6
5
LOS range (ρ ), m Azimuth angle (θ ), ◦ Relative velocity (ρ˙ ), m/s Azimuth angle rate (θ˙ ), ◦ /s
68
Case 1
69
Case 2
Initial state
Desired state
Initial state
Desired state
70
80 0 0 0
1 0, 90, 180, 270 0 0
80 0, −90, −180, −270 0 0
1
72
−90
73
0 0
74
10 11 12 13 14 15
77
Table 3 Rendezvous navigation sensor uncertainties [23]. Measurands
Standard deviation (σ )
LOS range (ρ ), cm Azimuth angle (θ ), ◦
1.518 0.002787
78 79 80 81 82
17 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46
83
the rendezvous is completed. This is because the LOS frame is centered at the CM of the target and there is a minimum value of ρ measured from the CM to external surface of the target.
49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66
84 85 86 87
5. Case studies Fig. 3. In-plane approaching trajectories for Case 1.
The advantages of the proposed LOS NMPC in the spacecraft rendezvous are demonstrated by dynamic simulation of two rendezvous scenarios. The rendezvous conditions are listed in Table 2. The chaser is assumed to navigate by optical sensing system with direct measurement of the range and azimuth angle to the target. The sensor uncertainties, listed in Table 3 as per [23], are taken as a disturbance to the state, where σρ and σθ are the standard deviations of the Gaussian distributions in terms of LOS range ρ and the azimuth angle θ , respectively. The target is assumed moving in an elliptic orbit with a semimajor axis of 10,000 km and a true anomaly of 0◦ . The initial position of the target is assumed at perigee with an initial relative distance between the chaser and the target spacecraft of 80 m. To explore the application limit of the proposed approach, the eccentricity of 0.3 is used in the case studies. The predictive horizon and control horizons are defined as N = 10 and N c = 1. The mass of the chaser is assumed as 500 kg and the maximum available control thrust is limited to U max = 1 N. The simulation duration is set to 600 s and the sampling time is T s = 0.1 s. Note that the simulation results are shown in the LOS frame and the relative states are transformed to the orbital frame to show the approaching trajectory. 5.1. Approach target in different directions
47 48
75 76
16 18
71
In this case, the coplanar rendezvous with the target from various directions is studied. At t = 0, the chaser is assumed to running in the same orbit as the target, but 80 m ahead with 0◦ azimuth angle. The control objective is to approach the target along x0 -axis and y 0 -axis directions by controlling the azimuth angle θ equal to 0◦ , 90◦ , 180◦ and 270◦ , respectively. For the purpose of collision avoidance and subsequent proximity operations, the chaser is required to stop at 1 m away from the target in all directions at the end of approach. The detailed rendezvous conditions are listed as the Case 1 in Table 2. The simulation results of Case 1 are shown in Figs. 3–5. Fig. 3 shows the coplanar approaching trajectories in the orbital frame along x0 -axis and y 0 -axis directions, respectively. The chaser maneuvers around the target first to align the azimuth angle to the desired ones and then approaches to the desired position nearly in a straight line in the final stage as required, leading to a successful rendezvous. Next, Fig. 4 displays the time histories of the LOS range, relative velocity, azimuth angle and azimuth angle rate, respectively. Since
the penalty term 1/ρ is relatively small at the beginning, the control on the azimuth angle rate is weak and the angle approaches to the desired value at a high rate initially. As the LOS range reduces and the penalty term 1/ρ increases, the azimuth angle rate is reduced and the change of azimuth is flattened, see Fig. 4(c). As a result, the chaser adjusts its trajectory with respect to the target first and then approaches in a straight line. This approaching mode is widely adopted in space rendezvous missions to avoid collision with the target [6,?]. The azimuth angle reaches its desired value at 100 s while the chaser reaches its final position at 600 s. The proposed control scheme yields smooth trajectories in all cases, which is evident in Fig. 4. It shows there is nearly no overshoot in the time histories of the LOS range and the azimuth angle. To further demonstrate the effectiveness of the proposed control law, the control force profiles in the LOS range and azimuth angle channels are plotted in Fig. 5. It is interesting to note that the control force F ρ in the LOS range channel reaches its maximum output magnitude, −1 N, in the first tens of seconds, see Fig. 5(a), in order to quickly decrease the LOS range to reduce the rendezvous time. The negative symbol means the corresponding thrust is against the motion. The control force F θ in the azimuth angle channel converges to zero at around 100 s after the azimuth angle θ reaches its desired value, see Fig. 5(b), reflecting the effectiveness of the proposed two-phased rendezvous strategy. Furthermore, the chattering in the control force profiles is very small, which indicates that the proposed control law is effective in rejecting the measurement uncertainties to reach the desired values of the LOS range and azimuth angle smoothly, see Figs. 3–4. 5.2. Approach target from different initial azimuth angles
88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120
In this case, the coplanar rendezvous with the target is studied by assuming the chaser starting from different initial azimuth angles but with the same distance, 80 m, from the target. The control objective is to drive the chaser approaching the target in the − y 0 -axis direction from 4 different initial azimuth angles θ0 : 0◦ , −90◦ , −180◦ and −270◦ , respectively. Similarly to Case 1, the chaser is required to stop at 1 m away from the target in the − y 0 -axis direction at the end of rendezvous for subsequent proximity operations. The detailed rendezvous conditions are listed as the Case 2 in Table 2. The simulation results are shown in Figs. 6–8. Fig. 6 plots the rendezvous trajectories in the orbital frame. In all cases, the chaser
121 122 123 124 125 126 127 128 129 130 131 132
JID:AESCTE
6
AID:4086 /FLA
[m5G; v1.218; Prn:30/06/2017; 16:23] P.6 (1-8)
P. Li, Z.H. Zhu / Aerospace Science and Technology ••• (••••) •••–•••
1
67
2
68
3
69
4
70
5
71
6
72
7
73
8
74
9
75
10
76
11
77
12
78
13
79
14
80
15
81
16
82
17
83
18
84
19
85
20
86
21
87
22
88
23
89
Fig. 5. Control input: (a) LOS range and (b) azimuth angle.
24
90
25
91
26
92
27
93
28
94
29
95
30
96
31
97
32
98
33
99
34
100
35
101
36
102
Fig. 6. In-plane approaching trajectories for Case 2.
37 38
104
39 40 41 42 43 44 45 46 47
Fig. 4. Case 1: time histories of (a) LOS range (b) relative velocity (c) azimuth angle and (d) azimuth angle rate.
48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66
103
successfully maneuvers to the vicinity of the target in the − y 0 -axis direction first and then decreases the LOS range to the designated position in a straight line as required. Fig. 7 shows the responses of the LOS range, relative velocity, azimuth angle and azimuth angle rate, respectively. Similarly, the azimuth angle rate changes quickly in the first 70 s since the chaser is far away from the target and the penalty term 1/ρ is relatively small. As the chaser getting closer to the target, the penalty 1/ρ increases hyperbolically. Accordingly, the azimuth angle rate approaches zero quickly which implies the azimuth angle reaches its desired value. The chaser adjusts its orientation with respect to the target prior to approach a straight trajectory in the final approach as per the rendezvous strategy. The trajectories are smooth as required. Finally, Fig. 8 plots the control forces in the LOS range and azimuth angle channels, which show the same trends as those in Fig. 5. The simulations were conducted based on Eq. (2) and Eq. (A.1) under the simulation scenarios 1 and 2, respectively. The computational times were 134.26 s and 136.11 s under the environment
of i7-4510 2.6 GHz central processing unit and 16.0 GB RAM. The computational time of LOS model is slightly less as expected. It should be noted that the direct outputs of relative navigation system are the LOS range and the azimuth angle. Our comparison did not include the time delay due to the transformation of navigation information between the LOS and orbital frames, which depends on individual converter and circuitry. Therefore, the current approach is advantageous in saving computational time if one considers the possibility of eliminating navigation information transformation process under LOS frame. It is also worthy pointing out that the control consumption and state evolution are almost identical in the above cases for difference orbital eccentricities ranging from 0 to 0.3 although only the results of e = 0.3 are reported. This is because the trajectory of the target is almost straight during the short rendezvous period no matter what the eccentricity is, for the given initial condition of rendezvous, the rendezvous location and elliptical orbital parameters ( f = 0◦ and a = 10,000 km) in the current study. This is an interesting finding that indicates the eccentricity of target orbit has little impact on the rendezvous due to the short rendezvous period and the short distance between the target and the chaser compared with the distance from the target to the center of Earth.
105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127
6. Conclusions
128 129
This paper develops a LOS based NMPC for the coplanar autonomous rendezvous with non-cooperative targets. The LOS formulation simplifies the equations of relative dynamics of space-
130 131 132
JID:AESCTE AID:4086 /FLA
[m5G; v1.218; Prn:30/06/2017; 16:23] P.7 (1-8)
P. Li, Z.H. Zhu / Aerospace Science and Technology ••• (••••) •••–•••
7
1
67
2
68
3
69
4
70
5
71
6
72
7
73
8
74
9
75
10
76
11
77
12
78
13
79
14
80
15
81
16
82
17
83
18
84
19
85
20
86
21
87
22
88
23
89
24
90
25
91
Fig. 8. Control input: (a) LOS range and (b) azimuth angle.
26
92
27
93
Acknowledgements
28 29
94
This work is supported by the Discovery Grant and Discovery Accelerator Supplement Grant of the Natural Sciences and Engineering Research Council of Canada (NSERC).
30 31 32 33
37 38 39 40 41
42 43 44
x¨ y¨
=
45 47
Fig. 7. Case 2: time histories of (a) LOS range (b) relative velocity (c) azimuth angle and (d) azimuth angle rate.
50 51 53 54 55 56 57 58 59 60 61
craft using the relative navigation information from the measurement directly. The NMPC model is derived into a state dependent representation and the control is optimized for the control error and control smoothness simultaneously to achieve a smooth transient response. Furthermore, the proposed methodology is robust for various eccentricities of the target’s orbit because the trajectory of the target is almost straight during the short rendezvous period. Numerical simulations demonstrate the newly proposed method is effective and capable of achieving the control objective and it is insensitive to navigation measurement errors.
62
Conflict of interest statement
65 66
None declared.
1 mc
−2ω
0 Fx Fy
x˙ y˙
+
ω2 −ω˙ ω˙ ω2
x y
+
μ R 3T
−x
2y
102 103 104 105 106 107 108 109 110 111
(A.1)
112 113
where {x, y } T , {˙x, y˙ } T are the relative position and velocity between the chaser and the target spacecraft, ω is first order time derivative of the true anomaly, R T is the distance from the center of the Earth to the target, and { F x , F y } T is the applied force. The relationship between the orbital frame and the LOS frame can be expressed as,
114
120
x y
=
cos θ 0
0 sin θ
ρ
(A.2)
115 116 117 118 119 121 122 123
Accordingly, the first and second order derivatives of the relative states {x, y } T can be expressed in terms of LOS frame, such that,
63 64
0 2ω
+
46
52
98
101
The derivation of LOS model is based on Ref. [19]. Assume the distance between two spacecraft is sufficiently small compared to the orbit radius but sufficiently large compared to the largest dimension of either the target or the chaser. Then, the relative motion of spacecraft can be described by the linearized Tschauner– Hempel equation in the orbital frame, such that,
36
49
97
100
35
48
96
99
Appendix A. Derivation of line-of-sight model
34
95
x˙ = ρ˙ cos θ − ρ θ˙ sin θ y˙ = ρ˙ sin θ + ρ θ˙ cos θ
x¨ = ρ¨ − ρ θ˙ 2 cos θ − ρ θ¨ + 2ρ˙ θ˙ sin θ y¨ = ρ¨ − ρ θ˙ 2 sin θ + ρ θ¨ + 2ρ˙ θ˙ cos θ
Substituting Eq. (A.3) and Eq. (A.4) into Eq. (A.1) yields,
124 125 126
(A.3)
127 128 129
(A.4)
130 131 132
JID:AESCTE
AID:4086 /FLA
[m5G; v1.218; Prn:30/06/2017; 16:23] P.8 (1-8)
P. Li, Z.H. Zhu / Aerospace Science and Technology ••• (••••) •••–•••
8
1 2 3 4 5 6 7 8 9
⎧ Fρ μ μ ⎪ ⎪ ρ¨ − ρ (θ˙ − ω)2 + 3 ρ − 3 1 − 3 sin2 θ sin θ = ⎪ ⎨ mc RC RT 3 ⎪ R μ F θ ⎪ ⎪ ρ θ¨ + 2ρ˙ (θ˙ − ω) − ωρ ˙ − 2 1 − 3T cos θ = ⎩ RT
(A.5)
mc
RC
where { F ρ , F θ } T is the applied force expressed in the LOS frame and can be transformed from the orbital frame as,
10 11
Fρ Fθ
= A −1
Fx Fy
,
A −1 =
and
cos θ − sin θ
sin θ cos θ
12
Equation (A.5) can be further linearized by approximating the term
13
R 3T
14
R 3C
with first order Taylor expansion, such as,
15 16 17 18
R 3T R 3C
19
22 23 24 25 26 27
≈ 1+ ≈1−
20 21
30 31 32 33 34 35 36 37 38 39 40 41 42 43
RT
− 32
≈1−
3 2
·
2y RT
+
15 8
·
2y RT
2 + ···
3y
(A.6)
RT
and,
⎧ μ Fρ ⎪ 2 ⎪ ˙ 2 + 2ωρ θ˙ + 1 − 3 sin2 θ ¨ ρ − ρ − ω ρ= θ ⎪ ⎨ 3 mc RT ⎪ μ F θ ⎪ ⎪ ρ θ¨ + 2ρ˙ θ˙ − 2ωρ˙ − ω˙ + 3 sin θ cos θ ρ = ⎩ 3 RT
28 29
2y
(A.7)
mc
References [1] J.C. Liou, N.L. Johnson, A sensitivity study of the effectiveness of active debris removal in LEO, Acta Astronaut. 64 (2009) 236–243, http://dx.doi.org/10.1016/ j.actaastro.2008.07.009. [2] J.C. Liou, N.L. Johnson, N.M. Hill, Suboptimal power-limited rendezvous with fixed docking direction and collision avoidance, J. Guid. Control Dyn. 66 (2010) 648–653, http://dx.doi.org/10.1016/j.actaastro.2009.08.005. [3] W. Hao, Z.H. Zhu, D. Jin, H. Hu, Model predictive control with output feedback for a deorbiting electrodynamic tether system, J. Guid. Control Dyn. (2016), http://dx.doi.org/10.2514/1.G000535, in press. [4] N.-H. Chau, I. Sharf, Adaptive reactionless motion and parameter identification in postcapture of space debris, J. Guid. Control Dyn. 36 (2013) 404–414, http://dx.doi.org/10.2514/1.57856. [5] P. Lu, X. Liu, Autonomous trajectory planning for rendezvous and proximity operations by conic optimization, J. Guid. Control Dyn. 36 (2013) 375–389, http://dx.doi.org/10.2514/1.58436.
[6] G. Boyarko, O. Yakimenko, M. Romano, Optimal rendezvous trajectories of a controlled spacecraft and a tumbling object, J. Guid. Control Dyn. 34 (2011) 1239–1252, http://dx.doi.org/10.2514/1.47645. [7] F. Gavilan, R. Vazquez, E.F. Camacho, Chance-constrained model predictive control for spacecraft rendezvous with disturbance estimation, Control Eng. Pract. 20 (2012) 111–122, http://dx.doi.org/10.1016/j.conengprac.2011.09.006. [8] H.J. Peng, J. Xin, B.S. Chen, Optimal nonlinear feedback control of spacecraft rendezvous with finite low thrust between libration orbits, Nonlinear Dyn. 76 (2014) 1611–1632, http://dx.doi.org/10.1007/s11071-013-1233-9. [9] Z. Ma, O. Ma, N. Shashikanth, Optimal approach to and alignment with a rotating rigid body for capture, J. Astronaut. Sci. 55 (2007) 407–409, http://dx.doi. org/10.1007/BF03256532. [10] H.G. Zhang, H.L. Yan, D.R. Liu, Neural-network-based near-optimal control for a class of discrete-time affine nonlinear systems with control constraints, IEEE Trans. Neural Netw. 20 (2009) 1490–1503, http://dx.doi.org/10.1109/TNN.2009. 2027233. [11] E.F. Camacho, C. Bordons, Model Predictive Control, Springer-Verlag, London, 2004. [12] J.M. Maciejowski, Predictive Control with Constraints, Pearson Education, Harlow, UK, 2002. [13] E.N. Hartley, P.A. Trodden, A.G. Richards, J.M. Maciejowski, Model predictive control system design and implementation for spacecraft rendezvous, Control Eng. Pract. 20 (2012) 695–713, http://dx.doi.org/10.1016/j.conengprac.2012. 03.009. [14] L. Breger, J.P. How, Gauss’s variational equation-based dynamics and control for formation flying spacecraft, J. Guid. Control Dyn. 30 (2007) 437–448, http://dx.doi.org/10.2514/1.22649. [15] S.D. Cairano, H. Park, I. Kolmanovsky, Model predictive control approach for guidance of spacecraft rendezvous and proximity maneuvering, Int. J. Robust Nonlinear Control 22 (2012) 1398–1427, http://dx.doi.org/10.1002/rnc.2827.
67
[16] M. Leomanni, E. Rogers, S.B. Gabriel, Explicit model predictive control approach for low-thrust spacecraft proximity operations, J. Guid. Control Dyn. 37 (2014) 1780–1790, http://dx.doi.org/10.2514/1.G000477. [17] S.H. Yu, Range-rate control algorithms and space rendezvous schemes, J. Guid. Control Dyn. 20 (1997) 206–208, http://dx.doi.org/10.2514/2.4024. [18] S.H. Yu, Autonomous rendezvous in elliptical orbits, Acta Astronaut. 41 (1997) 95–101, http://dx.doi.org/10.1016/S0094-5765(97)00204-X. [19] T. Liu, J. Zhao, Spacecraft Dynamics, Harbin Institute of Technology Press, Harbin, China, 2003, pp. 86–87 (in Chinese). [20] A. Weiss, M. Baldwin, S.R. Erwin, I. Kolmanovsky, Model predictive control for spacecraft rendezvous and docking: strategies for handling constraints and case studies, IEEE Trans. Control Syst. Technol. 23 (2015) 1638–1646, http://dx.doi. org/10.1109/TCST.2014.2379639. [21] P. Li, Z.H. Zhu, S.A. Meguid, State dependent model predictive control for orbital rendezvous using pulse-width pulse-frequency modulated thrusters, Adv. Space Res. 58 (2016) 64–73, http://dx.doi.org/10.1016/j.asr.2016.04.022. [22] C.D. Karlgaard, Robust rendezvous navigation in elliptical orbit, J. Guid. Control Dyn. 29 (2006) 459–499, http://dx.doi.org/10.2514/1.19148. [23] L.C. Ma, X.Y. Meng, Z.Z. Liu, L.F. Du, Suboptimal power-limited rendezvous with fixed docking direction and collision avoidance, J. Guid. Control Dyn. 36 (2012) 229–239, http://dx.doi.org/10.2514/1.56449.
92
68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109
44
110
45
111
46
112
47
113
48
114
49
115
50
116
51
117
52
118
53
119
54
120
55
121
56
122
57
123
58
124
59
125
60
126
61
127
62
128
63
129
64
130
65
131
66
132