# Real-Time Optimal Approach and Capture of ENVISAT Based on Neural Networks

Authors: Hongjue Li; Yunfeng Dong; Peiyun Li
Journal: International Journal of Aerospace Engineering
 (2020)

---

## Abstract

A neural network-based controller is developed to enable a chaser spacecraft to approach and capture a disabled Environmental Satellite (ENVISAT). This task is conventionally tackled by framing it as an optimal control problem. However, the optimization of such a problem is computationally expensive and not suitable for onboard implementation. In this work, a learning-based approach is used to rapidly generate the control outputs of the controller based on a series of training samples. These training samples are generated by solving multiple optimal control problems with successive iterations. Then, Radial Basis Function (RBF) neural networks are designed to mimic this optimal control strategy from the generated data. Compared with a traditional controller, the neural network controller is able to generate real-time high-quality control policies by simply passing the input through the feedforward neural network.


---

## Body

## 1. Introduction

ENVISAT is one of the largest satellites launched by the European Space Agency (ESA) which lost communication with the ground station in April 2012 [1]. A collision between ENVISAT and an object of mass 10 kg could produce an extremely large cloud of debris, initiating a self-sustaining chain reaction of collisions and fragmentation with the production of new debris, a phenomenon known as the Kessler syndrome [2].The capture of ENVISAT is a complex space mission complicated by its uncertain rotation and the coupling motion between the chaser spacecraft and its manipulators. Many efforts were made by researchers to determine the attitude of ENVISAT. Sommer et al. [3] estimated the rotation vector and rotation speed of ENVISAT through ISAR images from 2011 to 2017, and he indicated that the rotational axis is not fixed and ENVISAT performs a nutation. Koshkin et al. [4] determined the spin motion of ENVISAT using photometry and provided an estimation formula to calculate the rotation period. Pittet et al. [5] suggested that ENVISAT performed no precession and a possible but limited nutation. Theoretically, these ground observations pave the way of estimating the rotation of ENVISAT by multiple equations. Unfortunately, in practice, there are still some inconsistencies within these observations that may suffer from different measurement information, measurement equipment, and target motion data. Such inaccuracy affects the performance of the capture mission.The capture mission is usually divided into two stages, i.e., the approach stage and the capture stage. In the approach stage, only 6-Degree-Of-Freedom (DOF) motion planning and control are considered since the manipulators are locked. In general, motion planning can be executed before the implementation of the tracking control. Boyarko et al. [6] employed Gauss Pseudospectral Method (GPM) to plan 6-DOF time and fuel optimal motion trajectories. Chu et al. [7] further considered an obstacle avoidance constraint during GPM planning. However, [6, 7] require accurate measurement of three-dimensional orbit position and relative distance information. In contrast, Zhang et al. [8] introduced an angles-only offline trajectory planning method which is tracked by a nonlinear model predictive control (NMPC). Offline planning has the advantage to optimize trajectories under various performance indices and constraints. However, onboard implementation is not efficient since these algorithms require accurate measurements and high computational expense. Lee and Vukovich [9] established a robust adaptive terminal sliding mode control scheme on the Lie group SE(3) that ensures the finite-time convergence of trajectory tracking errors. Sun and Zheng [10] extended the adaptive control to consider thrust misalignment in spacecraft rendezvous and proximity mission. Zhang et al. [11] further applied an adaptive finite-time fast terminal sliding mode controller to achieve spacecraft formation flying control. These approaches can obtain finite-time convergence control in the presence of model uncertainties and external disturbances.New approaches have been developed that can yield rapid or even real-time, near-optimal, trajectory generation by taking simplifications [12, 13] or approximations [14]. Other methods include planar optimal energy planning [15], inverse methods [16], modified Rapidly exploring Random Tree (RRT∗) [17], modified Fast Marching Tree (FMT∗) [18], convex planning [19], artificial potential field-based guidance [20], and polynomial planning [21]. Convex planning is suitable for problems subject to convex constraints such as descent and landing [19], unmanned path planning [22], and rendezvous and proximity operations [23, 24]. On the contrary, time-varying nonconvex motion constraint can be considered by combining artificial potential field with backstepping methodology [20]. Inverse dynamics based polynomial planning methods consisted of parameterizing nominal trajectories via high-order polynomials, converting the original optimal control problem into a nonlinear programming problem [21, 25]. These studies mainly focused on the rapid guidance of the chaser body, ending with final conditions suitable for docking or capture, while neglecting the dynamic coupling between the manipulator and the base in the capture stage. The orbit [26] and/or attitude [27] of the space robots are generally assumed to be uncontrolled so that the focus is placed on the motion planning and control of manipulators [26–28]; real-time near-optimal [29] or optimal [30] manipulator motion trajectories can be obtained.Significant advancements have been made with regard to both relative motion control and manipulator motion control; however, the reconciliation between the two fields remains to be solved. The division of the mission into approach stage and capture stage can reduce the complexity of the problem [15], but the sum of the two optimal results of the two stages is not necessarily the optimal solution of the whole mission. Considering the highly dynamic nature of the capture mission, the narrow capture window, and the limited energy, the objective is to find optimal trajectories and controls in real-time, during the whole approach and capture procedure. Benefiting from machine learning technology [31, 32], coordination between the optimal control and real-time control becomes possible. Izzo et al. [33] performed a survey on the artificial intelligence trends in spacecraft guidance dynamics and control, which involved a review of the intersection of the guidance and control with evolutionary optimization, tree searches, and machine learning. Sánchez-Sánchez and Izzo [34] applied deep learning methods to represent the optimal control action during a pinpoint landing. Furthermore, Izzo et al. [35] investigated the use of machine learning and evolutionary techniques in interplanetary trajectory design. Biggs and Fournier [36] trained a multilayer perceptron (MLP) to achieve real-time optimal attitude control of a 12U CubeSat. All these inspiring works attempted to train deep neural networks to learn from traditional optimal control solutions and later directly provide control commands by inputting the measured variables or even images.Offline solutions of optimal control problems require heavy computational costs to solve an iterative optimization that restricts their applications in real time. Meanwhile, online guidance and control only consider designing control policy in a single scenario, e.g., either 6-DOF guidance control or manipulator motion control. In this work, a more comprehensive implementation of the dynamic system that covers the spacecraft, the manipulators, and their inherent interactions is considered. While the system is more complicated, they can be efficiently solved via our proposed neural network-based approach. Three neural networks are trained to mimic the optimal control policy in a data-driven manner. Then, the optimal relative motion guidance and control for 6-DOF relative motion spacecraft can be implemented in real-time through neural network inference. The main contributions are stated as follows:(i)
Machine learning technology has been successfully applied to solving complicated optimal control problems. This problem involves considerably more complicated dynamics, DOFs of state variables and control variables compared with landing problems or interplanetary transfer problems
(ii)
Multiple neural networks are jointly designed to collaboratively learn optimal solutions in a synchronized framework. It is difficult to train a single network with a small number of samples, while the generation of huge optimal control solutions as training samples is costly. Simulation results indicate that multiple networks perform well through carefully designed inputs and outputs even when the number of samples are small
(iii)
The proposed control strategy exhibits comparable capture performance as tracking of optimal solutions in ideal cases and can be implemented in real timeThe remainder of this article is organized as follows. Section2 is devoted to the problem formulation, where the motion dynamic model of the approach and capture mission is presented. In Section 3, control strategies applied to generate optimal relative motion control and manipulator motion control solutions are presented. In Section 4, a series of Radial Basis Function (RBF) neural networks are designed and their inputs and outputs are discussed. As an illustration, numerical simulations are carried out in Section 5. Finally, brief concluding remarks in Section 6 end this study.


## 2. Problem Formulation

According to the measurement of Kucharski et al. [1], ENVISAT is presently flying in an approximate sun-synchronous orbit with a side rotating attitude. According to the ESA, the main body of ENVISAT can be regarded as a 5m×10m×5m cube rotating approximately around its z axis at an angular velocity of 0.8631°/s, as obtained by the approximate calculation of the formula proposed by Koshkin et al. [4]. The +z surface of this body is full of sensors and devices, the –z surface is obscured by antennas, the ±x surfaces are obscured by the Advanced Synthetic Aperture Radar (ASAR), and the –y surface is obscured by the solar array of the satellite. Consequently, only the +y surface is suitable for capture. A chaser spacecraft (hereinafter referred to as the chaser) installed with two 7-DOF manipulators is designed to approach and capture the ENVISAT. The chaser is a small spacecraft launched from a mother spacecraft, which reaches up to a dozen meters from the ENVISAT after a Lambert transfer. Subsequently, the spacecraft starts to approach and capture the ENVISAT, which is the focus of this paper.The overall mission is depicted in Figure1. The inertial coordinate system is termed Si (Oxiyizi). The centroid orbital coordinate system is termed So (Oxoyozo). The x axis of this system is along the flight velocity direction, z axis points to the earth, y axis can be determined by right-hand rule, and origin is at the centroid of the ENVISAT. The body coordinate system of the ENVISAT is termed Stb (Oxtbytbztb), and the body coordinate system of the chaser is termed Scb (Oxcbycbzcb). The x axes of these systems are along the spacecraft origin forward directions, and their z axes point to the earth when both the spacecraft are in the normal earth pointing flight mode.Figure 1
Schematic diagram of the approach and capture mission. The chaser spacecraft is equipped with two 7-DOF manipulators.


### 2.1. Dynamical Equations

