Real-Time Optimal Control for Irregular Asteroid Landings Using Deep Neural Networks

TL;DR

This paper introduces a deep neural network-based control system for real-time, precise asteroid landings, significantly improving computational efficiency and robustness compared to traditional methods.

Contribution

It develops a novel DNN-based approach that approximates asteroid gravity fields and learns optimal control actions, enabling real-time, robust asteroid landing control.

Findings

DNNs accurately approximate irregular asteroid gravity fields

The control system achieves real-time optimal landings in simulations

The approach demonstrates high robustness to initial conditions

Abstract

Precise soft landings on asteroids are central to many deep space missions for surface exploration and resource exploitation. To improve the autonomy and intelligence of landing control, a real-time optimal control approach is proposed using deep neural networks (DNN) for asteroid landing problems wherein the developed DNN-based landing controller is capable of steering the lander to a preselected landing site with high robustness to initial conditions. First, to significantly reduce the time consumption of gravity calculation, DNNs are used to approximate the irregular gravitational field of the asteroid based on the samples from a polyhedral method. Then, an approximate indirect method is presented to solve the time-optimal landing problems with high computational efficiency by taking advantage of the designed gravity approximation method and a homotopy technique. Furthermore, five DNNs are developed to learn the functional relationship between the state and optimal actions obtained by the approximate indirect method, and the resulting DNNs can generate the optimal control instructions in real time because there is no longer need to solve the optimal landing problems onboard. Finally, a DNN-based landing controller composed of these five DNNs is devised to achieve the real-time optimal control for asteroid landings. Simulation results of the time-optimal landing for Eros are given to substantiate the effectiveness of these techniques and illustrate the real-time performance, control optimality, and robustness of the developed DNN-based optimal landing controller.

Figures

Tables

Table 1. Table 1: The parameter setting of the gravity DNN for 443 Eros

	Hidden Layer	Output Layer
Activation Function	ReLU	Linear
Size	6 Layers/256 Units
LR	0.0001
N	1000

Table 2. Table 2: Parameters setting of the control DNNs

		Activation Function	Size	Loss
$Ne{t}_{\alpha }$	Hidden Layer	ReLU	6 Layer/512 Units	MAE
	Output Layer	Tanh
$Ne{t}_{t_f}$	Hidden Layer	ReLU	6 Layer/256 Units	MSE
	Output Layer	Softplus
$Ne{t}_{\Delta m}$	Hidden Layer	ReLU	6 Layer/256 Units	MSE
	Output Layer	Softplus
LR	0.0001
N	1000

Table 3. Table 3: Time cost comparison

	Polyhedron (C)	DNN (Python)
Computation Efficiency Ratio	36	1
Gravity Calculation	41.3 ms	0.8 ms
Trajectory Optimization without Homotopy	295.2 min	5.7 min
Trajectory Optimization with Homotopy	48.5 min	1.1 min

Table 4. Table 4: Robustness analysis

	Indirect Method				DNN-Based Control
	Position Error (km)		Velocity Error (km/s)		Position Error (km)		Velocity Error (km/s)
	Mean	Max	Mean	Max	Mean	Max	Mean	Max
$-100\%$	3.37	5.28	4.98E-3	6.43E-3	1.25E-1	1.37E-1	2.33E-4	3.96E-4
$-50\%$	1.75	2.78	2.66E-3	3.52E-3	7.72E-2	8.62E-2	6.24E-5	1.88E-4
$-10\%$	3.61E-1	5.81E-1	5.59E-4	7.51E-4	1.81E-2	2.73E-2	1.08E-4	1.37E-4
$0$	0	0	0	0	3.29E-3	1.10E-2	1.51E-5	3.93E-5
$+10\%$	3.66E-1	5.91E-1	5.53E-4	7.38E-4	1.52E-2	1.97E-2	9.71E-5	1.21E-4
$+50\%$	1.735	2.80	1.66E-3	2.45E-3	8.36E-2	9.11E-2	5.08E-4	5.58E-4
$+100\%$	3.405	5.71	3.28E-3	5.20E-3	1.61E-1	1.77E-1	6.32E-4	7.02E-4

Equations

$${\min }_{\mathbf{\mu },\mathbf{\alpha }}{t}_f={\int }_{t_0}^{t_f}1dt$$

$$\begin{array}[]{l}\dot{\bm{r}}=\bm{v}\\ \dot{\bm{v}}+2\bm{\Omega}\times\bm{v}+\bm{\Omega}\times(\bm{\Omega}\times\bm{r})=\bm{G}(\bm{r})+\frac{T_{\max}\mu\bm{\alpha}}{m}\\ \dot{m}=-\frac{T_{\max}\mu}{I_{\mathrm{sp}}g_{0}}\end{array}$$

$$\mathbf{r}(0)={\mathbf{r}}_0,\mathbf{v}(0)={\mathbf{v}}_0,m(0)=m_0$$

$$\mathbf{r}\left(t_f\right)={\mathbf{r}}_f,\mathbf{v}\left(t_f\right)=0$$

$$\dot{\mathbf{v}}+2\mathbf{\Omega }\times \mathbf{v}+\mathbf{\Omega }\times (\mathbf{\Omega }\times \mathbf{r})=\frac{T_{\max }\mu \mathbf{\alpha }}{m}$$

$$\mathbf{G}=-G_u\rho \sum _{e\in edges}{L}_e{\mathbf{E}}_e\cdot {\mathbf{r}}_e+G_u\rho \sum _{f\in faces}{\omega }_f{\mathbf{F}}_f\cdot {\mathbf{r}}_f$$

$$L(\mathbf{\omega })=\frac{1}{N}\sum _{i=1}^N{\Vert Ne{t}_{\mathrm{G}}({\mathbf{r}}_i|\mathbf{\omega })-{\mathbf{G}}_i\Vert }^2$$

$${\mathbf{G}}_{\mathrm{axis}}^{\mathrm{error}}=\frac{Ne{t}_{\mathrm{G}}(\mathbf{r}|\mathbf{\omega })-\mathbf{G}}{{\Vert \mathbf{G}\Vert }_2}\times 100\%$$

$${\mathbf{G}}_{\mathrm{all}}^{\mathrm{error}}=\frac{{\Vert Ne{t}_{\mathrm{G}}(\mathbf{r}|\mathbf{\omega })-\mathbf{G}\Vert }_2}{{\Vert \mathbf{G}\Vert }_2}\times 100\%$$

