Abstract
Guidance systems are important to autonomous rendezvous with uncooperative targets such as an active debris removal (ADR) mission. A novel guidance frame is established in rotating line-of-sight (LOS) coordinates, which resolves the coupling effect between pitch and yaw planes in a general 3D scenario. The guidance law is named augmented proportional navigation (APN) by applying nonlinear control along LOS and classical proportional navigation normal to LOS. As saving time is a critical factor in space rescue and on-orbit service, the finite time convergence APN (FTCAPN) is further proposed which proves to possess convergence and high robustness. This paper builds on previous efforts in polynomial chaos expansion (PCE) to develop an efficient analysis technique for guidance algorithms. A large scope of uncertainty sources are considered to make state evaluation trustworthy and provide precise prediction of trajectory bias. The simulation results show that the accuracy of the proposed method is compatible with Monte Carlo simulation which requires extensive computational effort.
1. Introduction
A critical challenge for missions such as on-orbit servicing and active debris removal is to perform autonomous rendezvous with uncooperative target which provides limited measurement information and even evades with unexpected motion [1]. Generally, the Clohessy-Wiltshire (CW) equations [2] or Tschauner-Hempel (TH) equations [3] are used as relative dynamic equations for cooperative rendezvous guidance, with coordinate origin set on the target. However, for an autonomous rendezvous mission such as active debris removal (ADR), the target location and center of mass are unknown. Further, there is collision risk in an impulsive trajectory using CW equations. The measurement information is confined to be relative distance and line-of-sight (LOS) angles, provided by a laser range finder and visible sensor [4]. Besides, the unmodeled perturbation factors and unknown target motion require that the guidance law should be robust with uncertainty input, so it is necessary to pursue a novel and robust frame for guidance. And the high-precision relative motion model contributes to improving efficiency of the navigation system.
Augmented proportional navigation (APN) is initially proposed in our peer-review paper [5] to construct a three-dimensional dynamic equation set in a rotating coordinate system. The rotating coordinate system is established based on differential geometric theory. The expression is more clear compared with the study of Meng and Chen [6] in which the motion model is established in inertial coordinates based on differential geometric theory: where indicates unit vector along LOS and is the unit vector of LOS velocity; denoting , is a normal vector of LOS, and then, a rotating coordinate named “LOS rotation coordinate” is established with three orthonormal axes, , , and . The plane composed of and is the instantaneous rotation plane of LOS (IRPL); is the LOS rotation rate, and . is the IRPL rotation rate or the angular velocity of a fly-around plane, and .
The rotating LOS coordinate system and relative motion are shown in Figure 1.
(a) Coordinate systems
(b) 3D rendezvous geometry in inertial coordinates
Figure 1 illustrates the relationship between three coordinate systems where subscript represents inertial frame, subscript represents ordinary LOS coordinate, and denotes the rotating LOS coordinate. and are the elevation angle and azimuth angle, respectively.
APN implements feedback control along LOS, and classical proportional navigation is applied in transverse direction to impel the LOS rate to zero which guides the chaser to the rendezvous trajectory. With feedback control along LOS, the dynamic model is decoupled and linearized to obtain analytical solutions of the LOS rate and then a covariance analysis (CA) technique can be applied for performance assessment of the APN guidance system [5].
In this paper, nonlinear APN is further developed by introducing finite time convergence sliding mode control along LOS, which is named finite time convergence augmented proportional navigation (FTCAPN) for the first time. Its aim is ensuring that the state variable converges to the sliding surface in finite time; as a result, the chaser finishes the midrange rendezvous mission in finite time. FTCAPN is a highly nonlinear guidance system, and consequently, the linear CA technique is no longer proper for performance assessment. For a nonlinear system, the output is likely to be non-Gaussian even with Gaussian input. So the mean and covariance are not sufficient to represent the probability density function (PDF) of a state variable, which is ignored in current references [7–10]. So, higher moments are essential to characterize state variable distribution, i.e., skewness and kurtosis. In references [11, 12], non-Gaussian nature is studied accounting for space debris long period propagation. They provide incentive for our study, whereas there is no close-loop control in the above state propagation.
One important part of performance assessment is to determine the effect of input uncertainty on output result, which is also called uncertainty propagation or sensitivity analysis. A proper propagation method could make a precise prediction of trajectory bias, through which it is possible to evaluate performance of the guidance system and provide safety analysis result for trajectory programming and following operations.
Many techniques [7–12] have been proposed to tackle the above issue, e.g., Monte Carlo, covariance analysis, and some other quantified evaluation methods. However, some challenges occurred when applying these methods to implement performance evaluation of the guidance system. (1)Cost: cost of Monte Carlo (MC) is high with a number of propagations although it is commonly used in practice. An MC approach is based on random sampling. MC shows slow convergence, and a large number of executions are needed. (2) Linearization: the most popular nonsampling method is the covariance analysis method. It is a method for directly determining the statistical properties of solutions of linear systems with random inputs. Covariance analysis describing function technique (CADET) is developed to solve nonlinear systems. A nonlinear system is approximated by linear function using quasilinearization methods. The essence undying quasilinear analysis is the substitution of one or more approximate linear gains (describing functions) for each system nonlinearity. The analytic form of linearization is determined by the nonlinearity, number of state variables, and error criterion. Besides, the validity of CADET is established on the assumption that the state variable is Gaussian. Non-Gaussian nature of state variables after the nonlinearity leads to inaccurate CADET results. (3) Non-Gaussian: covariance realism is necessary while not a sufficient realism in describing uncertainty [5]. For non-Gaussian distribution, higher-order cumulants are needed to better represent the errors beyond the mean and covariance, i.e., skewness and kurtosis. Through data analysis, the non-Gaussian feature of state variables is significant because of a nonlinear dynamic model, large input uncertainty, etc.
This paper consists of 5 sections. In Section 2, the augmented proportional navigation (APN) is introduced. Then, its non-Gaussian nature is investigated and presented in Section 3. Section 4 depicts application of polynomial chaos expansion (PCE) in performance assessment, and a study of adaptive PCE is conducted to further reduce the computational cost. Section 5 gives the conclusions of this study after analyzing numerical simulation results.
2. Finite Time Convergence Augmented Proportional Navigation
Different navigation and control strategies should be adopted corresponding to specific measurement information in different rendezvous phases [13]. Our work focuses on guidance and control solely on midrange rendezvous which lasts from several kilometers to tens of meters [14]. In this paper, a three-dimensional (3D) dynamic equation set is constructed in a rotating coordinate system. In FTCAPN, finite time convergence sliding mode control is implemented along LOS, and classical proportional navigation is applied in transverse direction to impel the LOS rate to zero which guides the chaser to the rendezvous course. FTCAPN is a nonlinear guidance system, and it is introduced as follows.
2.1. Scheme of Finite Time Convergence Augmented Proportional Navigation
The 3D relative dynamic equations are established in a rotating LOS coordinate system which resolves the coupling effect between pitch and yaw planes in a general 3D scenario. Introduction of IRPL (instantaneous rotation plane of LOS) simplifies control of chaser; besides, it improves guidance efficiency with smaller miss distance and less fuel consumption. The applications of a rotating LOS coordinate system are also investigated in 3D missile intercept scenarios in references [15, 16].
A set of dynamic equations is derived based on equation (1): where and represent accelerations of the target and chaser, respectively. is the navigation constant, and it is usually set as a value between 3 and 5. “, , and ” denote projections along three axes of the rotating coordinate, respectively.
If the target is nonmaneuvering, the relative dynamic equation (1) could be written as
As mentioned above, the sliding mode control is introduced along the direction of LOS. The tracking error is defined as
Rendezvous occurs when the relative distance and relative velocity are equal to zero; thus, a sliding surface function is defined as , where is a constant value. Then, the control of the chaser is defined as follows: where is relative velocity, is disturbance, and . , , and are the suitable parameters, and is a sign function. The Lyapunov function candidate is defined:
Differentiating equation (7) with respect to time using equations (5) and (6), can be obtained:
So,
The stability of the control system along LOS can be guaranteed. Further, the convergence time can be calculated according to the following lemma.
Lemma. Suppose that there is a continuously differentiable function defined on , and there are real numbers and , besides is positive-definite on and also . Then, there is an area , and any from can reach in finite time. Moreover, if is the time required to reach , then where subscript 0 represents the initial values of the state variable. From equation (8), we can derive that Assuming and , then , , and equation (7) are substituted into equation (10):
To guarantee that the acceleration is within the limitation of maximum thrust, a saturation function is introduced:
In order to avoid chatter at the end of proximity, the sliding surface is modified into the following form: where is a small constant.
From the above derivation, we can see that the system will converge and the upper bound of converge time can be derived from initial conditions and control parameter settings. The midrange rendezvous guidance system, FTCAPN, is a highly nonlinear system. Its effectiveness is demonstrated by simulation results as follows. Its performance will be further analyzed in Section 4 to reveal the extent to which an output parameter’s value depends on an input parameter’s value, which is also named sensitivity analyses. The output, i.e., trajectory bias between the target and chaser, will be calculated directly through performance assessment or sensitivity analyses. Then, the bias can be used for evaluating guidance precision and proximity safety analysis.
2.2. Simulation of FTCAPN
FTCAPN proposed in our paper can be used for autonomous rendezvous with space debris to avoid risk to humans. As we know, debris remediation in LEO is more urgent compared with GEO and MEO, so we choose a target (NORAD catalog number: 39357) at 700 km altitude, and its orbital elements are listed in Table 1.
This paper focuses on the midrange proximity. The initial state of the chaser is provided by a relative state to the target as shown in Table 2.
According to Tables 1 and 2, the initial relative distance, velocity, and LOS rate can be calculated km, m/s, and rad/s, respectively; then, simulation results are shown below. To be noted, only the midrange rendezvous phase is considered in this paper. When the relative distance is 50 m, APN control stops considering proximity operation safety, visibility conditions, accuracy requirement, etc. Then, close-range navigation and control laws are adopted corresponding to sensors like LiDAR with more precise measurements. The control coefficients are set as , , , and .
Figure 2 shows the trajectory between the target and the chaser. This paper tries to achieve finite time rendezvous in an urgent situation, which is completely contrast to a conventional rendezvous mission in which fuel saving is given top priority. Under current scenario settings, the LOS rate keeps stable after 60 seconds and the chaser is able to achieve rendezvous with target in 85 seconds as shown in Figure 3. According to equation (12), the rendezvous time is calculated as , which demonstrates the convergence feature of FTCAPN.
When relative distance reaches 50 meters, simulation stops and the relative velocity is 5 m/s at the instant of stop (Figure 4). Total fuel consumption () is about 98.5 m/s. Initial relative velocity is quite high and decreases quickly (Figure 5). Consequently, the rendezvous mission can be accomplished faster, which can also be adjusted by different parameter settings.
Zero-Effort-Miss (ZEM) is defined as the closest distance between the chaser and target after control stops [17]. In this paper, it indicates the minimum relative distance after midrange control stops at the distance of 50 m. It is commonly used as the index to evaluate the guidance efficiency: where and are velocity components along LOS and perpendicular to LOS, respectively.
is transferred into a function of relative state. Based on the above equation, ZEM is calculated with a value of 0.003767 m. We can see that when the control on the chaser stops, the chaser is already on a rendezvous course with in Figure 3, so the chaser keeps approaching close to the target with constant bearing.
The deviation of equation (15) is calculated:
So, the uncertainty of ZEM could also be deduced by uncertainty in the LOS rate ().
3. Non-Gaussian Nature Analysis
Section 2 demonstrates the effectiveness of FTCAPN, and the robustness of FTCAPN with input uncertainty will be further investigated. More importantly, this section will explore the non-Gaussian feature of LOS rate uncertainty through MC analysis as the LOS rate determines the efficiency of the control system [17, 18]. Skewness and kurtosis are calculated to illustrate its non-Gaussian distribution based on the following equations:
Gaussian distribution is generally a good approximation for initial state uncertainty and measurement uncertainty. Propagating the uncertainty of the state variable through a nonlinear dynamic model transforms the state variable into a non-Gaussian distribution. Non-Gaussian behavior is determined by the following two terms: (1) nonlinearity of function and (2) magnitude of initial distribution uncertainty (standard deviation). Quantification of uncertainty sources is shown in Table 3.
3.1. Robustness Analysis of FTCAPN
Figures 6 and 7 show the time curves of the LOS rate mean and standard variance, which prove to be asymptotically stable. We can come to the conclusion that the FTCAPN still shows good performance with input uncertainty and measurement uncertainty.
Although the skewness value fluctuates around 0 in Figure 8, the value of initial kurtosis in Figure 9 is far from 3 which indicates that the state variable distribution is not Gaussian. However, the LOS rate distribution approaches Gaussian distribution gradually as the chaser approaches the target.
3.2. Non-Gaussian Nature Analysis of FTCAPN
This part will investigate the non-Gaussian nature of FTCAPN by the MC technique with larger magnitude of initial uncertainty. It takes about 1405 seconds to run 1000 times of MC simulation which is implemented single thread in MATLAB environment on a laptop with Intel Core i7 2.4 GHz and 8.0 GB of RAM. Propagations of the LOS rate are shown below. (1)Larger uncertainty in initial position error (Figures 10 and 11)
Figures 10 and 11 illustrate skewness and kurtosis of the LOS rate based on the dynamical model in Section 2. Compared with the simulation result in Section 3.1, a non-Gaussian feature is quite distinct at the early phase of proximity. When coming closer to a nonmaneuvering target, the values of skewness and kurtosis vary gradually and indicate that the distribution of the LOS rate turns to be Gaussian after reaching the sliding surface. (2)Larger uncertainty in initial position and velocity error (Figures 12 and 13)
With position and velocity initial uncertainty enlarged by 10 times, the non-Gaussian feature of variable distribution is significant as shown in Figures 12 and 13. Whereas at the end of the proximity phase, the state variable remains stable near the sliding surface and the variable of the LOS rate drifts towards Gaussian distribution. So, there will be large error if the LOS rate uncertainty is propagated solely by the mean and covariance.
Finally, we come to the conclusion that it is not accurate enough to study uncertainty propagation through the covariance analysis method which is based on propagation of the mean and covariance. We need to develop a more precise method.
4. Application of Polynomial Chaos Expansion
Robustness of APN is already demonstrated with input error in the initial state and the LOS rate. The efficiency with input of other types of uncertainty is yet to be discovered. In near-distance proximity, uncertainty sources include initial position error, initial velocity error, LOS measurement error, target maneuver acceleration and uncertainty, distance measurement error, commanded acceleration saturation, control delay, and nonlinearity of the control system. Consequently, the APN guidance system is a nonlinear and time-varying system with noisy input and other disturbances.
As mentioned in Introduction, existing performance assessment methods are constrained when dealing with a highly nonlinear system. In this section, PCE is introduced for performance assessment of the proposed guidance system, i.e., FTCAPN. In space community, PCE is already applied in orbit propagation of space debris by Jones et al. [11], entry descent and landing of a Mars vehicle by Prabhakar et al. [20], and other missions [21]. However, it is not yet used in performance assessment of a nonlinear guidance system.
4.1. General Introduction of Polynomial Chaos Expansion (PCE)
The term polynomial chaos (PC) is invented by Norbert Wiener and Hermite who used polynomials as an orthogonal basis to study the Gaussian stochastic processes [22, 23]. In polynomial chaos, the type of orthogonal polynomial is determined by the probability distribution of the random inputs.
An efficient numerical method is the stochastic collocation approach, which is a deterministic sampling method. The full tensor product of one-dimensional nodes is extensively used, and postprocessing is conducted to obtain the statistical feature from the solution ensemble. However, the number of total nodes surges exponentially along with the increasing quantities of uncertainty sources, which is known as curse of dimensionality. Sparse PCE, a subset of full tensor PCE, is gaining increasing interest for computational saving with higher uncertainty dimensions.
The approximation of function, , by PCE is the optimal linear approximation in a finite dimensional space spanned by certain multidimensional polynomials. PCE can be regarded as multidimensional polynomial surrogates to the original problem as shown in the following equation:
The basis functions are orthonormal polynomials. Distribution of determines the pattern, e.g., is tensor products of Hermite polynomials if is Gauss distributed. Higher accuracy of approximation is possible with larger . The total number of polynomial terms is derived: where indicates the order and is the number of dimensions:
Estimation of coefficients is the result of expectation operation, which is actually integral operation. There are two methods to calculate the integration, i.e., tensor product and sparse grid.
Once the coefficients are calculated, the statistical characteristics of can be derived by sampling or directly calculated by the coefficients. The mean value is
From the above equation, we can see that the mean value of equals the coefficient with the basis function .
The covariance matrix is
4.2. Full Tensor Polynomial Chaos
As mentioned in Section 4.1, the major method for the coefficient is to solve the integration by the collocation on certain quadrature nodes. Thus, equation (20) is solved by evaluating at given quadrature nodes using nonintrusive polynomial chaos expansion: where is the weighting coefficient corresponding to each collocation node.
The extension from univariate to multivariate is easily achieved by a tensor product of the univariate node set. Thus, the total number of evaluations is . If the number of nodes () in each dimension is equal, then , and the maximum order of the polynomial, i.e., , equals .
Figures 14 and 15 illustrate the propagation of mean value and standard deviation of the LOS rate, respectively, which are calculated through tensor PCE (TPCE) and Monte Carlo (MC) technique. In the simulation, the MC technique runs 1000 times. The level of TPCE is set as two which means the number of node in each dimension is 2.
The result curve of PCE matches with MC, which means they can achieve the same precision in performance analysis of the terminal guidance system. However, the computation time of MC and TPCE is 1400 and 182 seconds, respectively. So, the simulation results substantially demonstrate that TPCE outperform the MC method in computation efficiency.
4.3. Sparse Polynomial Chaos
Sparse polynomial chaos is proposed to reduce the number of nodes by the Smolyak algorithm. The reduction is achieved by eliminating high-order interactions. In the polynomials, where is the order of polynomials and is the level of sparse PCE (point in one dimension); however, regarding full tensor PCE in Section 4.2,
So, the SPCE achieves high computation efficiency with fewer nodes and expense of accuracy. To be noted, computational burden of SPCE is larger than TPCE in lower dimensions (generally, ).
In this case, the SPCE of level 2 is used to calculate the uncertainty propagation of the LOS rate. The model is evaluated for merely 15 times, and the computation time is reduced to 22 seconds. Compared with TPC, SPC achieves a high-precision result with less computation time as shown in Figures 16 and 17.
4.4. Adaptive Polynomial Chaos Expansions
Without mitigation techniques, the PCE method suffers from the curse of dimensionality. Computation time increases exponentially with respect to input dimensions , (). Further improvement is required to achieve the ZEM more efficiently in case of large . The adaptive PCE method would be a solution to solve the curse of dimensionality [24]. The adaptive polynomial chaos includes two parts: (1)Basis adaptive PCE. The maximum degree of polynomials is increased gradually to find the maximum degree which meets the requirement of accuracy(2)Sparse PCE. Only a subset of the most relevant polynomials is retained, while other coefficients are set to 0
The size of uncertainty dimensionality is assumed to be 9 which includes uncertainty in initial position , initial velocity (), LOS rate measurement (), distance measurement (), and target maneuver (), which are listed in Table 4. Target maneuver is projected into three directions along the axis of LOS coordinates. The simulation step is chosen to be 0.5 seconds.
Compared with former simulations, two more uncertainty sources are included, i.e., distance measurement error and target maneuver. Consequently, equation (16) is not proper to calculate the ZEM distribution. So, the programming codes are modified with the output set as ZEM directly. In the simulation, the range of polynomial degree is [2, 6] and truncation coefficient is set to be 0.5.
From Table 5, we can see that PCE is able to approximate the MC result with a high precision, whereas its time consumption is only six percent of MC. Besides, the distribution of ZEM is approximately Gaussian. It can also be verified by Section 2.2, because at the end of proximity, the relative motion is stable. However, during the midcourse of proximity, the non-Gaussian feature is significant and covariance realism is not sufficient to depict the distribution of LOS rate uncertainty.
5. Conclusion
This paper proposed the finite time convergence augmented proportional navigation (FTCAPN) which can be used in time-critical autonomous rendezvous missions such as active debris removal or space rescue. The time to converge is verified and worked out through the Lyapunov function and formula derivation. Robustness of FTCAPN is proved considering input uncertainty in the initial state and LOS rate measurement.
Further, this paper investigated non-Gaussian nature of the proposed guidance law (FTCAPN) which is not tackled in current literatures. Then, polynomial chaos expansion (PCE) is applied to performance assessment of FTCAPN. At last, adaptive PCE is proposed to reduce the computation burden when there is a large scope of random input sources. Simulation results demonstrate that the accuracy of the proposed method is compatible with Monte Carlo, whereas the latter requires extensive computational effort.
Data Availability
The data used to support the findings of this study are included within the article.
Conflicts of Interest
The authors declare that they have no conflicts of interest.
Acknowledgments
This work was supported by the Major Program of National Natural Science Foundation of China under Grant Numbers 61690210 and 61690213.