A series of equations are required to describe the motion regulations of the ENVISAT, chaser, and manipulators, namely, the orbit kinematics equations, orbit dynamics equations [37], attitude kinematics equations, attitude dynamics equations [38], multibody kinetics equations, and multibody dynamics equations [39].The orbit and attitude motion regulation of the target can be expressed as(1)ṙtar=vtar,v̇tar=−μr3rtar,q̇tar=12Ωωtarqtar,Itarω̇tar+ωtar×Itarωtar=0,where the subscription “tar” refers to “target,” rtar=rtarx,rtary,rtarz is the 3-axis position vector of the ENVISAT, vtar=vtarx,vtary,vtarz is the 3-axis velocity vector of the ENVISAT, Itar is the moment of inertia of the target, qtar=qtar0,qtar1,qtar2,qtar3 is the quaternion of the ENVISAT, ωtar=ωtarx,ωtary,ωtarz is the 3-axis angular velocity of the ENVISAT, and Ωωtar can be expressed as
(2)Ωωtar=0ωtarz−ωtaryωtarx−ωtarz0ωtarxωtaryωtary−ωtarx0ωtarz−ωtarx−ωtary−ωtarz0.The motion regulation of the chaser with its manipulators can be expressed as(3)HΘΘ¨+CΘ̇,Θ=U,where Θ=rcha,echa,θ is the generalized coordinate of the chaser including its 3-axis position vector rcha=rchax,rchay,rchaz, 3-axis Euler angle echa=echax,echay,echaz (which is replaced by chaser’s quaternion qcha=qcha0,qcha1,qcha2,qcha3 to avoid ambiguity in Equation (5)), and the joint angles of the two 7-DOF manipulators θ=θ1,θ2,⋯,θ14; HΘ is a symmetrically positive definite inertia matrix of the chaser with manipulators; CΘ̇,Θ is a nonlinear matrix including Coriolis and centrifugal forces; and U=F,M,TT is the control variables for the chaser including its 3-axis control forces F=Fx,Fy,Fz, 3-axis control moment M=Mx,My,Mz, and the 14 joint torques T=T1,T2,⋯,T14. The subscription “cha” refers to “chaser.” Detailed definition of HΘ and CΘ̇,Θ can be seen in [40] and is not discussed here. The dynamic coupling of the chaser spacecraft and the manipulators can be more obviously seen if Equation (3) is written in matrix form as follows:
(4)MEMr~0gTJTωMr~0gHωHωϕJTωTHωϕHϕrchaechaθ+CΘ̇,Θ=FMT.The control variablesF, M, and T are related to all the components of the generalized coordinate. For example, F is not only related to the position (rcha), velocity (ṙcha in Θ̇), and acceleration (r¨cha in Θ¨) of the chaser spacecraft but also related to echa and θ. The coupled nature can be more clearly reflected by applying the Newton-Euler method [41] due to its recursive form; i.e., a forward recursion is performed to propagate the link velocities and accelerations, followed by a backward recursion to propagate forces.Equations (1) and (3) define a 54-state system of differential equations governing the approach and capture dynamics. Combined into the state vector X, these states are
(5)X=rcha,rtar,vcha,vtar,qcha,qtar,ωcha,ωtar,θ,θ̇T,where vcha=vchax,vchay,vchaz is the 3-axis velocity vector of the chaser, ωcha=ωchax,ωchay,ωchaz is the 3-axis angular velocity of the chaser, and θ̇=θ̇1,θ̇2,⋯,θ̇14 denotes the 14 angular velocities of the joints.The governing dynamics assume 20 controls that can be used by the chaser and the manipulators to achieve the approach and capture mission:(6)U=F,M,TT,which has been fully defined after Equation (3).


### 2.2. Performance Indices

Three performance indices are considered to evaluate the control optimality of the chaser, namely, the fuel consumptionJf, time consumption Jt, and impact of the manipulator end effector Jv.ConsideringIsp as the specific impulse of the thrusters, Jf can be expressed as
(7)Jf=∫t0tfFIspdt.Consideringtf as the final capture time, Jt can be expressed as
(8)Jt=tf.Jv can be measured using the x component of the velocity of the manipulator end effector relative to the desired capture point Δvend as follows:
(9)Jv=Δvendx,where Δvendx is the x axis projection of Δvend in Stb.By defining the fuel consumption coefficient askf, time consumption coefficient as kt, and impact coefficient as kv, the overall performance index J can be obtained as follows:
(10)J=kf⋅Jf+kt⋅Jt+kv⋅Jv.


### 2.3. Constraints

Two main constraints are considered during the approach and capture procedure, namely, the collision avoidance constraintC1 and state variable constraint C2. The collision avoidance can be taken into consideration by imposing the following path constraint: the mass center of the chaser spacecraft must remain at a distance larger than the minimum distance (a “keep out” sphere with a radius Δrmin) from the mass center of the ENVISAT, which can be defined as
(11)C1:rtar−rcha≥Δrmin.The state variable constraints state that all the state variables and control variables cannot exceed their range of values. By definingXmin and Xmax as the minimum and maximum values of state vector X, respectively, and Umin and Umax as the minimum and maximum values of controls U, respectively, the state variable constraints can be expressed as
(12)C2:Xmin≤X≤Xmax,Umin≤U≤Umax.


## 2.1. Dynamical Equations

A series of equations are required to describe the motion regulations of the ENVISAT, chaser, and manipulators, namely, the orbit kinematics equations, orbit dynamics equations [37], attitude kinematics equations, attitude dynamics equations [38], multibody kinetics equations, and multibody dynamics equations [39].The orbit and attitude motion regulation of the target can be expressed as(1)ṙtar=vtar,v̇tar=−μr3rtar,q̇tar=12Ωωtarqtar,Itarω̇tar+ωtar×Itarωtar=0,where the subscription “tar” refers to “target,” rtar=rtarx,rtary,rtarz is the 3-axis position vector of the ENVISAT, vtar=vtarx,vtary,vtarz is the 3-axis velocity vector of the ENVISAT, Itar is the moment of inertia of the target, qtar=qtar0,qtar1,qtar2,qtar3 is the quaternion of the ENVISAT, ωtar=ωtarx,ωtary,ωtarz is the 3-axis angular velocity of the ENVISAT, and Ωωtar can be expressed as
(2)Ωωtar=0ωtarz−ωtaryωtarx−ωtarz0ωtarxωtaryωtary−ωtarx0ωtarz−ωtarx−ωtary−ωtarz0.The motion regulation of the chaser with its manipulators can be expressed as(3)HΘΘ¨+CΘ̇,Θ=U,where Θ=rcha,echa,θ is the generalized coordinate of the chaser including its 3-axis position vector rcha=rchax,rchay,rchaz, 3-axis Euler angle echa=echax,echay,echaz (which is replaced by chaser’s quaternion qcha=qcha0,qcha1,qcha2,qcha3 to avoid ambiguity in Equation (5)), and the joint angles of the two 7-DOF manipulators θ=θ1,θ2,⋯,θ14; HΘ is a symmetrically positive definite inertia matrix of the chaser with manipulators; CΘ̇,Θ is a nonlinear matrix including Coriolis and centrifugal forces; and U=F,M,TT is the control variables for the chaser including its 3-axis control forces F=Fx,Fy,Fz, 3-axis control moment M=Mx,My,Mz, and the 14 joint torques T=T1,T2,⋯,T14. The subscription “cha” refers to “chaser.” Detailed definition of HΘ and CΘ̇,Θ can be seen in [40] and is not discussed here. The dynamic coupling of the chaser spacecraft and the manipulators can be more obviously seen if Equation (3) is written in matrix form as follows:
(4)MEMr~0gTJTωMr~0gHωHωϕJTωTHωϕHϕrchaechaθ+CΘ̇,Θ=FMT.The control variablesF, M, and T are related to all the components of the generalized coordinate. For example, F is not only related to the position (rcha), velocity (ṙcha in Θ̇), and acceleration (r¨cha in Θ¨) of the chaser spacecraft but also related to echa and θ. The coupled nature can be more clearly reflected by applying the Newton-Euler method [41] due to its recursive form; i.e., a forward recursion is performed to propagate the link velocities and accelerations, followed by a backward recursion to propagate forces.Equations (1) and (3) define a 54-state system of differential equations governing the approach and capture dynamics. Combined into the state vector X, these states are
(5)X=rcha,rtar,vcha,vtar,qcha,qtar,ωcha,ωtar,θ,θ̇T,where vcha=vchax,vchay,vchaz is the 3-axis velocity vector of the chaser, ωcha=ωchax,ωchay,ωchaz is the 3-axis angular velocity of the chaser, and θ̇=θ̇1,θ̇2,⋯,θ̇14 denotes the 14 angular velocities of the joints.The governing dynamics assume 20 controls that can be used by the chaser and the manipulators to achieve the approach and capture mission:(6)U=F,M,TT,which has been fully defined after Equation (3).


## 2.2. Performance Indices

Three performance indices are considered to evaluate the control optimality of the chaser, namely, the fuel consumptionJf, time consumption Jt, and impact of the manipulator end effector Jv.ConsideringIsp as the specific impulse of the thrusters, Jf can be expressed as
(7)Jf=∫t0tfFIspdt.Consideringtf as the final capture time, Jt can be expressed as
(8)Jt=tf.Jv can be measured using the x component of the velocity of the manipulator end effector relative to the desired capture point Δvend as follows:
(9)Jv=Δvendx,where Δvendx is the x axis projection of Δvend in Stb.By defining the fuel consumption coefficient askf, time consumption coefficient as kt, and impact coefficient as kv, the overall performance index J can be obtained as follows:
(10)J=kf⋅Jf+kt⋅Jt+kv⋅Jv.


## 2.3. Constraints

Two main constraints are considered during the approach and capture procedure, namely, the collision avoidance constraintC1 and state variable constraint C2. The collision avoidance can be taken into consideration by imposing the following path constraint: the mass center of the chaser spacecraft must remain at a distance larger than the minimum distance (a “keep out” sphere with a radius Δrmin) from the mass center of the ENVISAT, which can be defined as
(11)C1:rtar−rcha≥Δrmin.The state variable constraints state that all the state variables and control variables cannot exceed their range of values. By definingXmin and Xmax as the minimum and maximum values of state vector X, respectively, and Umin and Umax as the minimum and maximum values of controls U, respectively, the state variable constraints can be expressed as
(12)C2:Xmin≤X≤Xmax,Umin≤U≤Umax.


## 3. Real-Time Optimal Control Strategy

It is difficult to directly optimize the entire problem due to the complex coupling motion of the orbit, attitude, and manipulator. Thus, it is necessary to reduce the complexity of the problem step by step. The coupling motion of the orbit and attitude is first considered and optimized using the GPM. The manipulator deployment, point, and capture control are then considered, and a feedforward compensation control algorithm is designed to eliminate the coupling among the manipulators, orbit, and attitude. Finally, an iterative optimization is established to obtain the optimal control strategy of the manipulators, orbit, and attitude.LetW be the range of initial values of the state variables. A set of initial values X0 can be randomly selected in W. Let Φ be a neighborhood sufficiently close to X0. A series of value sets Xi can be randomly selected in Φ, and each value set is known as a neighborhood member. For each neighborhood member, the optimal control sequence Ut and state sequence Xt can be obtained by using the methods described in the following sections.


### 3.1. Relative Motion Control