$$H={\mathbf{\lambda }}_r\cdot \mathbf{v}+{\mathbf{\lambda }}_{\mathbf{v}}\cdot \left(-2\mathbf{\Omega }\times \mathbf{v}-\mathbf{\Omega }\times (\mathbf{\Omega }\times \mathbf{r})+Ne{t}_{\mathrm{G}}(\mathbf{r}|\mathbf{\omega }\right)+\frac{T_{max}\mu }{m}\mathbf{\alpha })+{\lambda }_m\cdot \left(-\frac{T_{max}\mu }{I_{sp}{g}_0}\right)+{\lambda }_0$$

$$\Vert \mathbf{\lambda }(t_0){\Vert }_2=1$$

$$\begin{array}[]{l}\bm{\alpha}^{*}=-\frac{\bm{\lambda_{v}}}{\|\bm{\lambda_{v}}\|_{2}}\\ \mu^{*}=1\end{array}$$

$$\begin{array}[]{l}\dot{\bm{\lambda}}_{r}=-\frac{\partial{H}}{\partial{\bm{r}}}=\bm{\Omega}\times(\bm{\Omega}\times\bm{\lambda_{v}})-\frac{\partial{Net_{\mathrm{G}}(\bm{r}|\bm{\omega})}}{\partial{\bm{r}}}\cdot\bm{\lambda_{v}}\\ \dot{\bm{\lambda}}_{v}=-\frac{\partial{H}}{\partial{\bm{v}}}==-\bm{\lambda_{r}}-2\bm{\Omega}\times\bm{\lambda_{v}}\\ \dot{\bm{\lambda}}_{m}=-\frac{T_{max}\mu\|\bm{\lambda_{v}}\|}{m^{2}}\end{array}$$

$$H(t_f)=0$$

$${\lambda }_m(t_f)=0$$

$$\Phi (\mathbf{z})=[\mathbf{r}(t_f)-{\mathbf{r}}_f,\mathbf{v}(t_f)-{\mathbf{v}}_f,{\lambda }_m(t_f),\Vert \mathbf{\lambda }(t_0){\Vert }_2-1,H(t_f)]=0$$

$$\begin{array}[]{l}\dot{\bm{v}}+2\bm{\Omega}\times\bm{v}+\bm{\Omega}\times(\bm{\Omega}\times\bm{r})=\varepsilon\cdot Net_{\mathrm{G}}(\bm{r}|\bm{\omega})+\frac{T_{\max}\mu\bm{\alpha}}{m}\end{array}$$

$$L(\mathbf{\omega })=\frac{1}{N}\sum _{i=1}^N|Ne{t}_{\alpha }({\mathbf{x}}_i|{\mathbf{\omega }}^{\alpha })-{\alpha }_i^{\ast }|$$

Full text

1 Introduction

Small bodies in the solar system such as asteroids and comets contain high scientific values, and various exploration missions have been carried out for exploring these small bodies [1]. Landing a spacecraft on an asteroid for collecting high-quality resolution data and soil samples is a crucial step in the asteroid sample return missions, such as the finished NEAR [2] and Hayabusa 1 [3] missions, and the ongoing Hayabusa 2 [4] and OSIRIS-REx [4] missions. However, there is a big challenge for asteroid landing due to the remote distance, irregular shapes, rapidly variable gravitational fields, and low thrust of the vehicle [5]. Considering the effects of these disturbances and uncertainties, the landing guidance and control approaches with high autonomy and reliability should be studied.

The majority of the existing powered landing systems use two separate and independently systems for guidance and control, wherein a guidance system generates optimal trajectories that connect the lander’s current state with its target state, and a control system is in charge of tracking these trajectories by determining the thrust magnitudes and directions [6]. Referring to this nominal trajectory tracking pattern, there are different derivative versions which are mainly focused on the improvement of trajectory optimization and trajectory tracking techniques. In generally, the asteroid landing problem has been formulated as an optimal control problem (OCP) and traditionally solved by the direct and indirect methods. Benefiting from the development of computers in recent years, pseudo-spectral method [7] and convex optimization techniques [8, 9] are two effective direct methods that transform the original trajectory optimization problems into nonlinear programming problems (NLPs). Based on the study in [8], the convex optimization method was improved in [9] for solving the time-optimal trajectories of asteroids with irregular gravity. Meanwhile, indirect methods are also effective approaches for trajectory optimization with guaranteed solution accuracy and optimality [10]. Some practical techniques are employed to further improve the solution performance of the indirect methods for asteroid landing trajectory optimization, such as the homotopy approach [11] and the costate normalization [12]. Additionally, many studies have been focused on developing more precise and robust trajectory tracking techniques with closed-loop control methods, such as linear sliding-mode control [13], terminal sliding-mode control [14], convex optimization-based control [15], and model predictive control[16].

However, multiple issues need to be resolved for the trajectory tracking approach in asteroid landing missions [6]. First, both the direct and indirect methods suffer some drawbacks, including unpredictable iterative process and computational time, no guarantee of convergence to the optimal or even feasible solutions, and requirement of good initial guesses [17]. Second, high-fidelity gravity calculations are computationally intensive which extends the time-consumption of trajectory optimization by hundreds of times [5]. For this reason, researchers usually have to design a reference trajectory in advance at the expense of losing flight autonomy. Third, the separate design of the guidance and control systems can easily lead to actuator saturation [6]. Moreover, the control system also affects the finial performance of the guidance system which increases the complexity of system optimization. Accordingly, it is worth exploiting more effective landing approaches with higher autonomy and reliability and guaranteed flight optimality and real-time performance.