When the manipulators are locked, only the orbit and attitude motions of the chaser are needed to be considered. The state variables and control variables can be written as(13)x=Δr,Δv,qcha,qtar,ωcha,ωtarT,u=F,MT,where Δr and Δv denote the relative position and velocity vectors between the chaser and the ENVISAT, respectively.The relative orbit dynamics of the chaser and ENVISAT in the orbit frame centered at ENVISAT can be described using the Hill–Clohessy–Wiltshire (C–W) equations [42]. By replacing the orbit dynamics equations with the C–W equations and removing the multibody kinetic and dynamic equations, a series of explicit differential equations to govern the law of motion of the approaching stage can be obtained as follows:
(14)ẋ=fx,u.Since the manipulators are not considered, the value ofkv in the overall performance index J is zero. Thus, the approach mission can be transferred into an optimal control problem. That is, for a system described by Equation (14), the control u is calculated to transfer the state vector x from the initial state xi to the terminal state xf under constraints C1 and C2 while optimizing the performance index J.In the optimal control problem, besides path constraintsC1 and C2, certain extra constraints at the final time tf (often referred to as the terminal constraints) are needed to define xf and make the optimization problem solvable. The terminal constraints include the position constraint T1, velocity constraint T2, and attitude constraint T3 as follows:
(15)T1:Δrd=Δrtf≥rmin,T2:Δvtf=Loiωchatf×Δrtf+V,T3:Lob⋅1,0,0T=ΔrtfΔrtf,where Loi is the rotation matrix between Si and So and V is the extra relative velocity. If V=0, the chaser remains relatively static with the ENVISAT; otherwise, the chaser flies over the ENVISAT at Δrd with a relative speed of V. Lob is the rotation matrix between Scb and So and is determined by qcha.The above optimal control problem can be solved using the well-developed numerical toolkit Gauss Pseudospectral OPtimization Software (GPOPS) [43]. The output is the optimal state variable sequence xt and the optimal control variable sequence ut, t=t0,⋯,tf. Boyarko et al. [6] proved that a pseudospectral solver provided a solution remarkably close to the truly optimal one by using the results obtained by applying the Minimum Principle (MP). In their research, the chaser started to approach at a distance of 5 m, and the MP solution returned a value of a minimum control performance index J=0.11341185 compared with J=0.1133obtained by the pseudospectral solver. The GPOPS has also been employed to compute the optimal trajectories and generate the appropriate training set in a lunar landing mission [44]. In this case, the training set was used to train a deep Convolutional Neural Network (CNN) followed by a Recurrent Neural Network connecting with a Long Short-Term Memory (RNN-LSTM) and is effective in the planar landing problem.Given the size of the ENVISAT, a similar approach scenario is considered, in which the chaser starts at a distance of approximately 15 m from the mass center of the ENVISAT. The ENVISAT has an initial small angular velocity around itsx body axis and an angular velocity of 0.1484 rad/s around its z body axis. The minimum safe distance Δrmin is set to be 6.8 m, which means that the chaser approaches the ENVISAT to a distance of 1.8 m since the distance of the mass center of the ENVISAT to its +y surface is 5 m. The detailed ranges of all the initial values are presented in Table 1 to form a neighborhood, from which 388 numbers of neighborhood members (initial conditions) are randomly selected.Table 1
Neighborhood.
Variable
Value range
Δrx (m)
[-0.5, 0.5]
Δry (m)
[14.0, 16.0]
Δrz (m)
[-0.5, 0.5]
ΔvxΔvyΔvz (m·s-1)
0.0
ωtarx (rad·s-1)
[0.0, 0.0349]
ωtary (rad·s-1)
0.0
ωtarz (rad·s-1)
0.1484
ωchaxωchayωchaz (rad·s-1)
0.0
etarxetaryetarz (rad)
0.0
echaxechayechaz (rad)
0.0Researchers determine the coefficients in the performance index according to the mission requirements. If the chaser is required to approach the ENVISAT as fast as possible,kf is set as zero to ensure that the optimal control problem becomes a minimum time control problem. On the contrary, if the chaser is required to approach the ENVISAT with the least fuel consumption, the problem becomes a minimum control problem, and kt is set as zero. Here, kf is set to be 2/3 and kt is set to be 1, which means that both the fuel consumption and approach time are considered.For each neighborhood member, an optimal control sequence and the corresponding optimal trajectory are determined using the GPM. Each trajectory consists of 40 discrete points. The 388 optimal relative orbit trajectories are shown in Figure2. A sample relative orbit trajectory (m=12) along with the attitude of the chaser is shown in Figure 3.Figure 2
Optimal trajectories of all 388 neighborhood members. The radius of the sphere on the left isΔrmin. The blue cuboid on the left is the ENVISAT with an attitude of t0. The red, blue, and green lines perpendicular to one another denote the +x, +y, and +z axes, respectively, of the orbital coordinate system So.Figure 3
Optimal trajectory of the 12th neighborhood member. The light green trajectory is the track of the “subsatellite point” of the chaser on the ENVISAT. The dark green cuboids denote the chaser body.Note that not all the optimal trajectories and corresponding control sequences are used to train the networks. 26 trajectories of bad quality are eliminated so that only 362 trajectories are used to train the networks.


### 3.2. Manipulator Control

The manipulator installed on the chaser spacecraft is shown in Figure4. The base coordinate system is termed Sm (more specifically termed Slm/Srm for the left and right manipulators). The modified Denavit–Hartenberg (D–H) method is used to establish the multibody kinetics of the manipulator, and the D–H parameters are listed in Table 2. The initial joint angles of the folded manipulators are listed in the last column of Table 2.Figure 4
Definition of joints and links of the 7-DOF manipulator.Table 2
D–H parameters of the 7-DOF manipulator.
Joint (link)i
αi (°)
ai (m)
di (m)
θi (°)
1
0
0
0
15.57335
2
-90
0
0
120.40453
3
90
0
l2+l3
28.84101
4
-90
0
0
99.4784
5
90
0
l4+l5
0
6
-90
0
0
46.70365
7
90
0
0
0A common deployment strategy is designed by using the following defined desired joint angular accelerations:(16)θ¨it=−ΔθiΔt/22,t<Δt2,ΔθiΔt/22,t≥Δt2,i=1,2,⋯,14,where Δt=ti−te is the total deployment time.By integrating Equation (16), the joint angular velocities can be obtained, and by integrating it again, the joint angles can be obtained. By determining θi,θ̇i,θ¨i,vc,v̇c,vt,v̇t,ωc,ω̇c, the joint torques T at any time during the deployment can be calculated by using Equation (3).After deployment, a Proportional-Differential (PD) control algorithm is designed to control the joints of the manipulators. The calculation of the desired joint angles can be divided into the following 4 steps:Step 1.
Letθ3=θ3const, and calculate the angles of the other joints as 1θi.
The desired angle of joints 4, 2, and 1 can be calculated as follows:(17)θ4=±π−cos−12l2−pdx2−pdy2−pdz22l2,θ2=sin−1pdzA2+B2−C,π−sin−1pdzA2+B2+C,θ1=tan−1pdxny−pdynx/pdxmy−ppymxmxny−mynx/mynx−mxny,where A, B, C, mx, nx, my, and ny are
(18)A=−l2sinθ4cosθ3,B=l1+l2cosθ4,C=−tan−1BA,mx=−lsinθ4sinθ3,nx=lsinθ4cosθ2cosθ3+l+lcosθ4sinθ2,my=nx,ny=lsinθ4sinθ3.
The rotation matrixM04 can be derived from the abovementioned desired angles of joints 1, 2, 3, and 4. M07 can be obtained according to the desired position and attitude of the manipulator end effector. Thus, M47 (transformation matrix between the coordinate systems of links 4 and 7) can be calculated by performing the inverse product of M04 and M07 as follows:
(19)M47=M04−1M07=g11g12g13dxg21g22g23dyg31g32g33dz0001.
Subsequently, the desired angles of joints 5, 6, and 7 can be calculated using the components ofM47 as follows:
(20)θ5=tan−1g33g13,θ6=cos−1−g23,θ7=tan−1g22g21,or
(21)θ5=tan−1−g33−g13,θ6=cos−1−g23,θ7=tan−1−g22−g21.
For the attitude of the end effector of a certain manipulator, 2 groups of values ofθ5, θ6, and θ7 are calculated using Equations (20) and (21). For the position of the end effector of a certain manipulator, 4 groups of values of θ1, θ2, θ3, and θ4 are calculated using Equation (17). Together, 8 groups of desired joint angles are available. One group of joint angles is chosen as the desired joint angles to minimize the sum of the rotation angles from the current angles to the desired angles.Step 2.
Letθ3=θ3const+Δθ3, and calculate the angles of the other joints again as 2θi using Equations (17), (20), and (21).Step 3.
Linearize the equation for the manipulator in intervals1θi,2θi. The variation of the angle of joint 3 Δθ3best can be optimized by using linear programming.Step 4.
Letθ3=θ3+Δθ3best, and obtain the other joint angles by solving Equations (17), (20), and (21).
The desired angular velocity of each joint can be calculated by using the Moore-Penrose inverse of the Jacobian matrixJθ as follows:
(22)θ̇=Jθ+⋅ve,ωeT,Jθ+=JθTJθJθT−1,where ve is the desired velocity of the manipulator end effector and ωe is the angular velocity of the manipulator end effector.
Since the desired angle and angular velocity of each joint have been calculated, the PD control law is employed to solve for the joint torque as follows:(23)Ti=kp⋅θi+kd⋅θ̇i.Twelve neighborhood members are randomly selected from the 388 neighborhood members. Each member of the manipulator neighborhood is calculated, and all the trajectories of the joint angle varying with time are obtained. Taking joint 2 of the left manipulator as an example, the optimal joint angle trajectories of all the 12 neighborhood members are determined, as shown in Figure5.Figure 5
Twelve optimal joint angle trajectories for joint 2 of the left manipulator.


### 3.3. Feedforward Compensation and Iterative Optimization

A recursive form of the multirigid body dynamic model can be established using the Newton-Euler method as follows:(24)Γi=hθi,θ̇i,θ¨i,v̇cha,v̇tar,i=1,2,.⋯,14,where Γi=fi,miT denotes the force and moment applied to each joint, mi denotes the transferring torque of each joint, and fi denotes the transferring force of each joint.Wheni=0, Γ0=f0,m0T denotes the transferring force and moment between the manipulator base and the chaser body. By calculating and applying these values inversely as the feedforward to the chaser orbit and attitude control, the control variables of the chaser can be updated as follows:
(25)Ut=ut+Γ0t.Since the disturbance is estimated, it cannot be precisely eliminated by using the feedforward compensation usingut. The possible influences may include the deviation of the optimal trajectory, which makes ut not optimal, as shown in Figure 6. The orbit and attitude motion trajectory is assumed as xt0,⋯,xti,⋯,xtf (solid line) by assigning the terminal state xtf when the manipulators do not move. When the manipulators move and finish the capturing, the orbit and attitude motion trajectory is modified to x∗t0,⋯,x∗ti,⋯,x∗tf (dash line) with the terminal state x∗tf. One or more iterations are performed to update ut to make x∗tf closer to xtf. Considering that Δx∗tf=x∗tf−xtf, the iteration stops when Δx∗tf is less than a preset constant value, and the optimal orbit and attitude trajectory x′tf (dot-dash line) is obtained. The corresponding optimal control sequence is Ut=u′t+Γ0t.Figure 6
Iterative optimization process.