In addition to improving the performance of the indirect and direct methods, new control paradigms were also explored to enhance the autonomy and intelligence of the control system. In recent years, theoretical advancement in machine learning and improvement in artificial intelligence provide the possibility of real-time optimal control [18, 19]. According to the learning mode, machine learning algorithms are divided into three categories: supervised learning, unsupervised learning, and reinforcement learning [20]. An improved reinforcement learning algorithm was proposed for the trajectory optimization of interplanetary solar sail spacecraft where an evolutionary algorithm was used to optimize the neural networks [21]. A new actor-critic-identifier architecture was presented for approximate optimal control of certain nonlinear systems [22]. Recently, an integrated guidance and control algorithm was developed for planetary powered descent and landing by applying the principles of reinforcement learning theory, and the trained policy could directly predict the thrust commands based on the lander’s estimated state [19]. However, the current reinforcement learning algorithms are faced with two shortcomings: poor convergence property and hard reward function design. To address these issues, some researchers have begun to focus on combining the advantages of traditional trajectory optimization methods with the deep neural network (DNN) technology. Neural networks were used to approximate the value functions of the Hamilton$-$Jacobi$-$Bellman equations in [23, 24, 25]. In addition, neural networks were also trained with supervised signals to approximate the optimal control of some simple motion tasks in [26, 26]. For example, the DNNs were trained in a supervised manner in [12] using the optimal state-action pairs obtained via an indirect method, and the trained DNNs could drive the on-board landing systems with optimal control. Recurrent neural networks were also explored in a similar context in [6]. Recently, an interactive network learning algorithm was proposed to achieve real-time optimal control for solar sail orbit transfer missions with high terminal guidance accuracy [27].

Different from the above mentioned works, a novel real-time optimal control approach to asteroid landings is investigated in this study by developing a DNN-based gravity model and proposing a new DNN learning architecture for landing control. The contributions of this paper are threefold: First, DNNs are developed to approximate the irregular gravitational field of the asteroid, and the well-trained DNNs can calculate the high-fidelity gravity with a much higher speed comparing to the polyhedral method [28]. Second, incorporating the DNN-based gravity model in dynamics, an approximate indirect method is utilized to obtain optimal trajectories for landing, and a homotopy method is employed to improve the computational efficiency of the indirect method by supplying good initial guesses [10]. With the help of the gravity approximation technique and the homotopy method, this approximate indirect method can quickly generate the optimal trajectories for asteroid landing missions with irregular gravitational fields. Third, new DNNs are developed to approximate the optimal control actions obtained by the approximate indirect method, and the resulting DNN-based landing controller can autonomously drive the lander from a random initial condition to the target point in a real-time optimal manner.

This paper is organized as follows: In Sec. 2, the time-optimal asteroid landing problem is formulated as a two-point boundary value problem, and an indirect-method-based network learning architecture is presented to achieve real-time optimal control for this asteroid landing mission. In Sec. 3, a DNN is designed to approximate the gravitational field of the asteroid. Using this DNN-based gravity model in the dynamics, an approximate indirect method with a homotopy technique is presented based on the Pontryagin’s Minimum Principle in Sec. 4, and then a DNN-based landing controller is developed through learning the optimal actions obtained by the approximate indirect method. In Sec. 5, numerical simulations are given to evaluate the effectiveness of the proposed techniques and illustrate the performance of the designed DNN-based landing controller. The main work of this study is summarized in Sec. 6.

2 Problem Formulation and Solution Architecture

2.1 Problem Formulation

In this paper, we consider a time-optimal trajectory design problem for powered descent and landing on the surface of an asteroid. To achieve a soft pinpoint landing, a 3-dimensional (3D) time-optimal trajectory control problem is considered to find the thrust program that minimizes the following cost function

[TABLE]

subject to the equations of motion [9]

[TABLE]

and the boundary conditions

[TABLE]

where $\bm{r}=[x,y,z]^{T}$ and $\bm{v}=[v_{x},v_{y},v_{z}]^{T}$ are the lander’s position and velocity respectively with respect to an asteroid-fixed coordinate frame, and $\bm{\Omega}=[0,0,\omega]^{T}$ denotes the angular velocity vector of the asteroid. The notation $\bm{G}$ represents the force vector due to the asteroid’s gravitational field, $T_{max}$ is the maximum magnitude of the thrust, $\mu\in[0,1]$ is the thrust ratio, $\bm{\alpha}=[\alpha_{1},\alpha_{2},\alpha_{3}]$ is the direction vector of the thrust, and $m$ is the mass of the lander. The symbol $I_{\mathrm{sp}}$ denotes the thrust specific impulse, and $g_{0}=9.80665{\rm m/s^{2}}$ is the standard acceleration of gravity at the Earth sea level. The landing site is assumed to be preselected, and the landing velocity is assumed to zero in the asteroid-fixed frame for the purpose of soft landing. Different from the traditional methods, the initial state variables $[\bm{r}_{0},\bm{v}_{0},m_{0}]^{T}$ of the lander are assumed to be unknown beforehand. That is to say, the landing controller developed in this study is expected to have autonomous landing capability.

The considered time-optimal problem of asteroid landing is an optimal control problem that was traditionally solved by the indirect and direct methods [5]. Indirect and direct methods show excellent performance on constraints satisfaction, optimality, and convergence property on common trajectory optimization tasks; however, it is extremely difficult for on-board generation of asteroid optimal landing trajectories [6], due to the following reasons: First, solving the TPBVPs by indirect methods and the NLPs by direct methods needs good initial guesses [29]. Second, there is no guarantee of convergence to an optimal or even feasible solution [11]. Third, because extremely irregular and diverse asteroid shapes are involved, complicated gravitational field needs to be addressed. However, an accurate gravity calculation is computationally intensive, which leads to high time cost with hundreds of trajectory iterations [5]. In this study, DNNs accompanied by an improved indirect method are used to address these issues.

Considering the intensive computation time cased by the irregular gravitational filed, the time-optimal control problem is connected with a gravity-free time-optimal control problem with the same performance index and constraints, but the gravity $\bm{G}(r)$ is zero in the dynamical equations, that is

[TABLE]

Since the time-consuming gravity calculation is no longer needed, the relatively easier problem for gravity-free time-optimal control can be quickly solved and its solutions could serve as good initial guesses for the original time-optimal control problem with the help of a homotopy technique [10, 11].

2.2 Overall Network Learning Architecture

Taking advantage of high computational efficiency of the indirect methods and strong learning ability of the DNN technique, an indirect-method-based network learning architecture is proposed and depicted in Fig. 1 to achieve real-time optimal control for irregular asteroid landings. In this architecture, an indirect method is employed to solve the time-optimal landing problem and generate optimal state-action pairs. The control DNNs are trained off-line to approximate the functional relationship between the state and optimal actions, and can be used to generate optimal instructions for on-line landing control. Considering that the gravity calculation is usually computationally intensive [5], a DNN is developed to learn the irregular gravitational field based on samples generated from a polyhedral method. The gravity approximation can significantly speed up the solution process of the indirect methods in generating optimal trajectories. Additionally, a homotopy approach in terms of the gravity is used to connect the time-optimal control problem with the gravity-free time-optimal control problem. Since the gravity-free time-optimal control problem does not contain the gravity $\bm{G}(\bm{r})$, it could be solved quickly. Thereafter, the original time-optimal control problem could be solved with the good initial costates provided by the gravity-free time-optimal control problem. With the gravity approximation and homotopy techniques, the resulting approximate indirect method can efficiently generate optimal landing trajectories which contain optimal actions corresponding to different flight states. Through learning the samples extracted from the obtained optimal trajectories, the DNN-based controller could establish a functional relationship between the states and optimal actions and enable optimal action predictions. Since it no longer needs to solve the time-consuming trajectory optimization problem onboard, the DNN-based landing controller enjoys the ability of real-time optimal control for asteroid landing missions.

3 DNN Approximation of Gravity

To address the issue that traditional gravity models are usually computationally intensive, a DNN is developed in this section to approximate the irregular gravitational field of the asteroid. This DNN-based gravity model is aimed to improve the calculation speed of gravity prediction with guaranteed accuracy.

3.1 Gravitational Attraction Model

The accurate irregular gravity is traditionally computed by two methods: the spherical harmonic expansion method and the polyhedral method [28]. Although the spherical harmonic expansion method is easy to use, its drawback is that severe divergence appears close to the asteroid’s surface. In contrast, considering the exact potential of the polyhedron anywhere of an asteroid in space, the polyhedron potential method is a high-fidelity approximation of the gravitational field in the region close the surface on an asteroid [28]. Thus, the polyhedral method is employed to generate samples for DNN training in this study.

The exterior gravitational potential of a constant-density polyhedron can be derived analytically as follows [28]

[TABLE]

where $G_{u}$ is the gravitational constant and $\rho$ is the density of the asteroid. $\bm{r}_{e}$ is a vector from the field point to an arbitrary point on each edge, $\bm{E}_{e}$ is a dyad defined in terms of the face and edge normal vectors associated with each edge, $L_{e}$ is a logarithmic term expressing the potential of a one-dimenstional straight wire, $\bm{r}_{f}$ is a vector from the field point to an arbitrary point on each face, $\bm{F}_{f}$ is the outer product of the face normal vectors, $\omega_{f}$ is the solid angle subtended by a face when viewed from the field point.

As we can see that the formula in Eq. (6) involves the computation of sums over all the edges and faces of the polyhedron. Since an accurate shape model requires a large number of faces and edges, this approach is computationally intensive, particularly when the indirect method runs tens of thousands of times to generate optimal trajectories. To address this issue, a DNN is developed to learn the gravitational force filed. Then, the trained network can replace the polyhedral item in the dynamical model, and the computational speed of predicting the gravity can be significantly improved. In the subsequent section, we will present the detailed design of the DNN.

3.2 DNN Approximation

In this subsection, the detailed implementations of the gravity DNN are shown which includes the training data generation, data normalization, network architecture design, and optimization algorithm selection for network training. These algorithm implementations involve a series of scheme comparisons and parameter optimization, and due to the limited space, only the optimized results are provided here.

The first step is to generate the training data of the gravitational field for DNN learning. This study is focused on the asteroid landing missions, hence, the developed gravity DNN mainly approximates the position-based gravity around the target asteroid. Taking the near-Earth asteroid 443 Eros as an example, it is a typical irregular-shaped asteroid with estimated diameters of $34.4\times 11.2\times 11.2$ km. Considering that the lander flies near and outside the target asteroid, the polyhedral method is used to calculate the gravity force vector in a spherical shell space where the radius is greater than 3 km and less than 60 km, as displayed in Fig. 2. The position expressed in the spherical coordinate system is randomly chosen for the gravity calculation, and a total of 1,000,000 samples (position v.s. gravity force) are generated. According to the practical requirement of training, the samples are further randomly divided into two subsets with a ratio of 8:2 between the training dataset and test dataset. As their names imply, the training dataset is used to train the neural network, and the test dataset is used to evaluate the performance of network training and adjust the training parameters accordingly. Additionally, since the magnitude of the gravitational force is very small, we convert the input and output of the generated samples into the range of $[-1,1]$. Simulation experiments have proved that this normalization technique can improve the data identification ability of the DNN and thus improve the learning effect of gravity approximation. In addition, there is a data anti-normalization process to obtain the actual gravity prediction from the trained DNN.

A DNN is developed to learn the nonlinear functional relationship between the position and the gravitational force in this work. Since the gravity is dependent on the current position, deep feed-forward fully connected networks are employed. The design of the DNN architecture includes its activation function selection and size determination. The activation functions should be selected for both the output layers and the hidden layers. The activation functions of the output layers are determined by the value range of the sample output. Considering that the output of samples has no specific range limitation, thus, the gravity DNN adopts the Linear $[-\infty,\infty]$ activation function. Meanwhile, ReLU $[0,1]$ and Tanh $[-1,1]$ are two common activation functions for the hidden layers. In this approximation task, the ReLU activation function is found to be competent for the hidden layers and is employed in this study. The network size is another important user-defined parameter that affects the leaning effect, and both the depth and width should be carefully determined. Because the neuron numbers of the input and output layers are determined by the dimensions of the sample input and output, we only need to optimize the size of the hidden layers. As displayed in Table 1, 6 layers and 256 units are found to be competent for gravity approximation of the asteroid 443 Eros.

The parameter vector $\bm{\omega}$ of the DNN comprises the weights and the biases of the neurons. By adjusting $\bm{\omega}$ in a way that the mean square error (MSE) between the DNN’s predictions and the correct outputs from the training datasets is minimized, the DNN can accurately learn and approximate the gravitational field. The MSE is computed as follows

[TABLE]

where $Net_{\mathrm{G}}$ represents the gravity DNN and $\bm{G}_{i}$ is the correct gravity vector calculated by the polyhedral method for the ith sample. The symbol $N$ defines the number of samples that are randomly chosen from the training data to calculate the MSE in each training iteration. The algorithm Adam [30] is employed to minimize the MSE. As displayed in Table 1, the batch number $N$ is chosen as 1000 and the learning rate (LR) is set to be 0.0001.