## 3.1. Relative Motion Control

When the manipulators are locked, only the orbit and attitude motions of the chaser are needed to be considered. The state variables and control variables can be written as(13)x=Δr,Δv,qcha,qtar,ωcha,ωtarT,u=F,MT,where Δr and Δv denote the relative position and velocity vectors between the chaser and the ENVISAT, respectively.The relative orbit dynamics of the chaser and ENVISAT in the orbit frame centered at ENVISAT can be described using the Hill–Clohessy–Wiltshire (C–W) equations [42]. By replacing the orbit dynamics equations with the C–W equations and removing the multibody kinetic and dynamic equations, a series of explicit differential equations to govern the law of motion of the approaching stage can be obtained as follows:
(14)ẋ=fx,u.Since the manipulators are not considered, the value ofkv in the overall performance index J is zero. Thus, the approach mission can be transferred into an optimal control problem. That is, for a system described by Equation (14), the control u is calculated to transfer the state vector x from the initial state xi to the terminal state xf under constraints C1 and C2 while optimizing the performance index J.In the optimal control problem, besides path constraintsC1 and C2, certain extra constraints at the final time tf (often referred to as the terminal constraints) are needed to define xf and make the optimization problem solvable. The terminal constraints include the position constraint T1, velocity constraint T2, and attitude constraint T3 as follows:
(15)T1:Δrd=Δrtf≥rmin,T2:Δvtf=Loiωchatf×Δrtf+V,T3:Lob⋅1,0,0T=ΔrtfΔrtf,where Loi is the rotation matrix between Si and So and V is the extra relative velocity. If V=0, the chaser remains relatively static with the ENVISAT; otherwise, the chaser flies over the ENVISAT at Δrd with a relative speed of V. Lob is the rotation matrix between Scb and So and is determined by qcha.The above optimal control problem can be solved using the well-developed numerical toolkit Gauss Pseudospectral OPtimization Software (GPOPS) [43]. The output is the optimal state variable sequence xt and the optimal control variable sequence ut, t=t0,⋯,tf. Boyarko et al. [6] proved that a pseudospectral solver provided a solution remarkably close to the truly optimal one by using the results obtained by applying the Minimum Principle (MP). In their research, the chaser started to approach at a distance of 5 m, and the MP solution returned a value of a minimum control performance index J=0.11341185 compared with J=0.1133obtained by the pseudospectral solver. The GPOPS has also been employed to compute the optimal trajectories and generate the appropriate training set in a lunar landing mission [44]. In this case, the training set was used to train a deep Convolutional Neural Network (CNN) followed by a Recurrent Neural Network connecting with a Long Short-Term Memory (RNN-LSTM) and is effective in the planar landing problem.Given the size of the ENVISAT, a similar approach scenario is considered, in which the chaser starts at a distance of approximately 15 m from the mass center of the ENVISAT. The ENVISAT has an initial small angular velocity around itsx body axis and an angular velocity of 0.1484 rad/s around its z body axis. The minimum safe distance Δrmin is set to be 6.8 m, which means that the chaser approaches the ENVISAT to a distance of 1.8 m since the distance of the mass center of the ENVISAT to its +y surface is 5 m. The detailed ranges of all the initial values are presented in Table 1 to form a neighborhood, from which 388 numbers of neighborhood members (initial conditions) are randomly selected.Table 1
Neighborhood.
Variable
Value range
Δrx (m)
[-0.5, 0.5]
Δry (m)
[14.0, 16.0]
Δrz (m)
[-0.5, 0.5]
ΔvxΔvyΔvz (m·s-1)
0.0
ωtarx (rad·s-1)
[0.0, 0.0349]
ωtary (rad·s-1)
0.0
ωtarz (rad·s-1)
0.1484
ωchaxωchayωchaz (rad·s-1)
0.0
etarxetaryetarz (rad)
0.0
echaxechayechaz (rad)
0.0Researchers determine the coefficients in the performance index according to the mission requirements. If the chaser is required to approach the ENVISAT as fast as possible,kf is set as zero to ensure that the optimal control problem becomes a minimum time control problem. On the contrary, if the chaser is required to approach the ENVISAT with the least fuel consumption, the problem becomes a minimum control problem, and kt is set as zero. Here, kf is set to be 2/3 and kt is set to be 1, which means that both the fuel consumption and approach time are considered.For each neighborhood member, an optimal control sequence and the corresponding optimal trajectory are determined using the GPM. Each trajectory consists of 40 discrete points. The 388 optimal relative orbit trajectories are shown in Figure2. A sample relative orbit trajectory (m=12) along with the attitude of the chaser is shown in Figure 3.Figure 2
Optimal trajectories of all 388 neighborhood members. The radius of the sphere on the left isΔrmin. The blue cuboid on the left is the ENVISAT with an attitude of t0. The red, blue, and green lines perpendicular to one another denote the +x, +y, and +z axes, respectively, of the orbital coordinate system So.Figure 3
Optimal trajectory of the 12th neighborhood member. The light green trajectory is the track of the “subsatellite point” of the chaser on the ENVISAT. The dark green cuboids denote the chaser body.Note that not all the optimal trajectories and corresponding control sequences are used to train the networks. 26 trajectories of bad quality are eliminated so that only 362 trajectories are used to train the networks.


## 3.2. Manipulator Control

The manipulator installed on the chaser spacecraft is shown in Figure4. The base coordinate system is termed Sm (more specifically termed Slm/Srm for the left and right manipulators). The modified Denavit–Hartenberg (D–H) method is used to establish the multibody kinetics of the manipulator, and the D–H parameters are listed in Table 2. The initial joint angles of the folded manipulators are listed in the last column of Table 2.Figure 4
Definition of joints and links of the 7-DOF manipulator.Table 2
D–H parameters of the 7-DOF manipulator.
Joint (link)i
αi (°)
ai (m)
di (m)
θi (°)
1
0
0
0
15.57335
2
-90
0
0
120.40453
3
90
0
l2+l3
28.84101
4
-90
0
0
99.4784
5
90
0
l4+l5
0
6
-90
0
0
46.70365
7
90
0
0
0A common deployment strategy is designed by using the following defined desired joint angular accelerations:(16)θ¨it=−ΔθiΔt/22,t<Δt2,ΔθiΔt/22,t≥Δt2,i=1,2,⋯,14,where Δt=ti−te is the total deployment time.By integrating Equation (16), the joint angular velocities can be obtained, and by integrating it again, the joint angles can be obtained. By determining θi,θ̇i,θ¨i,vc,v̇c,vt,v̇t,ωc,ω̇c, the joint torques T at any time during the deployment can be calculated by using Equation (3).After deployment, a Proportional-Differential (PD) control algorithm is designed to control the joints of the manipulators. The calculation of the desired joint angles can be divided into the following 4 steps:Step 1.
Letθ3=θ3const, and calculate the angles of the other joints as 1θi.
The desired angle of joints 4, 2, and 1 can be calculated as follows:(17)θ4=±π−cos−12l2−pdx2−pdy2−pdz22l2,θ2=sin−1pdzA2+B2−C,π−sin−1pdzA2+B2+C,θ1=tan−1pdxny−pdynx/pdxmy−ppymxmxny−mynx/mynx−mxny,where A, B, C, mx, nx, my, and ny are
(18)A=−l2sinθ4cosθ3,B=l1+l2cosθ4,C=−tan−1BA,mx=−lsinθ4sinθ3,nx=lsinθ4cosθ2cosθ3+l+lcosθ4sinθ2,my=nx,ny=lsinθ4sinθ3.
The rotation matrixM04 can be derived from the abovementioned desired angles of joints 1, 2, 3, and 4. M07 can be obtained according to the desired position and attitude of the manipulator end effector. Thus, M47 (transformation matrix between the coordinate systems of links 4 and 7) can be calculated by performing the inverse product of M04 and M07 as follows:
(19)M47=M04−1M07=g11g12g13dxg21g22g23dyg31g32g33dz0001.
Subsequently, the desired angles of joints 5, 6, and 7 can be calculated using the components ofM47 as follows:
(20)θ5=tan−1g33g13,θ6=cos−1−g23,θ7=tan−1g22g21,or
(21)θ5=tan−1−g33−g13,θ6=cos−1−g23,θ7=tan−1−g22−g21.
For the attitude of the end effector of a certain manipulator, 2 groups of values ofθ5, θ6, and θ7 are calculated using Equations (20) and (21). For the position of the end effector of a certain manipulator, 4 groups of values of θ1, θ2, θ3, and θ4 are calculated using Equation (17). Together, 8 groups of desired joint angles are available. One group of joint angles is chosen as the desired joint angles to minimize the sum of the rotation angles from the current angles to the desired angles.Step 2.
Letθ3=θ3const+Δθ3, and calculate the angles of the other joints again as 2θi using Equations (17), (20), and (21).Step 3.
Linearize the equation for the manipulator in intervals1θi,2θi. The variation of the angle of joint 3 Δθ3best can be optimized by using linear programming.Step 4.
Letθ3=θ3+Δθ3best, and obtain the other joint angles by solving Equations (17), (20), and (21).
The desired angular velocity of each joint can be calculated by using the Moore-Penrose inverse of the Jacobian matrixJθ as follows:
(22)θ̇=Jθ+⋅ve,ωeT,Jθ+=JθTJθJθT−1,where ve is the desired velocity of the manipulator end effector and ωe is the angular velocity of the manipulator end effector.
Since the desired angle and angular velocity of each joint have been calculated, the PD control law is employed to solve for the joint torque as follows:(23)Ti=kp⋅θi+kd⋅θ̇i.Twelve neighborhood members are randomly selected from the 388 neighborhood members. Each member of the manipulator neighborhood is calculated, and all the trajectories of the joint angle varying with time are obtained. Taking joint 2 of the left manipulator as an example, the optimal joint angle trajectories of all the 12 neighborhood members are determined, as shown in Figure5.Figure 5
Twelve optimal joint angle trajectories for joint 2 of the left manipulator.