Finally, the pseudo code of the training algorithm is summarized in Algorithm 1, which is implemented in Python environment on a typical desktop using Tensorflow [31], and the symbol $M$ denotes the maximum iteration number for the network training. To summarize, a specific DNN is developed in this subsection to approximate the gravitational field of our selected asteroid, and we will evaluate the approximation performance of the network $Net_{\mathrm{G}}$ through multiple experiments in the following subsection.

3.3 Accuracy Evaluation of DNN-based Gravity Model

Next, our attention turns to evaluate accuracy of the designed DNN for the irregular gravitational filed of the asteroid. The first experiment is to quantitatively analyze the approximation error of the gravity DNN where the relative errors projected onto the three axes and the synthetic relative error are defined as

[TABLE]

and

[TABLE]

respectively, where the subscript “${\mathrm{axis}}$” represents one of the three axes x, y, and z, and $\|\cdot\|_{2}$ denotes the Euclidean norm of a vector. For the 443 Eros example, Fig. 3 shows these four errors of 10,000 samples from the training dataset and test dataset using two box plots which display the distribution of four errors based on five number summary: minimum, first quartile, median, third quartile, and maximum. As can be seen from this figure, the relative errors along these three axes stay within the range of $\pm 0.75\%$, and most of their absolute values are less than $0.2\%$. Meanwhile, as a combination of the three axial errors, the synthetic relative error is relatively bigger with a median of about $0.3\%$ and a maximum of less than $1.6\%$. More than $3/4$ of the errors are less than $0.75\%$. Through this experiment, we can conclude that the developed gravity DNN can approximate the gravitational field of an asteroid with small prediction errors. However, considering that the landing is a long-lasting process, we should further evaluate the cumulative effect of the approximation error.

Since the polyhedral method is also an approximate gravity model, and there must be approximation errors due to the effects of uneven distribution of mass, density error, polyhedron calculation error, etc. In our second experiment, we compare the orbit cumulative difference between the assumed deviation of the polyhedral model and the approximation error of the DNN. Taking the 443 Eros as an example, Fig. 4 displays four trajectories of an unpowered lander that are driven by the ideal polyhedral model, the polyhedral model with $\pm 1\%$ deviation, and the gravity DNN. The left subgraph is an enlarged view of the end of these four trajectories in the right subgraph. As it is shown in the figure that the trajectories forced by the gravity with errors would gradually deviate from the nominal trajectory. Since there are numerical difference between the assumed deviation of the polyhedral model and the approximation error of the DNN, the deviations of the trajectories are also different. We can see from the figure that the trajectory deviation due to the DNN approximation error is much smaller than the ones due to $\pm 1\%$ shifting of the polyhedral model. Considering that $\pm 1\%$ approximation error is usually acceptable in traditional gravitational models, we can conclude that the gravity DNN has outstanding approximation accuracy.

The third experiment is to show the approximation results of the gravity DNN using graphs in the Cartesian coordinates. The gravitational force at one position of an asteroid is a 3D vector. Figure 5 (a) displays the vector-end curve of the gravity in a randomly chosen unpowered flight trajectory. As can be seen from the figure, the predictions of the trained DNN are almost the same as the results of the polyhedral method. Additionally, another advantage of the DNN-based gravity model is that the DNN is an explicit function $Net_{\mathrm{G}}(\bm{r}|\bm{\omega})$, and the partial derivatives of the gravity prediction (DNN output) with respect to the position (DNN input) can be obtained explicitly. And also, the subsequent indirect method would use these partial derivatives. Considering that the polyhedral method needs extra calculation of the partial derivatives, the DNN-based gravity model has another advantage in reducing the computational time. Subgraphs 5 (b)-(d) compare the vector-end curves of the partial derivatives $\partial\bm{G_{x}}/\partial\bm{r}$,$\partial\bm{G_{y}}/\partial\bm{r}$,$\partial\bm{G_{z}}/\partial\bm{r}$ under these two models. As can be seen from the subgraphs that the gravity DNN behaves good approximation ability to the partial derivatives of the gravity.

Through learning the samples of gravity versus position obtained by the polyhedral method, we obtain a new gravity model represented by DNNs. This DNN-based gravity model can quickly and accurately predict the gravity forces and the gradient of the gravitational field with respect to the position. The approximation accuracy of this DNN-based gravity model has been verified based on the above three experiments, and its real-time performance will be discussed in Section 5. The design DNN-based gravity model plays an important role in generating optimal landing trajectories in the subsequent section, and the trajectory generation efficiency will be significantly improved due to the excellent real-time performance of this model.

4 DNN-Based Real-Time Optimal Landing Control

Based on the DNN-based gravity model in Section 3, an approximate indirect method is developed in this section to generate optimal trajectories for asteroid landings with irregular gravitational fields. Then, DNNs are designed to approximate the optimal control actions with respect to the flight states through off-line learning samples extracted from the generated optimal trajectories. Since there is no need to solve optimal control problems onboard, the resulting DNN-based landing controller is capable of optimal control with great potential for real-time asteroid landing applications.

4.1 Approximate Indirect Method for Optimal Landing Trajectory

In this subsection, an approximate indirect method is developed for the time-optimal landing problem based on the above DNN-based gravity model. Applying the Pontryagin’s Minimum Principle, the Hamiltonian of the time-optimal problem is formulated as

[TABLE]

where $\bm{G}(r)$ is replaced with the trained network $Net_{\mathrm{G}}(\bm{r}|\bm{\omega})$. By introducing a positive costate normalization factor $\lambda_{0}$, the initial costate vector $\bm{\lambda}=[\lambda_{0},\bm{\lambda}_{r},\bm{\lambda_{v}},\lambda_{m}]$ can be restricted on the surface of a unit sphere and satisfies

[TABLE]

To minimize the Hamiltonian $H$ in Eq. (10), the optimal control law of thrust direction and magnitude can be obtained as

[TABLE]

where $\|\cdot\|$ represents the Euclidean norm and the superscript ∗ denotes the optimal control. It has been proved that the thrust of the lander is always maximized for the time-optimal landing missions [32].

The costate differential equations corresponding to the Euler-Lagrange equations are derived by $\dot{\bm{\lambda}}=-\partial{H}/\partial{\bm{x}}$ and shown below

[TABLE]

where ${\partial{Net_{\mathrm{G}}(\bm{r}|\bm{\omega})}}/{\partial{\bm{r}}}$ approximates the gradient vector of the gravitational force with respect to the corresponding position, and its value is calculated by a back-propagation algorithm in Tensorflow.

According to the static condition and the transversality condition, the optimal arrival time $t_{f}$ and the final mass costate are determined as

[TABLE]

Combining the above derivation, a two-point boundary value problem (TPBVP) is obtained wherein the corresponding shooting function is

[TABLE]

and the corresponding designed variable vector is $\bm{z}=[\lambda_{0},\bm{\lambda}_{r},\bm{\lambda}_{v},\lambda_{m},t_{f}]$. This TPBVP is a nonlinear multi-variable root-finding problem and MINPACK’s hybrd and hybrj routines (modified Powell method) [33] are found to be competent to solve this problem.

Although the DNN-based gravity model consumes much smaller time compared with the polyhedral model, high-fidelity integration of Eq. (2) is still time consuming. Moreover, the convergence of the solution method for the TPBVP depends on a good initial costate vector, and the computational time increases when a large number of iterations are required. Since substantial optimal trajectories are needed to train the landing controller, a fast approach to providing good initial costates for the approximate indirect method is necessary. To this end, a homotopy approach in terms of the gravity is used to connect the time-optimal landing problem with the gravity-free time-optimal landing problem, which can be solved much more quickly because the time-consuming calculation of gravitational force is not needed as shown in Eq.(5). Using a homotopy factor $\varepsilon$, we can rewrite the differential equation of the velocity as

[TABLE]

where the parameter $\varepsilon$ links the gravity-free time-optimal landing problem ($\varepsilon=0$) with the time-optimal landing problem ($\varepsilon=1$). After obtaining the solution to the gravity-free time-optimal landing problem, we will gradually increase the parameter $\varepsilon$ following a certain series and take the solution obtained in the current iteration as an initial guess for the next iteration. Finally, the solution to the original time-optimal landing problem can be obtained. Simulations show that because the gravity forces of asteroids are usually very small, the solutions to the gravity-free time-optimal landing problem can always provide good initial costates for the time-optimal landing problems, i.e., the homotopy takes one step in most scenarios.

4.2 Generation of Training Data

In traditional guidance and control approaches to asteroid landings, the lander’s initial states are assumed to be known and the time-consuming trajectory planning can be conducted in advance. However, this assumption is almost impossible due to the remoteness of the asteroids, rapidly variable gravitational fields, and uncertain force such as the solar wind and the aberration of light. In this study, we aim to develop a DNN-based intelligent landing controller which can steer the lander from an arbitrary initial state to the target point with high autonomy and reliability. To achieve this goal, we should get enough samples about the flight states versus optimal control actions and develop control DNNs to approximate the functional mapping relationship. In accordance with Bellman’s optimality principle, knowing the optimized vector $\bm{z}$ satisfying $\Phi(\bm{z})$ is equivalent to knowing the whole optimal trajectory (referring to Eq. (16)). Consequently, any state-action solution pair $[\bm{x},\bm{u}^{*}]$ along the optimal trajectory satisfies the first-order optimality condition, where $\bm{x}=[\bm{r},\bm{v}]$ denotes the flight state and $\bm{u}^{*}$ consists of the corresponding optimal control $[\bm{\alpha}^{*},\mu^{*}]$, time-to-go $t_{f}^{*}$, and expected mass consumption $\Delta m$. The trajectories generated by the approximate indirect method should cover the whole space that the landing process may happen. Taking the Eros as an example, we assume that the lander initiates its landing mission from a 3-dimensional cubic space as shown in Fig. 6 (a). The side length of this cube is 10 km, and the lander’s initial position is randomly generated in this cube. Meanwhile, the lander’s three initial velocity components are randomly chosen in the range of $[-0.01/1000,0.01/1000]$ km/s. A total of 10,000 optimal trajectories (Fig. 6 (b)) randomly starting from this cubic space are obtained via the approximate indirect method, and for every trajectory, 200 optimal state-action pairs are extracted with equal time intervals. Then, the resulting 2,000,000 optimal state-action pairs are divided into two subsets for the DNNs: the training dataset with 80% of the data and the test dataset with 20% of the data.

4.3 Network Optimization and Training

Five DNNs, $Net_{\alpha 1}$, $Net_{\alpha 2}$, $Net_{\alpha 3}$, $Net_{t_{f}}$, and $Net_{\Delta m}$, are designed to approximate the state-based optimal actions, which correspond to the three thrust directions ($\alpha_{1}^{*}$ , $\alpha_{2}^{*}$, and $\alpha_{3}^{*}$), time-to-go $t_{f}^{*}$, and mass consumption $\Delta m^{*}$. In order to avoid the approximation interference between control instructions, three independent DNNs, $Net_{\alpha 1}$, $Net_{\alpha 2}$, and $Net_{\alpha 3}$ with the same structure are used to approximate the thrust directions. Meanwhile, there is no need to learn the thrust ratio $\mu^{*}$ since it is known to be a constant of one according to Eq. (12). Table 2 exhibits the structures and parameters of the five DNNs, including the selection of activation functions, network size determination, loss definition, and user-defined parameters. With regard to the selection of activation functions, all five DNNs adopt ReLU $[0,\infty]$ for the hidden layers. Meanwhile, the activation functions for the output layers are determined based on their output value ranges. Specifically, Tanh $[-1,1]$ is used for the three DNNs for the thrust directions, while Softplus $[0,\infty]$ is adopted in the DNNs of $Net_{t_{f}}$ and $Net_{\Delta m}$ since the time-to-go $t_{f}^{*}$ and mass consumption $\Delta m^{*}$ are always greater than zero. Loss functions are used to quantify the difference between the predictions of the DNN and the target outputs from the training data. In this regard, the DNNs $Net_{t_{f}}$ and $Net_{\Delta m}$ adopt MSE as the loss function which is similar to the gravity DNN in Eq. (7). Meanwhile, the DNN $Net_{\alpha}$ uses the mean absolute error (MAE) as the loss function that is defined as

[TABLE]

where $\alpha$ represents one of the thrust directions $\alpha_{1}$ , $\alpha_{2}$, and $\alpha_{3}$, and the symbol $\bm{w}^{\alpha}$ represents the parameter vector of the networks. Similar to the gravity DNN developed in Section 3, the algorithm Adam is employed to minimize the MSE and MAE of the five networks wherein the learning rate (LR) is set to 0.0001 and the batch number $N$ is 1000.