## 3.3. Feedforward Compensation and Iterative Optimization

A recursive form of the multirigid body dynamic model can be established using the Newton-Euler method as follows:(24)Γi=hθi,θ̇i,θ¨i,v̇cha,v̇tar,i=1,2,.⋯,14,where Γi=fi,miT denotes the force and moment applied to each joint, mi denotes the transferring torque of each joint, and fi denotes the transferring force of each joint.Wheni=0, Γ0=f0,m0T denotes the transferring force and moment between the manipulator base and the chaser body. By calculating and applying these values inversely as the feedforward to the chaser orbit and attitude control, the control variables of the chaser can be updated as follows:
(25)Ut=ut+Γ0t.Since the disturbance is estimated, it cannot be precisely eliminated by using the feedforward compensation usingut. The possible influences may include the deviation of the optimal trajectory, which makes ut not optimal, as shown in Figure 6. The orbit and attitude motion trajectory is assumed as xt0,⋯,xti,⋯,xtf (solid line) by assigning the terminal state xtf when the manipulators do not move. When the manipulators move and finish the capturing, the orbit and attitude motion trajectory is modified to x∗t0,⋯,x∗ti,⋯,x∗tf (dash line) with the terminal state x∗tf. One or more iterations are performed to update ut to make x∗tf closer to xtf. Considering that Δx∗tf=x∗tf−xtf, the iteration stops when Δx∗tf is less than a preset constant value, and the optimal orbit and attitude trajectory x′tf (dot-dash line) is obtained. The corresponding optimal control sequence is Ut=u′t+Γ0t.Figure 6
Iterative optimization process.


## 4. Learning the Optimal Control

A series of RBF networks are chosen to reduce the optimization time consumption while maintaining the solution to be optimal or at least suboptimal. The network input and output are obtained from the optimal solutions and their initial values of the state variables. By assuming the number of neighborhood members asm and the line number of state sequence X of each neighborhood member as n, the number of lines of all the neighborhood members can be defined asm·n. Considering xt and ut as the tth lines of Xt and Ut, respectively, xt and ut together form a sample, where xt is the input and ut is the output.It is unrealistic to directly provide all thext and ut as samples to train a single large neural network, and decoupling the input and output according to the motion characteristics of the chaser is necessary. In this case, the motion of the chaser can be divided into 3 parts: the orbit motion, attitude motion, and manipulator motion. Each part can be further divided according to the dimensions.


### 4.1. Coordinate Transformation

The optimal trajectories and the corresponding optimal controls form the so called “state-action pairs” and are supposed to be used to train networks. However, the RBF networks cannot be directly trained by using these state-action pairs because the state variables are anisotropic, which means that the same relative trajectory has different relative state variables in different directions. To eliminate the anisotropy, it is necessary to transform the state variables and control variables to a new body-based coordinate system.Here, a training coordinate systemSn is defined as shown by the blue line in Figure 7. The x axis of Sn coincides with the x axis of So. The y axis of Sn is the projection of the rotation axis of the ENVISAT on the yoz plane of So. The z axis can be determined using the right-hand rule. The rotation matrix between the training coordinate system and the orbital coordinate system is Lno. By transferring the inputs and outputs to Sn, the anisotropy can be eliminated, and the regulations between the inputs and outputs can be more easily learned.Figure 7
Training coordinate system. The subscript “o” denotes the orbital coordinate systemSo, and the subscript “n” denotes the training coordinate system Sn.


### 4.2. Orbit Control Networks

Three RBF neural networksNFx, NFy, and NFz are designed according to the three-axis motion of the chaser orbit. These networks are supposed to provide the three-axis control forces by providing certain inputs. Therefore, it is vital to determine the factors that affect the decision to implement the control forces. The instant relative position Δrt and relative velocity Δvt play a dominant role in the decision, while the terminal relative position Δrtf and relative velocity Δvtf may influence the decision, especially when the chaser is close to the terminal position. The maximum relative velocity Δvmax affects the decision of the acceleration or deceleration. Hence, rt, vt, rtf, vtf, and Δvmax are selected to form n×5 matrices IFi (i=x,y,z) as the input, and F from ut is selected to form n×1 matrices OFi (i=x,y,z) as the output, which can be formulated as
(26)IFi=Δrit,Δvit,Δritf,Δvitf,ΔvmaxT,i=x,y,z,OFi=FitT,i=x,y,z.The training procedure for the orbit control networks is shown in Figure8. First, IFi and OFi are transferred to the training coordinate system while neglecting the orbital angular velocity. Next, both IFi and OFi are divided into a training dataset and a testing dataset after normalization. The networks are trained using the training dataset and tested using the testing dataset. The trained networks onboard provide the three-axis forces by inputting the appropriate values from the sensors.Figure 8
Training process of the orbit control networkNFi (i=x,y,z).During the normalization procedure, the inputs and outputs are both normalized to the [-1, 1] interval.MFin~ is defined as the overall dataset of the input and output data after normalization, σFi∈0,1 as the test ratio, j=σFi·n·m as the line number of the training data, and k=1−j as the line number of the test data. j lines are randomly selected from MFin~ to form the training dataset MFin_train~. The remaining k lines are used to form the test dataset MFin_test~. The NFi network can be trained using MFin_train~ and tested using MFin_test~ after setting an appropriate spread speed and training goal.The output of the trained network when the input isIFin_test~ is OFin~. The deviation between OFin~ and the desired output OFin_test~ is ΔOFin. For every moment tk, the deviation can be written as ΔFFintk. Mean square error (MSE) is used to evaluate the training error of the network, which can be calculated as
(27)εsdFi=∑k=1jΔFin~tk2.By antinormalizingOFin~, the actual output OFi=Fi of the network NFi (i=x,y,z) can be obtained.


### 4.3. Attitude Control Networks

Similarly, the attitude motion can be controlled by three networksNMx, NMy, and NMz. The training process of NMi (i=x,y,z) is shown in Figure 9. The only difference between the orbit control network and the attitude control network is in terms of the input segment, which contains a process to transform qcha and qtar to relative quaternions Δq and relative angular velocities Δω.Figure 9
Training process of the attitude control networkNMi (i=x,y,z).Furthermore,Δq is transferred to the relative Euler angle Δe so that the dimension of xt can be reduced as follows:
(28)xt→Δe=fqtar,qchaΔω=ωtar−ωchax′t=rt,vt,Δet,ΔωtT.Similar to that for the orbit control networks, the selection of the input and output elements for the attitude control networks is determined by analyzing the factors affecting the output moments. Here, the elements selected fromx′t to form the inputs IMx, IMy, and IMz and the outputs OMx, OMy, and OMz are as follows:
(29)IMx=Δext,ΔωxtT,IMy=Δrzt,Δvzt,Δeyt,ΔωytT,IMz=Δryt,Δvyt,Δezt,ΔωztT,OMi=MitTi=x,y,z.Note that the attitude control is related not only to the relative attitude states, i.e.,Δeit and Δωit, but also to the relative position and velocity states, i.e., Δrit and Δvit in Equation (29), because the orbit and attitude are coupled in relative motion.The other procedures are similar to those of the orbit control networks and are not repeated here.


### 4.4. Manipulator Joint Control Networks

Joints 2 and 4 are selected to be controlled by the networks because they exert a considerable influence on the manipulator shape. For a 7-DOF manipulator, once a redundant joint (joint 2 or joint 4) is determined, the other joint angles can be calculated. The left and right manipulator control networks can be designed asNl4 and Nr4, respectively. The training process of Nl4 is shown in Figure 10 as an example.Figure 10
Training process of the manipulator control networkNl4.The inputsIl4 and Ir4 and outputs Ol4 and Or4 no longer need to be transferred to a certain coordinate system, and they are designed as follows:
(30)Il2=Il4=Ir2=Ir4=rxt,rzt,vxt,vzt,Δext,Δeyt,Δezt,Δωxt,Δωyt,ΔωztT,Ol2=θ2tT,Ol4=θ4tT,Or2=θ8tT,Or4=θ10tT.Note that the outputs of the manipulators are the joint anglesθi, and the joint torques T can be obtained using Equation (23) after determining θ̇i.The other procedures are similar to those of the orbit control networks and are not repeated here.


### 4.5. Training Process

The training parameters of the orbit, attitude, and manipulator control networks are listed in Table3. In our investigation, the training performances using the parameters different from those listed in Table 3 either match the training performance using the parameters in the table or exhibit a degradation. The training process of NFx with 50 iterations is shown in the left part of Figure 8, and the training process of the left manipulator joint 4 control network Nl4 with 22 iterations is shown in the right part of Figure 11.Table 3
Training parameters of the control networks.
Network
Spread speed
Training goal
NFx
0.1
0.005
NFy
0.1
0.008
NFz
0.1
0.012
NMx
0.1
0.003
NMy
0.1
0.005
NMz
0.1
0.003
Nl2
1
0.001
Nl4
1
0.001
Nr2
1
0.001
Nr4
1
0.001Figure 11
Training process for the RBF networks. The gray dashed line denotes the preset training goal. (a) Training process forNFx. (b) Training process for Nl4.
(a)
(b)After training, the MSEs of the orbit control networksNFx, NFy, and NFz are 0.499%, 0.788%, and 0.119%, respectively. The MSEs of the attitude control networks NMx, NMy, and NMz are 0.299%, 0.499%, and 0.299%, respectively. The MSEs of the manipulator control networks Nl2, Nl4, Nr2, and Nr4 are 0.084%, 0.099%, 0.085%, and 0.096%, respectively.The trained networks can be completely described by a series of parameters containing the input dimensionni, output dimension no, number of neurons nn, input layer weight Iw, output layer weight Lw, hidden layer weights b1 and b2, maximum value Oomax, and minimum value Oomin of the normalization range. Once the networks are well trained, the network parameters are loaded into the control computer of the chaser to control its orbit, attitude, and manipulators.


## 4.1. Coordinate Transformation

The optimal trajectories and the corresponding optimal controls form the so called “state-action pairs” and are supposed to be used to train networks. However, the RBF networks cannot be directly trained by using these state-action pairs because the state variables are anisotropic, which means that the same relative trajectory has different relative state variables in different directions. To eliminate the anisotropy, it is necessary to transform the state variables and control variables to a new body-based coordinate system.Here, a training coordinate systemSn is defined as shown by the blue line in Figure 7. The x axis of Sn coincides with the x axis of So. The y axis of Sn is the projection of the rotation axis of the ENVISAT on the yoz plane of So. The z axis can be determined using the right-hand rule. The rotation matrix between the training coordinate system and the orbital coordinate system is Lno. By transferring the inputs and outputs to Sn, the anisotropy can be eliminated, and the regulations between the inputs and outputs can be more easily learned.Figure 7
Training coordinate system. The subscript “o” denotes the orbital coordinate systemSo, and the subscript “n” denotes the training coordinate system Sn.


## 4.2. Orbit Control Networks

Three RBF neural networksNFx, NFy, and NFz are designed according to the three-axis motion of the chaser orbit. These networks are supposed to provide the three-axis control forces by providing certain inputs. Therefore, it is vital to determine the factors that affect the decision to implement the control forces. The instant relative position Δrt and relative velocity Δvt play a dominant role in the decision, while the terminal relative position Δrtf and relative velocity Δvtf may influence the decision, especially when the chaser is close to the terminal position. The maximum relative velocity Δvmax affects the decision of the acceleration or deceleration. Hence, rt, vt, rtf, vtf, and Δvmax are selected to form n×5 matrices IFi (i=x,y,z) as the input, and F from ut is selected to form n×1 matrices OFi (i=x,y,z) as the output, which can be formulated as
(26)IFi=Δrit,Δvit,Δritf,Δvitf,ΔvmaxT,i=x,y,z,OFi=FitT,i=x,y,z.The training procedure for the orbit control networks is shown in Figure8. First, IFi and OFi are transferred to the training coordinate system while neglecting the orbital angular velocity. Next, both IFi and OFi are divided into a training dataset and a testing dataset after normalization. The networks are trained using the training dataset and tested using the testing dataset. The trained networks onboard provide the three-axis forces by inputting the appropriate values from the sensors.Figure 8
Training process of the orbit control networkNFi (i=x,y,z).During the normalization procedure, the inputs and outputs are both normalized to the [-1, 1] interval.MFin~ is defined as the overall dataset of the input and output data after normalization, σFi∈0,1 as the test ratio, j=σFi·n·m as the line number of the training data, and k=1−j as the line number of the test data. j lines are randomly selected from MFin~ to form the training dataset MFin_train~. The remaining k lines are used to form the test dataset MFin_test~. The NFi network can be trained using MFin_train~ and tested using MFin_test~ after setting an appropriate spread speed and training goal.The output of the trained network when the input isIFin_test~ is OFin~. The deviation between OFin~ and the desired output OFin_test~ is ΔOFin. For every moment tk, the deviation can be written as ΔFFintk. Mean square error (MSE) is used to evaluate the training error of the network, which can be calculated as
(27)εsdFi=∑k=1jΔFin~tk2.By antinormalizingOFin~, the actual output OFi=Fi of the network NFi (i=x,y,z) can be obtained.


## 4.3. Attitude Control Networks

Similarly, the attitude motion can be controlled by three networksNMx, NMy, and NMz. The training process of NMi (i=x,y,z) is shown in Figure 9. The only difference between the orbit control network and the attitude control network is in terms of the input segment, which contains a process to transform qcha and qtar to relative quaternions Δq and relative angular velocities Δω.Figure 9
Training process of the attitude control networkNMi (i=x,y,z).Furthermore,Δq is transferred to the relative Euler angle Δe so that the dimension of xt can be reduced as follows:
(28)xt→Δe=fqtar,qchaΔω=ωtar−ωchax′t=rt,vt,Δet,ΔωtT.Similar to that for the orbit control networks, the selection of the input and output elements for the attitude control networks is determined by analyzing the factors affecting the output moments. Here, the elements selected fromx′t to form the inputs IMx, IMy, and IMz and the outputs OMx, OMy, and OMz are as follows:
(29)IMx=Δext,ΔωxtT,IMy=Δrzt,Δvzt,Δeyt,ΔωytT,IMz=Δryt,Δvyt,Δezt,ΔωztT,OMi=MitTi=x,y,z.Note that the attitude control is related not only to the relative attitude states, i.e.,Δeit and Δωit, but also to the relative position and velocity states, i.e., Δrit and Δvit in Equation (29), because the orbit and attitude are coupled in relative motion.The other procedures are similar to those of the orbit control networks and are not repeated here.


## 4.4. Manipulator Joint Control Networks

Joints 2 and 4 are selected to be controlled by the networks because they exert a considerable influence on the manipulator shape. For a 7-DOF manipulator, once a redundant joint (joint 2 or joint 4) is determined, the other joint angles can be calculated. The left and right manipulator control networks can be designed asNl4 and Nr4, respectively. The training process of Nl4 is shown in Figure 10 as an example.Figure 10
Training process of the manipulator control networkNl4.The inputsIl4 and Ir4 and outputs Ol4 and Or4 no longer need to be transferred to a certain coordinate system, and they are designed as follows:
(30)Il2=Il4=Ir2=Ir4=rxt,rzt,vxt,vzt,Δext,Δeyt,Δezt,Δωxt,Δωyt,ΔωztT,Ol2=θ2tT,Ol4=θ4tT,Or2=θ8tT,Or4=θ10tT.Note that the outputs of the manipulators are the joint anglesθi, and the joint torques T can be obtained using Equation (23) after determining θ̇i.The other procedures are similar to those of the orbit control networks and are not repeated here.


## 4.5. Training Process

The training parameters of the orbit, attitude, and manipulator control networks are listed in Table3. In our investigation, the training performances using the parameters different from those listed in Table 3 either match the training performance using the parameters in the table or exhibit a degradation. The training process of NFx with 50 iterations is shown in the left part of Figure 8, and the training process of the left manipulator joint 4 control network Nl4 with 22 iterations is shown in the right part of Figure 11.Table 3
Training parameters of the control networks.
Network
Spread speed
Training goal
NFx
0.1
0.005
NFy
0.1
0.008
NFz
0.1
0.012
NMx
0.1
0.003
NMy
0.1
0.005
NMz
0.1
0.003
Nl2
1
0.001
Nl4
1
0.001
Nr2
1
0.001
Nr4
1
0.001Figure 11
Training process for the RBF networks. The gray dashed line denotes the preset training goal. (a) Training process forNFx. (b) Training process for Nl4.
(a)
(b)After training, the MSEs of the orbit control networksNFx, NFy, and NFz are 0.499%, 0.788%, and 0.119%, respectively. The MSEs of the attitude control networks NMx, NMy, and NMz are 0.299%, 0.499%, and 0.299%, respectively. The MSEs of the manipulator control networks Nl2, Nl4, Nr2, and Nr4 are 0.084%, 0.099%, 0.085%, and 0.096%, respectively.The trained networks can be completely described by a series of parameters containing the input dimensionni, output dimension no, number of neurons nn, input layer weight Iw, output layer weight Lw, hidden layer weights b1 and b2, maximum value Oomax, and minimum value Oomin of the normalization range. Once the networks are well trained, the network parameters are loaded into the control computer of the chaser to control its orbit, attitude, and manipulators.


## 5. Results and Discussion

### 5.1. Simulation Case Design

As discussed above, the ENVISAT can be regarded as a5m×10m×5m cube by neglecting its solar array and other deployed devices. The geometric center of the ENVISAT, mass center of the ENVISAT, and origin of Stb are considered to coincide for simplicity. According to Gómez and Walker [45], the mass of the ENVISAT is 8211 kg, the moment of inertia is (124801.21, 16979.74, 129180.25) kg·m2, and the product of inertia is 0 kg·m2. The z axis is the maximum inertia principal axis of the satellite, and ENVISAT is assumed to rotate at an angular velocity of 0.8631°/s [4]. At some point, the instantaneous angular velocity of ENVISAT is assumed to be (0.15, 0, 0.85)°/s and the corresponding nutation angle is 9.6752°, which satisfies the conclusion of Pittet et al. [5]. The instantaneous attitude of ENVISAT is assumed to be (0, 78, 86)° as if it is “lying” on its side in the orbit. According to Section 2, only the +y surface of the ENVISAT is suitable for capture. Therefore, (0, 5, -0.2) m and (0, 5, 0.2) m expressed in Stb are selected as the two desired capture points on the +y surface of the ENVISAT.A small chaser is designed as a0.88m×1m×1m cuboid with a length of 1.5 m after installing the two unfolded manipulators. The mass of the chaser is 100 kg. The moment of inertia of the chaser is (100, 100, 100) kg·m2, and the product of inertia is 0 kg·m2. The geometric center of the chaser is fixed with the mass center of the chaser, and it coincides with the origin of Scb. The two manipulators are installed at (0.49, -0.37, 0) m and (0.49, 0.37, 0) m on the +x surface (front surface) of the chaser. The mass and inertia parameters of the manipulator are listed in Table 4 by referring to Papadopoulos and Dubowsky [46] who considered the manipulator as a cylinder for simplicity.Table 4
Mass and inertia parameters of the manipulator.
Linki
Length (m)
Mass (kg)
Moment of inertia (kg·m2)
0
0
0
(0, 0, 0)
1
0.045
0.08
6.35,6.35,10×10−5
2
0.075
0.3
3.28,3.28,3.75×10−4
3
0.58
2.32
6.65,6.65,2.90×10−2
4
0.58
2.32
6.65,6.65,2.90×10−2
5
0.075
0.3
3.28,3.28,3.75×10−4
6
0.045
0.08
6.35,6.35,10×10−5
7
0.1
0.4
5.83,5.83,5×10−4The maximum thrust force of the chaser is 15 N, the maximum control moment is 15 N·m, and the specific impulse is 300 m·s-1. The maximum joint torque is 40 Nm, and the maximum rotation angular velocity of each joint is 20 rpm.An overall of six cases listed in Table5 are designed to evaluate the performance and verify the effectiveness of the network controller. The first three cases are designed to evaluate the performance of the network controller compared to a tracking control strategy of the optimal trajectories without considering measurement uncertainty. The chaser is locating at three different positions expressed in Stb relative to the ENVISAT near (0, 15, 0) m at the initial time. The latter three cases are designed to verify the effectiveness of the network controller under different nutation situations with initial measurement uncertainty. The chaser is locating at the same position expressed in Stb relative to the ENVISAT. In the latter three cases, the ENVISAT has a slightly smaller absolute rotation angular velocity (0.8631°/s compared to 0.8668°/s in the previous three cases). Case 4 simulates the situation in which the ENVISAT has nutation, although it is not considered and measured. Cases 5 and 6 simulate the situations in which the ENVISAT has nutation, but it is not measured accurately at the initial time. The real value of ωt in the table is the instantaneous angular velocity of the ENVISAT relative to the inertial coordinate system expressed in Stb at the initial time.Table 5
The 6 designed simulation cases.
Cases
Realωt (°/s)
Initial measuredωt (°/s)
Initial relativer (m)
1
(0.17, 0.0, 0.85)
(0.17, 0.0, 0.85)
(0.124, 15.765, -0.417)
2
(0.347, 14.871, -0.308)
3
(-0.423, 15.364, -0.124)
4
(0.15, 0.0, 0.85)
(0.0, 0.0, 0.85)
(-0.218, 14.964, -1.024)
5
(0.13, 0.02, 0.85)
6
(0.17, -0.01, 0.85)In each case, two different control strategies are applied on the chaser: the traditional control strategy that plans at first and then tracks the optimal control result and the network control strategy proposed in this paper. The following performance indices are used to evaluate the mission performance, they are(i)
Jt: the overall mission duration, defined by Equation (8)
(ii)
Δrend=∣Δrend∣: the relative position between the manipulator end effector and the desired capture point within the capture plane at the moment of impact
(iii)
Jv: the relative velocity between the manipulator end effector and the desired capture point at the moment of impact, defined by Equation (9)
(iv)
Jf: fuel consumption defined by Equation (7)