The network size is an important user-defined parameter to guarantee the quality of network approximation. In this study, the network size is optimized based on its performance on the training and test datasets. Taking the size optimization of $Net_{\alpha_{1}}$ as an example, Fig. 7 displays the MAE comparison of DNNs with different layers and units. We can obtain the following conclusions from this figure. First, compared to the increment of the network width, the network depth is more conducive to improve the learning performance of the DNNs. Second, a too mall network size would lead to the under-fitting phenomenon of the DNN, which means that the DNN can not achieve enough accuracy on the training dataset. While a too big network size also causes the deterioration of the learning effect, and these will be an obvious gap between the training error and test error, which is termed as the over-fitting phenomenon. Based on the above analysis, the network size of $Net_{\alpha}$ is determined as 6 Layers/512 Units. By the similar means, the sizes of the DNNs $Net_{t_{f}}$ and $Net_{\Delta m}$ are set to be 6 Layer/256 Units.

The samples of flight state versus optimal control actions are obtained by the approximate indirect method with high computational efficiency. However, there are numerical differences between the elements in both the input and output vectors of the samples, because these elements have different physical meanings and magnitudes. In this study, the Scikit-learn package of machine learning in Python environment is used to linearly convert the values of sample elements into the range of [-1,1]. This data normalization can improve the data identification ability and thus boost the prediction accuracy of the DNNs. Additionally, there is a data anti-normalization module to calculate the actual action predictions from the normalized data. Finally, the pseudo code of the training algorithm for optimal action approximation is summarized as follows in Algorithm 2, which is actually implemented in Python Environment using Tensorflow.

4.4 DNN-Based Real-Time Optimal Landing

Through learning the training data generated by the proposed approximate indirect method, the DNNs are trained to approximate the state-based optimal actions. Figure 8 describes the real-time optimal control for asteroid landing driven by the trained DNNs. Considering the preprocessing of the training data, a normalization module is used for network inputs and an anti-normalization module is used for network predictions. Additionally, a multiplier of $1/\|\bm{\alpha}\|_{2}$ is introduced to normalize the thrust direction vector $\bm{\alpha}$. Different from the traditional indirect and direct methods, the well-trained DNNs can directly predict the optimal control actions according to the flight state and no longer needs to solve the original OCPs onboard. Consequently, the DNN-based landing controller is capable of real-time optimal control for asteroid landing missions. Multiple experiments will be carried out to evaluate the performance of this DNN-based landing controller in the next section.

5 Simulation Results

Numerical simulations are presented in this section to evaluate the effectiveness and performance of the proposed DNN-based landing control for asteroids with irregular gravitational fields. The evaluation of the DNN-based landing controller involves several aspects including the real-time performance of command generation, the accuracy of optimal action prediction, and the generalization capability with respect to the states. Furthermore, to quantitatively evaluate the proposed control process, the solution optimality, terminal landing accuracy, and control robustness under uncertainties are also investigated. Based on these considerations, four simulation experiments are conducted, and the performance analysis of the controller is discussed in detail. The proposed methodology has been implemented in C and Python on a desktop with Inter Core i7-7700 CPU @3.60GHz RAM 8.00GB and operating system of Windows 10. The near-Earth asteroid 443 Eros is used in the considered landing scenarios.

5.1 Real-Time Performance

One of the biggest advantages of the proposed method is its real-time performance for landing control using DNNs. Specifically, the gravity DNN developed in Section 3 is used to improve the speed of trajectory optimization and the control DNNs developed in Section 4 are developed to achieve real-time decision-makings for on-board landing control. Considering that the polyhedral method is implemented in C environment while our DNNs are developed in Tensorflow using Python environment, we have to evaluate the computational efficiency of these two simulation environments first.

Table 3 shows that the computational efficiency ratio between C and Python is 36:1, which is obtained through time-consumption analysis of a series of benchmark algorithms. The computational efficiency comparison of these two environments is complicated and may be influenced by the running algorithms, and the software and hardware of the computers. However, the running of algorithms in C is generally faster than in Python. The awareness of the above conclusion does not affect the results of the subsequent time comparisons, however, it will be beneficial for further studying the efficiency of the DNN-based controller.

First, the speed improvement in trajectory optimization benefited from the gravity DNN is investigated. A polyhedral model with 25,350 vertices and 49,152 faces is used for 433 Eros, and its calculation is time-consuming. As displayed in Table 3, the gravity calculation using the polyhedral method consumes about 41.3 ms in C. Meanwhile, the well-trained gravity DNN $Net_{\mathrm{G}}(\bm{r}|\bm{\omega})$ can predict the accurate value of gravity in 0.8 ms. As such, leveraging the DNN technique, the calculation speed of gravity is improved by 51 times (without considering the impact of simulation environment), and the time consumption of trajectory optimization is reduced from 295.2 min to 5.7 min. Moreover, the homotopy technique can further decrease the time consumption by improving the success rate of the approximate indirect method. Considering that the running time of algorithms in C is usually far less than that in Python, we anticipate that the implementing the DNN in C would bring much more speed increment to trajectory optimization.

In addition, the real-time performance is also observed in control command generation by the developed DNN-based landing controller, which is capable of directly predicting the optimal control instructions based on the flight state. Without the need to solve the OCPs onboard, the optimal control commands can be generated within 0.9 ms with a much higher updating frequency than any direct or indirect method. Additionally, the DNN-based trajectory controller does not suffer the convergence issues and can provide deterministic solutions to asteroid landing problems. The real-time performance and stable calculation process of the DNN-based landing controller show great potential for practical engineering applications [12].

5.2 Optimality Analysis

The network approximation accuracy to the optimal actions and the state-based generalization ability are two factors that determine the performance of the developed DNN-based landing controller and are needed to be investigated. Figure 9 displays the approximation results of the DNNs in learning the thrust directions $\bm{\alpha}$, time-to-go $t_{f}$, and mass consumption $\Delta m$. In this figure, all red solid lines represent the ideal outputs referring to the training dataset, while the blue dotted lines indicate the predictions of the trained DNNs. The left subgraph displays the approximation results to the three thrust direction components and the right one depicts the approximation results to the time-to-go $t_{f}$ and the mass consumption $\Delta m$. The right subgraph has two Y-axes, one of which indicates the time-to-go and the other one denotes the mass consumption. As we can see from the two subgraphs, the approximate optimal solutions obtained by the trained DNNs are very accurate. Even at the state with sharp changes in the instructions, the DNN results can still follow the ideal solutions very well. Additionally, the functional relationships of the time and mass consumption with respect to the states are relatively simple, thus, their approximations seem much better than the approximation of the thrust directions.

Then, we further quantitatively analyze the network approximation errors using the box plot. Figure 10 shows the statistical results of approximation errors of the five control DNNs. This box plot is drawn using 10,000 samples randomly chosen from the training dataset and the test dataset. The Errors are calculated according to Eq. (8). As we can see from the figure that the relative approximation errors of the three thrust directions stay within $\pm 0.25\%$ and the maxima of their absolute values are smaller than $1\%$. Due to the relatively simple functional relationships of the time and mass consumption with respect to the states, the approximation errors of the DNNs $Net_{t_{f}}$ and $Net_{\Delta m}$ are much smaller than other networks with the means within $\pm 0.15\%$ and the limitations within $\pm 0.6\%$. Based on these analysis, we can draw the conclusion that the trained control DNNs have high approximate accuracy to the optimal actions and are potential for on-board solution generations for asteroid landing missions.

In addition, the asteroid landing flight is a long-lasting control process, and small network fitting error may propagate and cause large flight error. However, large trajectory deviation can be reduced through continuous command corrections. Considering the existing approximation errors of DNNs, two flight cases driven by the DNN-based landing controller are carried out to test the controller’s performance on reducing the propagated errors in a closed-loop control manner. The lander’s initial state is randomly selected from the initial space illustrated in Fig. 6. The whole flight is completely driven by the DNN-based landing controller, which independently determines the control instructions according to the actual flight state. The simulation is terminated based on the predicted landing time from the network $Net_{t_{f}}$. The resulting flight profiles are compared to the ideal flight solutions obtained by the indirect method and shown in Figs. 11 and 12. In these two figures, the left subgraphs exhibit the trajectories from these two approaches and the right subgraphs compare the corresponding profiles of the three thrust directions $\alpha_{1}$, $\alpha_{2}$, and $\alpha_{3}$. As it can be seen from the figures that the developed DNN-based landing controller is capable of independently driving the lander to the target state in an optimal manner, and the driven trajectories coincide with the ideal trajectories designed off-line by the indirect method. The above two experiments further verify the optimality of the solution generated by the proposed DNN-based controller for asteroid landing missions.

5.3 Terminal Guidance Error

In addition to the real-time calculation process and the guaranteed solution optimality, the accuracy at the terminal landing moment is also critical to the success of the entire mission. In this subsection, we will discuss the terminal guidance accuracy of the asteroid landing flight driven by the proposed DNN-based landing controller.

Figure 13 displays 300 DNN-driven trajectories starting from random initial states. As it is shown in the figure that although the start conditions of the lander are randomly initialized, the DNN-based landing controller can autonomously decide the optimal control instructions according to the flight states and steer the lander to the target site accurately. Considering that the flight state is continuously changing and there must be states that do not exist in the training dataset, we can conclude that this DNN-based controller has excellent generalization ability for state-based decision-makings. Figure 14 illustrates the terminal positions and velocities of these 300 trajectories. We can see from the figure that the terminal guidance position errors are within the range of $\pm 0.01$km, and the terminal velocity errors are within the range of $\pm 4\times 10^{-5}$ km/s. Hence, we can conclude that the DNN-based landing controller is capable of steering the lander to the target point with high accuracy. The terminal guidance accuracy can be further improved using some auxiliary techniques, for example, a multi-scale network cooperation strategy was proposed in [27] to improve the terminal flight accuracy of solar sail orbit transfers. In addition, a compound control scheme was developed in [34] to improve the terminal landing accuracy for the Moon landing mission.

5.4 Robustness Analysis

In consideration of the inevitable dynamical uncertainties caused by the approximation error of gravitational force, the solar wind, and the aberration of light, we test the robustness of the DNN-based landing controller in this subsection. For preliminary analysis, we assume that the gravity calculation errors and all other force interferences are considered together as the deviation of the gravitational field. Figure 15 shown the trajectories under the gravity with $-50\%$ deviation. The red solid trajectory is the ideal solution without gravity deviations. Under the effect of gravity deviation, the black dotted trajectory is obtained by the open-loop indirect method and the blue dashed trajectory is the solution by the DNN-based controller. As it can be seen from the figure that the flight by the open-loop indirect method has big terminal state error. Traditionally, this error could be reduced by a closed-loop trajectory tracking system. Designing a closed-loop tracking scheme is not an easy task due to the difficulty of generating a reference trajectory in advance, the risk of trajectory tracking failure, and the lack of flight autonomy. In this study, the proposed DNN-based landing controller can adjust the control instruction on-line based on actual flight state in every guidance cycle. Consequently, the proposed DNN-based controller can significantly reduce the flight errors through rolling command corrections and achieve high model robustness and terminal landing accuracy, even when large deviations occur in the gravity field.

Table 4 presents the statistical results of terminal flight errors under different gravity deviations. As we can see from the table that the open-loop control scheme based on the indirect method has relatively large flight errors and behaves poor robustness when obvious uncertainties are incorporated in the dynamics. In contrast, benefiting from its rolling trajectory correction capability, the DNN-based landing controller is more robust and achieves more accurate terminal flight accuracy than the traditional indirect methods.

6 Conclusion

In this study, deep neural network technology is introduced to improve the real-time performance and flight autonomy for asteroid landing control. This work processes the following three contributions. First, a new DNN-based gravity model is developed to replace the time-consuming homogeneous polyhedral method and significantly improve the computational efficiency of gravity determination. Then, an improved approximate indirect method is proposed to solve the time-optimal landing problems where the gravity term in the dynamical system is replaced with the designed DNN-based gravity model and a homotopy method is employed to further reduce the computational cost by connecting the original time-optimal problem to a gravity-free time-optimal problem. Furthermore, five DNNs are developed to learn the state-based optimal actions obtained from the approximate indirect method, and a DNN-based landing controller composed of these DNNs is eventually designed to successfully achieve the real-time optimal control for asteroid landing missions. Simulation results of time-optimal landing problems of landing on 443 Eros are given to substantiate the improvement of real-time performance and illustrate the effectiveness of the developed DNN-based controller in solution optimality, terminal flight accuracy, and robustness under uncertainties.

Acknowledgments

This work is supported by the National Natural Science Foundation of China (Grants 11872223, 11672146 and 11432001).