### 5.2. Capture without Measurement Uncertainty

The simulation results of the first three cases are presented in Table6. The average Jt of the traditional control strategy is 34.849 s, and that of the network method is 38.828 s. The average Jv of the traditional control strategy is 0.036 m/s, and that of the network method is 0.027 m/s. The average Jf of the traditional control strategy is 0.517 kg, and that of the network method is 0.586 kg. Considering that kf is 2/3, kt is 1, and kv is 100 (considerably larger to show the importance of the manipulator end impact and to balance the small magnitude of Jv), the overall J of the tracking control strategy is 27.350, and that of the network controller is 29.171. Compared to the tracking control strategy, the performance loss of the network controller is approximately 6.6%.Table 6
Simulation results for the first three verification cases.
Index
Case 1
Case 2
Case 3
Tracking
Network
Tracking
Network
Tracking
Network
Jt (s)
33.361
38.833
35.593
40.129
35.953
37.521
LeftΔrend (m)
0.030
0.100
0.051
0.053
0.056
0.053
RightΔrend (m)
0.086
0.064
0.065
0.047
0.066
0.055
LeftJv (m/s)
0.070
0.008
0.043
0.029
0.010
0.027
RightJv (m/s)
0.015
0.001
0.019
0.056
0.057
0.043
Jf (kg)
0.466
0.594
0.534
0.602
0.552
0.562Since there is no measurement uncertainty in the first three cases, the optimal trajectories planned offline are the true optimal solutions. The tracking result of those optimal trajectories can be regarded as the standard to evaluate the performance of the network controller proposed in this paper. The first three cases indicate that there exist approximately 6.6% performance loss of the network controller. That is, the network controller provides a suboptimal control close to the true optimal control solution. The losses in performance occur due to the following three aspects:(1)
Samples. Samples of inferior quality lead to poor learned networks. Unlike other machine learning problems in which the samples are commonly selected from open source datasets, the samples of the approach and capture problem are generated using the GPM and designed algorithms. Thus, the appropriate initial conditions of the samples contribute to the GPM planning stage, and necessary elimination of samples of inferior quality can improve the training accuracy.
(2)
Training errors (commonly evaluated considering the MSE) that reflect the gaps between the output of the networks and the test datasets. The training error can be reduced by selecting the appropriate network parameters, which can be realized with experience by using the Monte Carlo method.
(3)
Generalization errors that reflect the prediction precision to the unknown flight conditions. The generalization errors are generally related to the network types and structures. The RBF network is a shallow network, and the number of neurons is determined automatically in the training process. The division of the networks according to the dimensions is based on human intelligence. Further research on different kinds of networks and different division methods is in progress.


### 5.3. Capture with Measurement Uncertainty

The simulation results of the latter three cases are presented in Table7. The average Jt of the traditional control strategy is 36.502 s, and that for the network methods is 41.686 s. The average Jv of the traditional control strategy is 0.124 m/s, and that for the network controller is 0.045 m/s. The average Jf of the traditional control strategy is 0.517 kg, and that for the network controller is 0.834 kg. Considering that kf is 2/3, kt is 1, and kv is 100, the overall J of the traditional control strategy is 37.252, and that of the network controller is 33.125.Table 7
Simulation results for the latter three verification cases.
Index
Case 4
Case 5
Case 6
Tracking
Network
Tracking
Network
Tracking
Network
Jt (s)
31.857
41.281
40.273
41.809
37.377
41.969
LeftΔrend (m)
0.088
0.061
0.055
0.060
0.043
0.025
RightΔrend (m)
0.057
0.053
0.128
0.056
0.054
0.053
LeftJv (m/s)
0.086
0.041
0.077
0.077
0.145
0.074
RightJv (m/s)
0.170
0.023
0.114
0.018
0.154
0.039
Jf (kg)
0.577
0.897
0.485
0.808
0.490
0.797Since there exist measurement errors ofωt at the initial time, when the nutation is not measured accurately, the tracking of the offline planned trajectories results in larger Δrend and Jv (bold numbers in Table 7). By using the network controller, the average Jv is 63.71% lower than that of the tracking control strategy, which means that the capture mission has a lower risk to produce debris. On the contrary, the network controller consumes more time and fuel than those controlled by the traditional control strategy, which can be seen from the average values of Jt and Jf.Taking the simulation results of case 5 as an example, the relative orbit varying with time expressed inStb is shown in Figure 12(a). The joint 2 and joint 4 angles of the left manipulator varying with time are shown in Figure 12(b).Figure 12
Relative orbit motion and manipulator joint motions in case 5. (a) Relative position expressed inStb. (b) Angles of joint 2 and joint 4 of the left manipulator.
(a)
(b)The vertical velocities of both the left and right manipulator end effectors relative to the desired capture points in all the latter three cases under different control strategies are shown in Figure13. When the nutation was neglected, the chaser’s manipulator end effectors hit the ENVISAT at a relatively high speed when controlled by the traditional control strategy, as depicted using the red lines in case 4. When the nutation is taken into account but the measurement is not accurate (cases 5 and 6), the chaser’s manipulator end effectors still heavily impact the ENVISAT when controlled by the traditional tracking control strategy (red lines), compared to when the chaser is controlled by the networks (blue lines). This indicates that the network controller can conduct a softer capture with less impact.Figure 13
Vertical velocities of both the left and right manipulator end effectors relative to the desired capture points. The solid lines represent the left manipulator, and the dotted lines denote the right manipulator. The red lines denote the chaser controlled by the traditional control strategy, and blue lines denote the chaser controlled by the networks. All the lines are terminated when the manipulators hit the target. (a) Case 4. (b) Case 5. (c) Case 6.
(a)
(b)
(c)The relative positions between both the manipulator end effectors and the desired capture points within the capture plane at the moment of impact are shown in Figure14. The dots closer to the origin indicated a smaller capture position deviation. It can be seen from the image that in all the latter three cases, the chaser controlled by the networks exhibits a smaller capture position deviation (blue dots) compared to that of the chaser controlled using the traditional control strategy (red dots). This indicates that the network controller can conduct a more accurate capture action.Figure 14
Relative positions between both the manipulator end effectors and the desired capture points within the capture plane at the moment of impact. The red dots denote the chasers controlled using the traditional control strategy, and the blue dots denote the chaser controlled by the networks. The circular dots denote the left manipulator end effectors, and the square dots denote the right manipulator end effectors. (a) Case 4. (b) Case 5. (c) Case 6.
(a)
(b)
(c)The total flight history of the chaser in case 5 is shown in Figure15. Multiple flight moments from the initial time to the final time are shown in Figure 15(a). At the initial time, the chaser is at the right of the ENVISAT, and the ENVISAT rotates approximately around its z axis (not exactly around the z axis because of nutation). The relative flight trajectory between the chaser and the ENVISAT is plotted in the white line. Figure 15(b) shows the deployment, pointing, and capturing procedure of the two manipulators.Figure 15
Flight history of the chaser in case 5. (a) The approach and capture process in the orbital coordinate system. (b) The attitude motion of the chaser and the deployment of the two manipulators from the negativez direction of the orbital coordinate system.
(a)
(b)


## 5.1. Simulation Case Design

As discussed above, the ENVISAT can be regarded as a5m×10m×5m cube by neglecting its solar array and other deployed devices. The geometric center of the ENVISAT, mass center of the ENVISAT, and origin of Stb are considered to coincide for simplicity. According to Gómez and Walker [45], the mass of the ENVISAT is 8211 kg, the moment of inertia is (124801.21, 16979.74, 129180.25) kg·m2, and the product of inertia is 0 kg·m2. The z axis is the maximum inertia principal axis of the satellite, and ENVISAT is assumed to rotate at an angular velocity of 0.8631°/s [4]. At some point, the instantaneous angular velocity of ENVISAT is assumed to be (0.15, 0, 0.85)°/s and the corresponding nutation angle is 9.6752°, which satisfies the conclusion of Pittet et al. [5]. The instantaneous attitude of ENVISAT is assumed to be (0, 78, 86)° as if it is “lying” on its side in the orbit. According to Section 2, only the +y surface of the ENVISAT is suitable for capture. Therefore, (0, 5, -0.2) m and (0, 5, 0.2) m expressed in Stb are selected as the two desired capture points on the +y surface of the ENVISAT.A small chaser is designed as a0.88m×1m×1m cuboid with a length of 1.5 m after installing the two unfolded manipulators. The mass of the chaser is 100 kg. The moment of inertia of the chaser is (100, 100, 100) kg·m2, and the product of inertia is 0 kg·m2. The geometric center of the chaser is fixed with the mass center of the chaser, and it coincides with the origin of Scb. The two manipulators are installed at (0.49, -0.37, 0) m and (0.49, 0.37, 0) m on the +x surface (front surface) of the chaser. The mass and inertia parameters of the manipulator are listed in Table 4 by referring to Papadopoulos and Dubowsky [46] who considered the manipulator as a cylinder for simplicity.Table 4
Mass and inertia parameters of the manipulator.
Linki
Length (m)
Mass (kg)
Moment of inertia (kg·m2)
0
0
0
(0, 0, 0)
1
0.045
0.08
6.35,6.35,10×10−5
2
0.075
0.3
3.28,3.28,3.75×10−4
3
0.58
2.32
6.65,6.65,2.90×10−2
4
0.58
2.32
6.65,6.65,2.90×10−2
5
0.075
0.3
3.28,3.28,3.75×10−4
6
0.045
0.08
6.35,6.35,10×10−5
7
0.1
0.4
5.83,5.83,5×10−4The maximum thrust force of the chaser is 15 N, the maximum control moment is 15 N·m, and the specific impulse is 300 m·s-1. The maximum joint torque is 40 Nm, and the maximum rotation angular velocity of each joint is 20 rpm.An overall of six cases listed in Table5 are designed to evaluate the performance and verify the effectiveness of the network controller. The first three cases are designed to evaluate the performance of the network controller compared to a tracking control strategy of the optimal trajectories without considering measurement uncertainty. The chaser is locating at three different positions expressed in Stb relative to the ENVISAT near (0, 15, 0) m at the initial time. The latter three cases are designed to verify the effectiveness of the network controller under different nutation situations with initial measurement uncertainty. The chaser is locating at the same position expressed in Stb relative to the ENVISAT. In the latter three cases, the ENVISAT has a slightly smaller absolute rotation angular velocity (0.8631°/s compared to 0.8668°/s in the previous three cases). Case 4 simulates the situation in which the ENVISAT has nutation, although it is not considered and measured. Cases 5 and 6 simulate the situations in which the ENVISAT has nutation, but it is not measured accurately at the initial time. The real value of ωt in the table is the instantaneous angular velocity of the ENVISAT relative to the inertial coordinate system expressed in Stb at the initial time.Table 5
The 6 designed simulation cases.
Cases
Realωt (°/s)
Initial measuredωt (°/s)
Initial relativer (m)
1
(0.17, 0.0, 0.85)
(0.17, 0.0, 0.85)
(0.124, 15.765, -0.417)
2
(0.347, 14.871, -0.308)
3
(-0.423, 15.364, -0.124)
4
(0.15, 0.0, 0.85)
(0.0, 0.0, 0.85)
(-0.218, 14.964, -1.024)
5
(0.13, 0.02, 0.85)
6
(0.17, -0.01, 0.85)In each case, two different control strategies are applied on the chaser: the traditional control strategy that plans at first and then tracks the optimal control result and the network control strategy proposed in this paper. The following performance indices are used to evaluate the mission performance, they are(i)
Jt: the overall mission duration, defined by Equation (8)
(ii)
Δrend=∣Δrend∣: the relative position between the manipulator end effector and the desired capture point within the capture plane at the moment of impact
(iii)
Jv: the relative velocity between the manipulator end effector and the desired capture point at the moment of impact, defined by Equation (9)
(iv)
Jf: fuel consumption defined by Equation (7)


## 5.2. Capture without Measurement Uncertainty

The simulation results of the first three cases are presented in Table6. The average Jt of the traditional control strategy is 34.849 s, and that of the network method is 38.828 s. The average Jv of the traditional control strategy is 0.036 m/s, and that of the network method is 0.027 m/s. The average Jf of the traditional control strategy is 0.517 kg, and that of the network method is 0.586 kg. Considering that kf is 2/3, kt is 1, and kv is 100 (considerably larger to show the importance of the manipulator end impact and to balance the small magnitude of Jv), the overall J of the tracking control strategy is 27.350, and that of the network controller is 29.171. Compared to the tracking control strategy, the performance loss of the network controller is approximately 6.6%.Table 6
Simulation results for the first three verification cases.
Index
Case 1
Case 2
Case 3
Tracking
Network
Tracking
Network
Tracking
Network
Jt (s)
33.361
38.833
35.593
40.129
35.953
37.521
LeftΔrend (m)
0.030
0.100
0.051
0.053
0.056
0.053
RightΔrend (m)
0.086
0.064
0.065
0.047
0.066
0.055
LeftJv (m/s)
0.070
0.008
0.043
0.029
0.010
0.027
RightJv (m/s)
0.015
0.001
0.019
0.056
0.057
0.043
Jf (kg)
0.466
0.594
0.534
0.602
0.552
0.562Since there is no measurement uncertainty in the first three cases, the optimal trajectories planned offline are the true optimal solutions. The tracking result of those optimal trajectories can be regarded as the standard to evaluate the performance of the network controller proposed in this paper. The first three cases indicate that there exist approximately 6.6% performance loss of the network controller. That is, the network controller provides a suboptimal control close to the true optimal control solution. The losses in performance occur due to the following three aspects:(1)
Samples. Samples of inferior quality lead to poor learned networks. Unlike other machine learning problems in which the samples are commonly selected from open source datasets, the samples of the approach and capture problem are generated using the GPM and designed algorithms. Thus, the appropriate initial conditions of the samples contribute to the GPM planning stage, and necessary elimination of samples of inferior quality can improve the training accuracy.
(2)
Training errors (commonly evaluated considering the MSE) that reflect the gaps between the output of the networks and the test datasets. The training error can be reduced by selecting the appropriate network parameters, which can be realized with experience by using the Monte Carlo method.
(3)
Generalization errors that reflect the prediction precision to the unknown flight conditions. The generalization errors are generally related to the network types and structures. The RBF network is a shallow network, and the number of neurons is determined automatically in the training process. The division of the networks according to the dimensions is based on human intelligence. Further research on different kinds of networks and different division methods is in progress.


## 5.3. Capture with Measurement Uncertainty

The simulation results of the latter three cases are presented in Table7. The average Jt of the traditional control strategy is 36.502 s, and that for the network methods is 41.686 s. The average Jv of the traditional control strategy is 0.124 m/s, and that for the network controller is 0.045 m/s. The average Jf of the traditional control strategy is 0.517 kg, and that for the network controller is 0.834 kg. Considering that kf is 2/3, kt is 1, and kv is 100, the overall J of the traditional control strategy is 37.252, and that of the network controller is 33.125.Table 7
Simulation results for the latter three verification cases.
Index
Case 4
Case 5
Case 6
Tracking
Network
Tracking
Network
Tracking
Network
Jt (s)
31.857
41.281
40.273
41.809
37.377
41.969
LeftΔrend (m)
0.088
0.061
0.055
0.060
0.043
0.025
RightΔrend (m)
0.057
0.053
0.128
0.056
0.054
0.053
LeftJv (m/s)
0.086
0.041
0.077
0.077
0.145
0.074
RightJv (m/s)
0.170
0.023
0.114
0.018
0.154
0.039
Jf (kg)
0.577
0.897
0.485
0.808
0.490
0.797Since there exist measurement errors ofωt at the initial time, when the nutation is not measured accurately, the tracking of the offline planned trajectories results in larger Δrend and Jv (bold numbers in Table 7). By using the network controller, the average Jv is 63.71% lower than that of the tracking control strategy, which means that the capture mission has a lower risk to produce debris. On the contrary, the network controller consumes more time and fuel than those controlled by the traditional control strategy, which can be seen from the average values of Jt and Jf.Taking the simulation results of case 5 as an example, the relative orbit varying with time expressed inStb is shown in Figure 12(a). The joint 2 and joint 4 angles of the left manipulator varying with time are shown in Figure 12(b).Figure 12
Relative orbit motion and manipulator joint motions in case 5. (a) Relative position expressed inStb. (b) Angles of joint 2 and joint 4 of the left manipulator.
(a)
(b)The vertical velocities of both the left and right manipulator end effectors relative to the desired capture points in all the latter three cases under different control strategies are shown in Figure13. When the nutation was neglected, the chaser’s manipulator end effectors hit the ENVISAT at a relatively high speed when controlled by the traditional control strategy, as depicted using the red lines in case 4. When the nutation is taken into account but the measurement is not accurate (cases 5 and 6), the chaser’s manipulator end effectors still heavily impact the ENVISAT when controlled by the traditional tracking control strategy (red lines), compared to when the chaser is controlled by the networks (blue lines). This indicates that the network controller can conduct a softer capture with less impact.Figure 13
Vertical velocities of both the left and right manipulator end effectors relative to the desired capture points. The solid lines represent the left manipulator, and the dotted lines denote the right manipulator. The red lines denote the chaser controlled by the traditional control strategy, and blue lines denote the chaser controlled by the networks. All the lines are terminated when the manipulators hit the target. (a) Case 4. (b) Case 5. (c) Case 6.
(a)
(b)
(c)The relative positions between both the manipulator end effectors and the desired capture points within the capture plane at the moment of impact are shown in Figure14. The dots closer to the origin indicated a smaller capture position deviation. It can be seen from the image that in all the latter three cases, the chaser controlled by the networks exhibits a smaller capture position deviation (blue dots) compared to that of the chaser controlled using the traditional control strategy (red dots). This indicates that the network controller can conduct a more accurate capture action.Figure 14
Relative positions between both the manipulator end effectors and the desired capture points within the capture plane at the moment of impact. The red dots denote the chasers controlled using the traditional control strategy, and the blue dots denote the chaser controlled by the networks. The circular dots denote the left manipulator end effectors, and the square dots denote the right manipulator end effectors. (a) Case 4. (b) Case 5. (c) Case 6.
(a)
(b)
(c)The total flight history of the chaser in case 5 is shown in Figure15. Multiple flight moments from the initial time to the final time are shown in Figure 15(a). At the initial time, the chaser is at the right of the ENVISAT, and the ENVISAT rotates approximately around its z axis (not exactly around the z axis because of nutation). The relative flight trajectory between the chaser and the ENVISAT is plotted in the white line. Figure 15(b) shows the deployment, pointing, and capturing procedure of the two manipulators.Figure 15
Flight history of the chaser in case 5. (a) The approach and capture process in the orbital coordinate system. (b) The attitude motion of the chaser and the deployment of the two manipulators from the negativez direction of the orbital coordinate system.
(a)
(b)


## 6. Conclusions

This paper proposes a real-time optimal control strategy to capture the ENVISAT in the presence of nutation and measurement uncertainty using artificial intelligence. From the simulation results, the following observations have been obtained. Firstly, in ideal cases, the capture performance of our method is comparable with the tracking control approach. Secondly, the averageJv (the relative velocity between the manipulator end effector) is 63.71% lower than that of the tracking control approach, implying the capture mission has a lower risk to produce debris. Thirdly, the solution of the tracking control problem requires an average computational cost of 1128 seconds to generate the policies while our approach is accomplishable in real time. It is noted that the inference time of the neural network is 12.97 milliseconds on ARM (STM32F417) that can be omitted in this case. As a result, our proposed control strategy is particularly suitable to solve a variety of other control tasks where performance and efficiency are both required.



---

*Source: 8165147-2020-01-30.xml*