Robust Trajectory Planning Of Gliding-Guided Projectiles With Weak Maneuverability

1. Introduction

Artillery, as a weapon with a long and influential history in warfare, is hailed as the god of war. In the modern combat environment, with the advancement of science and technology and the innovation of combat concepts, the capability for beyond-visual-range precision strikes has increasingly become one of the key indicators of military strength. Artillery systems have gradually evolved from purely area-damage weapons to multifunctional strike platforms. Gliding-guided projectiles, as a new type of artillery-launched weapon, combine the advantages of both tactical missiles and conventional projectiles. They not only offer low cost per shot, rapid response, and flexible use but also possess extended range and high strike accuracy, making them one of the focal directions in the current development of weaponry and munitions[1].

Among the many key technologies in the research of gliding-guided projectiles, trajectory planning has always been a core and hot issue. Limited by the stringent requirements of the gun-launched platform and the harsh dynamic environment during launch, gliding projectiles are small in size with limited control fins. The navigation and control devices that can be loaded, as well as the overall structural layout of the projectile, have special requirements. Moreover, the projectile usually undergoes unpowered flight during the controlled trajectory phase. Consequently, gliding projectiles possess a low lift-to-drag ratio and limited maneuverability and anti-interference capabilities, under which the rationality of the trajectory design becomes particularly crucial[2]. Addressing this issue, Shi conducted an in-depth study on the trajectory characteristics of gliding projectiles using the energy analysis method, forming a relatively complete theoretical system[3]. Building on this, Yi planned a range-extension trajectory based on the maximum lift-to-drag ratio as the performance indicator[4]. Chen addressed the matching issues of dynamic and static parameters during the flight of gliding projectiles and proposed a full-trajectory planning method aimed at minimizing control effort as the performance indicator[5]. Following this approach, Xu introduced the concept of a composite efficiency factor, taking into account both the gliding efficiency and maneuverability of the projectile, and used it as a performance indicator to conduct trajectory design works[6]. The evolution of modern warfare has given rise to an increasingly diverse array of combat mission requirements, prompting a surge in research directions related to trajectory planning. Yan proposed a three-dimensional trajectory planning method for cruise missiles' low-altitude penetration issues based on ant colony optimization and Bezier curve optimization, addressing the problems of excessive turning points and large turn angles common in traditional methods[7]. He employed swarm intelligence algorithms such as the grey wolf optimizer, firefly algorithm, and particle swarm optimization for maximizing the maneuver characteristics of cruise missiles, and compared their respective advantages and disadvantages[8]. Liu introduced a multi-stage rapid trajectory planning algorithm based on sequential convex optimization for the online planning of multi-stage trajectories[9]. Beyond planning for individual projectiles, Liu and Pinon separately researched collaborative multiple-projectile and 'autonomous actors-artillery collaboration' trajectory planning[10,11]. With the advancement of science and technology, cutting-edge technologies like neural networks[12] and artificial intelligence[13] are being progressively incorporated into trajectory planning processes, demonstrating strong developmental momentum.

The aforementioned literatures extensively explore the trajectory planning problems from various perspectives, focusing primarily on deterministic optimal control issues, however, they do not take into account the impact of uncertainties, leading to overly idealized trajectory designs. During the actual flight of a projectile, uncertainties exist in the dynamics model, aerodynamic and environmental parameters, and states at the control-start point. Planning variables, performance indicators, and constraints are closely related to these parameters, and results obtained from deterministic models may lead to exceeding variable boundaries, deviations in performance indicators, and breaches in constraint conditions. This issue can be mitigated through robust planning, which fundamentally aims to minimize the impact of uncertainties on the trajectory without eliminating or reducing the sources of uncertainty[14]. In the process of robust planning, it is critical to understand the propagation of stochastic disturbances along the trajectory under current control commands, necessitating research into quantitative methods for dynamic system uncertainties. The Linear Analysis of Covariance (LAC) method, which requires only a single integration of the error propagation equation to obtain statistical results akin to those from Monte Carlo Simulation (MCS), has been widely applied due to its efficiency in robust planning issues[15,16]. However, LAC necessitates linearized models, which can lead to significant approximation errors when dealing with highly nonlinear problems. Additionally, robust planning based on LAC is sensitive to initial conditions since it reshapes around a reference trajectory, often resulting in local optima that prematurely terminate the algorithm[17]. Given the limitations of LAC, scholars studying robust planning problems for Unmanned Aerial Vehicles (UAVs) have adopted methods based on nonlinear models such as Polynomial Chaos Expansion (PCE)[18,19,20,21,22] and Unscented Transform[23,24,25,26,27], enhancing the precision and efficiency of problem-solving.

The aforementioned studies pertain to open-loop robust planning methods, overlooking the coupling between the planned trajectory and the guidance control system. In actual conditions, deviations caused by various uncertainties are typically mitigated by the guidance control system. Consequently, the compatibility between the planned trajectory and the guidance control system significantly impacts the effectiveness of mission execution. For UAVs with ample maneuvering capability, the performance loss associated with open-loop planning does not pose a substantial issue. However, for the gliding projectile with limited rudder control capability and unpowered flight, it is critical to minimize losses and fully exploit maneuverability. Therefore, it is necessary to pursue research on robust trajectory planning methods that incorporate closed-loop feedback from the guidance control system. Nonetheless, as of now, research addressing this specific problem remains insufficient.

In summary, this paper takes a certain type of gliding-guided projectile as a case study to tackle the trajectory planning issue for striking fixed target. To reduce the trajectory's sensitivity to uncertainties and enhance the compatibility between the planned trajectory and the guidance control system, thereby maximizing its limited control capability, a closed-loop robust trajectory planning method based on the concepts of robust optimization and tracking guidance is proposed. Initially, the main uncertainties during the actual flight process of the gliding projectile are modeled, and the non-invasive polynomial chaos expansion (NIPCE) method is used to obtain quantified propagation laws for these uncertainties. Subsequently, based on these propagation results and aiming to minimize the terminal dispersion as the performance indicator, open-loop robust trajectory planning is carried out, considering the ballistic characteristics of the gliding projectile and related constraints. Furthermore, a simple yet effective Proportional-Integral-Differential (PID) controller is employed to track the robust optimal trajectory. Based on feedback from the guidance control system, the planning approach is extended to a closed-loop format, designing a closed-loop robust optimal trajectory that integrates both uncertainty factors and the 'scheme-guidance compatibility'. Finally, through numerical simulations in various scenarios and different conditions, the impact of uncertainty deviations on the performance of closed-loop guidance is explored.

2. Models and Methods

2.1. Dynamic Model of Gliding Ballistic

To mitigate the adverse effects of factors such as aerodynamic asymmetry and engine thrust misalignment, while reducing the cost of the guidance control system, the tail fins of gliding-guided projectiles are typically designed with a certain angle of inclination, allowing the projectile to spin slowly during flight. Consequently, only pitch and yaw control are required. In the initial design process, to facilitate the determination of the gliding projectile's flight trajectory and main flight characteristics, its rotational motion around the center of mass is often neglected, treating it as a maneuverable point mass for study[28]. This paper is based on the following simplified assumptions.

Assumption 

1. The gliding projectile rotates without inertia.

Assumption 

2. The control system of the gliding projectile works ideally, without errors or delays.

The aforementioned assumptions posit that when the control mechanisms of the gliding projectile deflect, the torque acting on the projectile is in equilibrium at each instant, hence it is termed the 'instantaneous equilibrium' assumption. Given its axisymmetric shape, when the total angle of attack is not excessively large, the gliding projectile exhibits linear aerodynamic characteristics.

$$\{\begin{array}{c}{\left(\frac{{\delta }_z}{\alpha }\right)}_{\mathrm{E}}=-\frac{m_z^{\alpha }}{m_z^{{\delta }_z}}\\ {\left(\frac{{\delta }_y}{\beta }\right)}_{\mathrm{E}}=-\frac{m_y^{\beta }}{m_y^{{\delta }_y}}\end{array}$$

(1)

In the equations, ${\delta }_z$ and ${\delta }_y$ represent the deflection angles of the pitch and yaw control surfaces, respectively. $\alpha$and $\beta$ denote the angle of attack and sideslip of the projectile, with subscript E indicating the equilibrium state. ​$m_z^{\alpha }$ and $m_y^{\beta }$ represent the pitch and yaw static stability derivatives, while $m_z^{{\delta }_z}$​​ and $m_y^{{\delta }_y}$​​ denote the control effectiveness of the pitch and yaw control surfaces, respectively.

Under the instantaneous equilibrium assumption, the control variables for the gliding trajectory can be approximated using the angle of attack $\alpha$ and sideslip $\beta$ at equilibrium. Gliding-guided projectiles undergo unpowered flight during the controlled phase, relying on the aerodynamic forces acting on the projectile to achieve maneuver of the center of mass.

$$\{\begin{array}{l}{F}_x=qS{C}_{x0}\left[1+k_{\mathrm{c}}\left({\alpha }^2+{\beta }^2\right)\right]\hfill \\ F_y=qS{C}_y^{\alpha }\alpha \hfill \\ F_z=qS{C}_z^{\beta }\beta \hfill \end{array}$$

(2)

In the equations, $F_x$, $F_y$, and $F_z$ represent the total drag, lift, and lateral force acting on the projectile, respectively. $q=\rho V^2/2$denotes the dynamic pressure, $\rho$is the air density, $V$is the velocity of the projectile, $S$is the reference area of the projectile, $C_{x0}$is the zero-lift drag coefficient, $k_{\mathrm{c}}$is the induced drag coefficient, $C_y^{\alpha }$and $C_z^{\beta }$ are the partial derivatives of the lift and lateral force coefficients with respect to angle of attack and sideslip, respectively. For axisymmetric shapes,$C_z^{\beta }=-C_y^{\alpha }$.

The centroid motion model of the gliding trajectory is

$$\{\begin{array}{l}-m\dot{V}=mg\sin \theta +F_x\hfill \\ mV\dot{\theta }=-mg\cos \theta +F_y\hfill \\ -\left(mV\cos \theta \right)\dot{\psi }=F_z\hfill \\ \dot{x}=V\cos \theta \cos \psi \hfill \\ \dot{y}=V\sin \theta \hfill \\ \dot{z}=-V\cos \theta \sin \psi \hfill \end{array}$$

(3)

In the equations, $m$is the mass of the projectile, $g$is the gravitational acceleration, $\theta$and $\psi$are the inclination and deflection angle, respectively. $x$, $y$and$z$are the range, altitude and lateral displacement of the projectile, respectively.

Considering that the impact time of the gliding trajectory is free, and the target position to be hit is fixed, the range of the projectile is monotonously increasing and bounded in the course of flight. In order to facilitate the follow-up analysis, the independent variable of the model is converted from time t to range x, so the motion model shown in equation (3) is converted to

$$\{\begin{array}{l}\dot{V}=-\frac{mg\sin \theta +F_x}{mV\cos \theta \cos \psi }\hfill \\ \dot{\theta }=\frac{-mg\cos \theta +F_y}{m{V}^2\cos \theta \cos \psi }\hfill \\ \dot{\psi }=\frac{-F_z}{m{V}^2{\cos }^2\theta \cos \psi }\hfill \\ \dot{t}=\frac{1}{V\cos \theta \cos \psi }\hfill \\ \dot{y}=\frac{\tan \theta }{\cos \psi }\hfill \\ \dot{z}=-\tan \psi \hfill \end{array}$$

(4)

Remarking the combination of state variables as $X={[V,\theta ,\psi ,t,y,z]}^{\mathrm{T}}$ and the combination of control variables as $U={[\alpha ,\beta ]}^{\mathrm{T}}$, the equation (4) can be expressed as

$$\dot{X}=f\left(x,X,U\right)$$

(5)

2.2. Model of Major Uncertainties

There exist all kinds of uncertain factors in the actual flight environment, and it is not realistic to take them all into consideration. Therefore, this paper selects the following representative uncertain factors for research, and models them as the deviation between the real value and the reference value.

2.2.1. Aerodynamic Parameter Deviation

Due to the machining errors of the projectile and variations in the ballistic environment, there are inevitably differences between the actual flight aerodynamic parameters and the reference values. This paper assumes that the deviations of the aerodynamic parameters follow mutually independent normal distributions.

$$\{\begin{array}{l}{\tilde{C}}_{x0}=\left(1+N_d{\sigma }_d\right)C_{x0}\hfill \\ {\tilde{C}}_y^{\alpha }=\left(1+N_l{\sigma }_l\right)C_y^{\alpha }\hfill \\ {\tilde{C}}_z^{\beta }=-\left(1+N_{l^{\prime }}{\sigma }_{l^{\prime }}\right)C_y^{\alpha }\hfill \end{array}$$

(6)

In the equations, ${\tilde{C}}_{x0}$, ${\tilde{C}}_y^{\alpha }$ and ${\tilde{C}}_z^{\beta }$ represent the true values of zero-lift drag coefficient, lift coefficient derivative, and lateral force coefficient derivative, respectively. Subscripts$d$, $l$and $l^{\prime }$denote the first letters of drag, lift, and lateral, respectively.$N_d$, $N_l$, and $N_{l^{\prime }}$ are random numbers following a standard normal distribution N(0,1), and${\sigma }_d$, ${\sigma }_l$, and ${\sigma }_{l^{\prime }}$ represent the degrees to which the above parameters deviate from their reference values.

2.2.2. Meteorological Environment Deviation

The differences between actual meteorology and standard meteorology are mainly manifested in the deviations of wind field, temperature, and atmospheric density. The impact of wind field on trajectory calculation can be considered as an unmodeled random force, which is compensated by the guidance control system. Temperature primarily affects the calculation of speed of sound, thereby influencing the flight Mach number, which in turn leads to variations in aerodynamic parameters. Therefore, the deviation in temperature can be attributed to the deviation in aerodynamic parameters. In summary, this paper models the deviation in meteorological environment as the deviation in atmospheric density. The statistical characteristics of the deviation of actual atmospheric density from the standard value are related to the altitude y, and the ratio of standard deviation to standard atmospheric density can be approximated as

$${\sigma }_{\rho }=0.003517{\mathrm{e}}^{\frac{y}{26629.77}}$$

(7)

Therefore, the actual atmospheric density can be expressed as

$$\tilde{\rho }=\left(1+N_{\rho }{\sigma }_{\rho }\right)\rho$$

(8)

In the formula, $N_{\rho }$is a random number obeying the standard normal distribution N(0,1).

2.2.3. State Deviations at the Control-Start Point

Before terminal guidance begins, the full trajectory of a gliding guided projectile can be divided into four stages: the launching phase, boosting phase, climbing phase, and gliding phase, with the first three stages being uncontrolled flight[29]. After the projectile is fired, under the influence of various uncertainties, the state values at the control-start point will inevitably deviate from the designed trajectory. This paper assumes that the state deviations follow mutually independent normal distributions. Remarking $X_0={[V_0,{\theta }_0,{\psi }_0,t_0,y_0,z_0]}^{\mathrm{T}}$ as the design value of the trajectory at the control-start point, then the actual state at the control-start point can be expressed as

$${\tilde{X}}_0=\left[\begin{array}{l}{\tilde{V}}_0\hfill \\ {\tilde{\theta }}_0\hfill \\ {\tilde{\psi }}_0\hfill \\ {\tilde{t}}_0\hfill \\ {\tilde{y}}_0\hfill \\ {\tilde{z}}_0\hfill \end{array}\right]=\left[\begin{array}{l}{V}_0+N_V{\sigma }_V\hfill \\ {\theta }_0+N_{\theta }{\sigma }_{\theta }\hfill \\ {\psi }_0+N_{\psi }{\sigma }_{\psi }\hfill \\ t_0+N_t{\sigma }_t\hfill \\ y_0+N_y{\sigma }_y\hfill \\ z_0+N_z{\sigma }_z\hfill \end{array}\right]$$

(9)

In the equations, $N_V$, $N_{\theta }$, $N_{\psi }$, $N_t$, $N_y$and $N_z$are random numbers obeying the standard normal distribution N(0,1), and${\sigma }_V$, ${\sigma }_{\theta }$, ${\sigma }_{\psi }$, ${\sigma }_t$, ${\sigma }_y$ and ${\sigma }_z$ represent the standard deviation of each state at the control-start point respectively.

2.3. Uncertainty Propagation of the Gliding Trajectory

2.3.1. Transformation of the Stochastic Dynamic Model

Based on the discussion in Section 2.2, the main uncertainties during the flight of a gliding projectile can be categorized into two types: parameter uncertainties and initial state uncertainties. The dimension of parameter uncertainties$d_P=4$, while the dimension of initial state uncertainties$d_X=6$. Therefore, the dimension of uncertainty factors for the entire stochastic dynamic system$d=10$. Remarking the random parameters as $P={[{\sigma }_d,{\sigma }_l,{\sigma }_{l^{\prime }},{\sigma }_{\rho }]}^{\mathrm{T}}$ and introduce them into the dynamic system shown in Equation (4), then the stochastic dynamic model can be expressed as

$$\dot{X}=f\left(x,X,U,P\right)$$

(10)

This paper adopts the NIPCE method based on Gaussian quadrature nodes sampling and Weighted Stochastic Response Surface Method (WSRSM)[30] to transform Equation (10) into a high-dimensional deterministic dynamic model with PCE coefficients as state variables. Setting the order of PCE as$p=2$, then the number of PCE coefficients corresponding to each state variable after expansion is

$$P+1=\frac{\left(p+d\right)!}{p!d!}=66$$

(11)

The number of sampling points used for regression is

$$n_{\mathrm{s}}=2\left(P+1\right)=132$$

(12)

According to the theory of generalized polynomial chaos, each state variable in formula (10) can be expanded in the form of a polynomial.

$$X_i\approx {\displaystyle \sum _{j=0}^P{S}_{ij}{\Phi }_j\left(\xi \right)},i=1,2,\dots ,d_X$$

(13)

In the equation, $\xi \in R^{10}$is a $d$-dimensional cell array, where each element consists of $n_{\mathrm{s}}$ Gaussian quadrature sampling points corresponding to the respective uncertain factors. Let the univariate orthogonal basis functions corresponding to each uncertain factor in each dimension be denoted by${\varphi }_j({\xi }_k),k=1,2,\dots ,d$, then ${\Phi }_j(\xi )$ represents the orthogonal polynomial function obtained by taking tensor products of${\varphi }_j({\xi }_k)$. $S_{ij}$denotes the j-th PCE coefficient of the i-th state variable, and the total number of coefficients is

$$n_{\mathrm{PCE}}=d_X\left(P+1\right)=396$$

(14)

For the random parameters which obey the normal distribution, the first-order, one-dimensional orthogonal polynomial model can accurately describe the uncertainty[31]. Therefore, each random parameter can be expanded in the form of the sum of 2 terms.

$$P_i=p_{i0}{\Phi }_0\left({\xi }_P\right)+p_{i1}{\Phi }_1\left({\xi }_P\right),i=1,2,\dots ,d_P$$

(15)

In the equation, ${\xi }_P\in R^4$is a $d_P$-dimensional cell array, where each element consists of $n_{\mathrm{s}}$ Gaussian quadrature sampling points corresponding to the respective uncertain factors.

To sum up, the stochastic system shown in formula (10) is transformed into

$$\Psi {\dot{S}}_{\mathrm{PCE}}=F\left(x,S_{\mathrm{PCE}},U,P_{\mathrm{PCE}}\right)$$

(16)

In the formula, $S_{\mathrm{PCE}}\in R^{396}$denotes the PCE coefficients of all state variables, $P_{\mathrm{PCE}}\in R^8$denotes the PCE coefficients of random parameters. $F\in R^{792}$is the actual response values of $n_{\mathrm{s}}$ sample points into the original stochastic system formula (10), and $\Psi \in R^{792\times 396}$denotes the coefficient matrix of the transformed deterministic system which is described in detail as follows.

Remarking$\Omega ({\xi }_i)=[{\Phi }_0({\xi }_i),{\Phi }_1({\xi }_i),\dots {\Phi }_P({\xi }_i)],i=1,2,\dots ,n_{\mathrm{s}}$，defining the intermediate matrix as

$$\Theta ({\xi }_i)\in R^{6\times 396}=I_{d_X}\otimes \Omega ({\xi }_i)$$

(17)

In the equation, the symbol $\otimes$ denotes the Cronecker product, therefore

$$\Psi ={[\Theta {({\xi }_1)}^{\mathrm{T}},\Theta {({\xi }_2)}^{\mathrm{T}},\dots \Theta {({\xi }_{n_{\mathrm{s}}})}^{\mathrm{T}}]}^{\mathrm{T}}$$

(18)

Considering the different importance of each sample point, the weight vector is

$${\omega }_{\mathrm{s}}\in R^{792}={\left[\underset{\mathrm{a}\mathrm{total}\mathrm{of}{d}_X\mathrm{items}}{\underbrace{{\omega }^{\mathrm{T}},{\omega }^{\mathrm{T}},\dots ,{\omega }^{\mathrm{T}}}}\right]}^{\mathrm{T}}$$

(19)

In the formula, $\omega$ denotes the weight corresponding to $n_{\mathrm{s}}$ sample points.

So far, according to WSRSM, the overdetermined equations shown in formula (16) can be solved by weighted least square regression.

$$\{\begin{array}{l}H={\left({\Psi }^{\mathrm{T}}{W}_{\mathrm{s}}\Psi \right)}^{-1}{\Psi }^{\mathrm{T}}{W}_{\mathrm{s}}\hfill \\ {\dot{S}}_{\mathrm{PCE}}=HF\left(x,S_{\mathrm{PCE}},U,P_{\mathrm{PCE}}\right)\hfill \end{array}$$

(20)

Where $H\in R^{396\times 792}$is the regression coefficient matrix and $W_{\mathrm{s}}\in R^{792\times 792}=\mathrm{diag}({\omega }_{\mathrm{s}})$ denotes the weight matrix.

By numerical integration of equation (20), the variation of all PCE coefficients with range $x$ can be obtained. The mean $\mu$ and standard deviation $\sigma$ of the state variables in the original stochastic system can be obtained directly from the PCE coefficients.

$$\{\begin{array}{l}\mu \left(X_i\right)=S_{i0}\hfill \\ \sigma \left(X_i\right)=\sqrt{{\displaystyle \sum _{j=1}^P{S}_{ij}^2\mathrm{E}\left[{\Phi }_j^2\left(\xi \right)\right]}}\hfill \end{array},i=1,2,\dots ,d_X$$

(21)

Where the operator E represents the mathematical expectation.

2.3.2. Truncation of Chaotic Polynomial Basis

According to the description in Section 2.3.1, the transformed deterministic system has a total of 396 state variables, and there are also 132 Gaussian quadrature sampling points used for regression. In this problem scale, the computation time for a single uncertainty propagation is acceptable. However, in subsequent robust trajectory planning of the gliding projectile, the uncertainty propagation module will be repeatedly called during the iterative process. Moreover, the convergence speed of optimization algorithms significantly decreases when the number of variables is large, leading to a significant reduction in computational efficiency.

It is known that as the dimension $d$ of the uncertainty factors increases, the total number of terms in the orthogonal polynomials grows exponentially. However, in many problems, not all terms in the polynomial expansion are equally important. Typically, only a few variables are involved in the important terms of the expansion. In most practical engineering systems, the output response is mainly influenced by the input variables and their low-order interaction terms, while the impact of higher-order interaction terms is relatively small.

Considering the above analysis, this paper adopts a maximum interaction restriction truncation strategy to directly remove high-order polynomial terms, thereby reducing the computational burden.

In the full-order PCE model (13), multivariate orthogonal polynomials are constructed by direct tensor product. When the PCE order is$p$, it is necessary that

$$A^{d,p}=\left\{\lambda \in N^d:{\Vert \lambda \Vert }_1\equiv {\displaystyle \sum _{j=1}^d{\lambda }_j⩽p}\right\}$$

(22)

In the formula, $\lambda =[{\lambda }_1,{\lambda }_2,\dots ,{\lambda }_d]$represents multiple indexes, $A^{d,p}$ is a set of multiple indexes, ${\lambda }_i$ is the order of the base of univariate orthogonal polynomials corresponding to the j-dimensional variables in multivariate orthogonal polynomials, and the operator ${\Vert \cdot \Vert }_1$ represents 1-norm.

The main idea of maximum interaction restriction truncation is to select a low rank subset $A^{d,p,r}$ from the multi-index set defined in formula (22) to ensure not only ${\Vert \lambda \Vert }_1⩽p$ but also${\Vert \lambda \Vert }_0⩽r$.

$$\{\begin{array}{l}{A}^{d,p,r}=\left\{\lambda \in A^{d,p}:\Vert \lambda \Vert \equiv {\Vert \lambda \Vert }_0⩽p\right\}\hfill \\ 1⩽r<\min \left(d,p\right)\hfill \end{array}$$

(23)

In the equations, ${\Vert \lambda \Vert }_0$represents the rank of multiple indexes, r represents the truncation coefficient.

For the case study in this paper, $d=10$and$p=2$. Therefore, based on Equation (23), the truncation coefficient is set to$r=1$. After removing the high-order interaction terms, the number of PCE coefficients for each dimension of random variables decreases from 66 at full order to 11. Consequently, the total number of state variables in the deterministic dynamic system after dimensionality reduction decreases from 396 to 66. The number of Gaussian quadrature sampling points used for regression reduces from 132 to 22, representing a reduction of 83.3%.

2.4. Robust Planning Model of Gliding Trajectories

2.4.1. Open-Loop Robust Planning

• Planning expectation

With the control variables determined, through the method introduced in Section 2.3 one can accurately and efficiently get the propagation laws of uncertainties along the gliding trajectory. This is manifested as deviations and dispersion of the trajectory envelope during the flight process and terminal states. By incorporating the results of uncertainty propagation into the objective function and constraints during trajectory planning, the robustness of the trajectory can be enhanced in the iterative process, reducing the sensitivity of the gliding trajectory to random disturbances.

Theoretically, statistical information about the terminal states of the trajectory can be obtained through the PCE agent model, allowing the covariance ellipse of the terminal states to be determined. Bringing back the control variables of the robust optimal trajectory to the original stochastic dynamic model for open-loop tracking enables the true distribution of the gliding trajectory's terminal position to be obtained. By observing the coverage of the theoretical covariance ellipse on the actual terminal position, the effectiveness of robust planning can be validated, with the expected results depicted in Figure 1.

2.

Dynamics model

The dynamic model is the deterministic system shown in equation (16) after dimensional expansion, the state variables are all PCE coefficients, and the control variables are consistent with the original stochastic system.

3.

Constraints

After considering the influence of uncertain factors, there are numerous results (bounded) of flight trajectory under the same open-loop control input. Therefore, the focus of open-loop robust planning is not the specific values of the states, but their statistical laws, and the statistical moments of state variables should also be considered in the constraints (this paper mainly focuses on the first and second order statistical moments, which are mean and variance).

• Boundary constraints

According to the analysis in Section 2.2.3, the state deviations at the control-start point obey the independent normal distribution. Therefore, the constraints of the mean of initial states are

$$\{\begin{array}{l}\mathrm{E}\left[V\left(x_0\right)\right]=V_0\hfill \\ \mathrm{E}\left[\theta \left(x_0\right)\right]={\theta }_0\hfill \\ \mathrm{E}\left[\psi \left(x_0\right)\right]={\psi }_0\hfill \\ \mathrm{E}\left[t\left(x_0\right)\right]=t_0\hfill \\ \mathrm{E}\left[y\left(x_0\right)\right]=y_0\hfill \\ \mathrm{E}\left[z\left(x_0\right)\right]=z_0\hfill \end{array}$$

(24)

Where $x_0=0$ is the $x$ coordinate at the control-start point.

In order to ensure the final strike accuracy and strike effect, the constraints on the mean values of the terminal states are

$$\{\begin{array}{l}\mathrm{E}\left[V\left(x_{\mathrm{f}}\right)\right]⩾V_{\mathrm{f}\min }\hfill \\ \mathrm{E}\left[y\left(x_{\mathrm{f}}\right)\right]=y_{\mathrm{T}}\hfill \\ \mathrm{E}\left[z\left(x_{\mathrm{f}}\right)\right]=z_{\mathrm{T}}\hfill \end{array}$$

(25)

In the formulas, the subscript f represents the first letter of final, and the subscript T represents the first letter of Target, $x_{\mathrm{f}}=x_{\mathrm{T}}$and $(x_{\mathrm{T}},y_{\mathrm{T}},z_{\mathrm{T}})$represents the fixed target position.

• Path constraints

Considering the limited rudder control capability of the gliding projectile, the control constraints during flight are

$$\{\begin{array}{l}\left|\alpha \right|⩽{\alpha }_{\max }\hfill \\ \left|\beta \right|⩽{\beta }_{\max }\hfill \end{array}$$

(26)

4.

Objective function

For the planned trajectory of a gliding projectile, the dispersion at the terminal position is usually used as a criterion to measure its sensitivity to uncertainty, while less attention is paid to the dispersion of other terminal states and the state envelope during the flight. Therefore, the objective function of the open-loop robust planning is designed as

$$\min J_{\mathrm{o}}=\mathrm{D}\left(y_{\mathrm{f}}\right)+\mathrm{D}\left(z_{\mathrm{f}}\right)$$

(27)

In the formula, the subscript o represents the first letter of open-loop, and the operator D represents the variance.

2.4.2. Closed-Loop Robust Planning

• Planning expectation

Closed-loop robust planning builds upon open-loop planning by incorporating a tracking guidance module. It assumes that all uncertainties are biased towards the semi-extreme (which means the $\pm 1.5\sigma$ value of each random parameter in this paper). The tracking of robust optimal trajectories provides closed-loop feedback for control variables, thus obtaining control effort consumption during the projectile's flight (based on simulation experiences and measured data, the condition where all uncertainties are biased to the semi-extreme is quite stringent, making trajectory tracking relatively challenging and control effort consumption higher. Adopting such a conservative approach during the planning process can generally ensure the smooth progress of closed-loop guidance during actual flight). By incorporating control effort consumption as a penalty function into the planning objective function, the iterative process can then consider the coupling and matching characteristics between the planned trajectory and the guidance system, enabling the planning of a closed-loop robust optimal trajectory that balances robustness and feasibility.

2.

Closed-loop guidance based on PID control

Due to the repeated invocation of the guidance module during the iterative process of closed-loop planning algorithms, a sufficiently simple and effective trajectory tracking algorithm is necessary. In the existing literature, some scholars employed the Linear Quadratic Regulator (LQR) algorithm for tracking programmed trajectories[17]. This method involves linearizing the model near the reference trajectory, which can lead to significant tracking errors when deviations are large. Furthermore, the LQR algorithm requires solving the Riccati equation, increasing the complexity of the solution. Based on the aforementioned analysis, this paper utilizes a PID controller with an analytic expression that does not rely on the model for tracking the planned trajectory. Let the altitude and lateral displacement at the range x on the reference trajectory be denoted as $\widehat{y}$ and $\widehat{z}$ respectively, and the control variables as $\widehat{U}$, then the position error of the projectile relative to the reference trajectory can be expressed as

$$e=\left[\begin{array}{l}{e}_y\hfill \\ e_z\hfill \end{array}\right]=\left[\begin{array}{l}\widehat{y}-y\hfill \\ \widehat{z}-z\hfill \end{array}\right]$$

(28)

The closed-loop control commands after feedback modification can be expressed as

$$U_{\mathrm{c}}=\left[\begin{array}{l}{\alpha }_{\mathrm{c}}\hfill \\ {\beta }_{\mathrm{c}}\hfill \end{array}\right]=\widehat{U}+\Delta U=\widehat{U}+\left[\begin{array}{l}\Delta \alpha \hfill \\ -\Delta \beta \hfill \end{array}\right]=\widehat{U}+k_{\mathrm{P}}e+k_{\mathrm{I}}{\displaystyle \underset{x_0}{\overset{x}{\int }}e\mathrm{d}\tau }+k_{\mathrm{D}}\dot{e}$$

(29)

In the formula, the subscript c represents the first letter of closed-loop.$k_{\mathrm{P}}$, $k_{\mathrm{I}}$and $k_{\mathrm{D}}$ represent proportional parameter, integral parameter and differential parameter, respectively. The minus sign in front of $\Delta \beta$ item indicates the direction.

3.

Dynamic model and constraints

The dynamic model and related constraints of closed-loop robust planning are the same as those in open-loop planning, which need not be repeated.

4.

Objective function

Defining the closed-loop control effort consumption of the projectile during flight as

$$E_{\mathrm{c}}={\displaystyle \underset{x_0}{\overset{x_{\mathrm{f}}}{\int }}\left({\alpha }_{\mathrm{c}}^2+{\beta }_{\mathrm{c}}^2\right)\mathrm{d}x}$$

(30)

Then the objective function of closed-loop robust planning can be designed as

$$\min J_{\mathrm{c}}={\omega }_{\mathrm{D}}{J}_{\mathrm{o}}+E_{\mathrm{c}}$$

(31)

Where ${\omega }_{\mathrm{D}}$ is the weight of the terminal position dispersion.

The closed-loop robust trajectory planning is a multi-objective optimal control problem. Considering that the order of magnitude difference between $J_{\mathrm{o}}$ and $E_{\mathrm{c}}$ in the calculation example is large, in order to achieve the balance of planning emphasis, this paper uses the deterministic energy-optimal trajectory without considering random disturbance as the benchmark trajectory to make the objective function dimensionless. Remarking the terminal position dispersion of the benchmark trajectory after uncertainty propagation as $\overline{J}$, the control effort consumption as $\overline{E}$, then the objective function is converted into

$$\min J_{\mathrm{c}}={\omega }_{\mathrm{D}}\frac{J_{\mathrm{o}}}{\overline{J}}+\frac{E_{\mathrm{c}}}{\overline{E}}$$

(32)

2.4.3. Flow of Robust Gliding Trajectory Planning

In summary, the open-loop robust planning problem of gliding trajectory can be described as follows: considering the main uncertain factors discussed in Section 2.2 and the high-dimensional deterministic dynamic model shown in equation (16), planning the open-loop control commands $U_{\mathrm{o}}$within a known range $[x_0,x_{\mathrm{f}}]$ based on the current battlefield environment and task requirements. After uncertain factors propagate along the gliding trajectory under the action of $U_{\mathrm{o}}$, the flight scheme can minimize the objective function shown in equation (27) while satisfying the constraint equations (24) to (26). The closed-loop robust planning problem can be described as follows: taking the open-loop robust optimal trajectory as the reference trajectory, considering the case where all uncertain factors are shifted to semi-extreme, obtaining the feedback-modified closed-loop control commands $U_{\mathrm{c}}$ through the PID controller shown in equation (29) and calculating the control effort consumption $E_{\mathrm{c}}$, incorporating $E_{\mathrm{c}}$ into the planning objective in the form of penalty functions. Using the deterministic energy-optimal trajectory as the benchmark, non-dimensionalizing the objective function as shown in equation (32), obtaining the closed-loop robust optimal trajectory through a bi-level planning format. The flow chart of gliding trajectory robust planning proposed in this paper is shown as Figure 2.

3. Simulation Results and Analysis

Due to the large number of states in the deterministic system after dimension-expansion and the presence of many zero elements in the PCE coefficients, the problem scale is large. Discretizing state variables will cause unnecessary computational time consumption. Meanwhile, the uncertainty propagation module only needs control input information to obtain the statistical laws of the system states that meet the accuracy requirements. Therefore, this paper adopts a direct shooting method with only discrete control variables to solve the robust planning problem of gliding trajectory, with the number of discrete points N set to 65. The relevant parameters of the projectile and target are shown in Table 1, and the parameters related to uncertainty are shown in Table 2. The parameters of the PID controller are adjusted and adapted during the simulation process, set to$k_{\mathrm{P}}=1.0\times {10}^{-2}$,$k_{\mathrm{I}}=1.5\times {10}^{-10}$ and$k_{\mathrm{D}}=1.2\times {10}^{-1}$.

3.1. Uncertainty Propagation of Gliding Trajectory under Specified Control Commands

When uncertainty factors are not considered, conducting trajectory planning based on minimizing control effort consumption as the objective function can yield the deterministic optimal control commands $\widehat{U}$, which is not detailed in this paper. To verify the effectiveness and superiority of the NIPCE method proposed in this paper, according to the discussion in Section 2.3, substituting $\widehat{U}$ into the high-dimensional system equation (16) for uncertainty propagation, the simulation results are shown in Figure 3 (where the symbols $\mu (\cdot )$ and $\sigma (\cdot )$ represent the mean value and standard deviation of the state variables, the black dotted lines in the figures represent the constraint thresholds of the state variables, and in the legend, FOPCE and COPCE respectively denote the full-order PCE and cutoff-order PCE).

From the results in Figure 3(a), (i), and (k), it is obvious that under the influence of various uncertain factors, the mean value of state at the end of the trajectory violates the constraint conditions indicated in equation (25). From Figure 3(j) and (l), it is observed that there is a significant dispersion in the terminal position of the trajectory. Therefore, robust planning of the gliding trajectory is required to modify the deviations in state and reduce the impact of uncertainty factors. Furthermore, by observing all subplots in Figure 3, it is noted that the results of the FOPCE method are highly consistent with the MCS results. When using the basis-truncation strategy (corresponding to the COPCE method) to eliminate high-order polynomial terms, the standard deviations of the terms$V$, $\theta$ and $\psi$ with relatively high nonlinearity (Figure 3(b), (d), (f)) in Equation (4) show some deviation from the MCS results, but within an acceptable range of accuracy. The deviation occurs only during intermediate processes, with the statistical characteristics at the end of the trajectory remaining unaffected. This analysis illustrates that the NIPCE method proposed in this paper is competent for quantifying the uncertainties of gliding trajectories, and in subsequent robust planning processes, a truncated PCE model can be used instead of a full-order PCE model to enhance computational efficiency.

The simulation in this paper is based on a personal computer with 4 cores CPU E3-1230 V2 3.30GHz. The calculation times of single uncertainty propagation of the gliding trajectory using different methods on the MATLAB 2022b platform are shown in Table 3.

The results from Table 3 indicate that compared to the MCS method, the FOPCE method can significantly improve computational efficiency while ensuring accuracy. The FOPCE method with single-core computation takes less time than MCS with multi-core parallel computation, achieving a reduction of 79.5% in time consumption. Employing a basis-truncation strategy to remove unnecessary high-order cross-terms can further reduce problem scale and computational time, facilitating the optimization process of robust planning. Under single-core computation the same, the computation time of COPCE is reduced by 84.2% compared to FOPCE.

3.2. Open-Loop Robust Planning

Based on the discussion in section 2.4.1, using the deterministic effort-optimal control commands $\widehat{U}$ as the initial value for iteration, the open-loop robust planning results of the gliding trajectory are shown in Figure 4 (in the legend, DEO-Trj represents the results of uncertainty propagation along the trajectory under the influence of $\widehat{U}$, and ORO-Trj represents the uncertainty propagation results of the open-loop robust optimal control commands $U_{\mathrm{o}}$).

Based on Figure 4(a), (i) and (k), the open-loop robust planning corrects the bias in the mean value of terminal states of the trajectory, ensuring compliance with the constraints outlined in Equation (25). Figure 4(a), (g), (i) and (k) show that the curvature in the middle section of the robust optimal trajectory increases compared to the reference trajectory, resulting in reduced projectile velocity and a longer flight time. This implies that increasing the mid-trajectory curvature enhances the robustness of the planned trajectory, aligning with the conclusions drawn using the LAC method in reference[2]. Figure 4(j) and (l) reveal that the open-loop robust planning significantly reduces the dispersion at the trajectory's terminal point, with standard deviations of altitude and lateral deviation decreasing by 23.6% and 35.3%, respectively. This effectively minimizes the sensitivity of the planned trajectory to uncertainties. However, as observed in Figure 4(m) and (n), the control commands of the robust optimal trajectory significantly increase compared to the reference trajectory, indicating that the trajectory's robustness is achieved at the expense of additional control effort consumption. Moreover, the results demonstrate that the projectile remains at control saturation for a considerable period, posing significant challenges to subsequent guidance control system design. If the deviation of uncertainties is substantial, a gliding projectile with limited maneuverability may fail to track the planned trajectory, leading to mission failure. Therefore, considering the alignment between the planned trajectory and the guidance control system during planning and conducting closed-loop robust trajectory planning research holds significant practical importance.

To further validate the effectiveness of the open-loop robust planning, the robust optimal control commands $U_{\mathrm{o}}$ are applied to the original stochastic dynamics model expressed in Equation (4) for 10,000 Monte Carlo simulations. The comparison between the actual distribution of the terminal positions of the gliding trajectory and the theoretical covariance ellipse is illustrated in Figure 5.

As from the results in Figure 5, it can be indicated that by the theoretical covariance ellipse adequately encompasses the actual terminal positions the of the trajectory, with the target located at the center of the ellipse, results consistent with the planning expectations in Section 2.4.1 (It should be pointed out that there is a situation where the altitude y is negative at the terminal position of the trajectory in the figure, which is not an erroneous result, but a statistical law obtained by the model through uncertainty propagation when the independent variable is x. Through subsequent closed-loop robust planning, the projectile will precisely hit the target). However, due to the limited control capabilities of the projectile and the absence of closed-loop feedback from a guidance control system, the open-loop robust planning still exhibits some dispersion at the terminal phase, even with minimized objective function. This analysis underscores the necessity of pursuing research in closed-loop robust planning.

3.3. Closed-Loop Robust Planning

Based on the discussion in Section 2.4.2, considering the scenario where uncertainties are skewed to the semi-extreme values, a PID controller is used to track the open-loop robust optimal trajectory, providing closed-loop feedback for control commands. The closed-loop control effort consumption is calculated and incorporated into the objective function. This paper employs the weight of terminal position dispersion ${\omega }_{\mathrm{D}}$ in Equation (32) as the basis for distinguishing the following three scenarios.

Scenario 1：${\omega }_{\mathrm{D}}=0$，which indicates that the terminal dispersion is allowed, only the state deviations are modified, and the control effort consumption is minimized.

Scenario 2：${\omega }_{\mathrm{D}}=1$，which means that both terminal dispersion and control effort consumption are taken into consideration, and both are equally important.

Scenario 3：${\omega }_{\mathrm{D}}=5$，which means taking both the terminal dispersion and the control effort consumption into consideration, but the former is more important.

The results of closed-loop robust planning of gliding trajectory in different scenarios are shown in Figure 6.

From the results in Figure 6, it is evident that with the introduction of a guidance control module, the closed-loop robust planning outcome is no longer a statistical law. The projectile, under the influence of the guidance law, precisely hits the target, eliminating terminal dispersion. Even under relatively stringent deviation conditions, the projectile still manages to track the reference trajectory with minimal control saturation (It can be imagined that with smaller deviations, the projectile's control effort consumption would be lower). This proves the effectiveness and superiority of closed-loop planning.

As the terminal dispersion weight ${\omega }_{\mathrm{D}}$ increases, the curvature in the middle section of the closed-loop robust optimal trajectory correspondingly increases, gradually approaching the open-loop robust optimal trajectory (as shown in Figure 4(i) and 4(k)), and the control commands also increase accordingly. It can be inferred that when ${\omega }_{\mathrm{D}}$ is sufficiently large, the proportion of control effort consumption in the objective function Equation (32) will approach zero. The closed-loop robust planning will degenerate into open-loop robust planning, and the closed-loop robust optimal trajectory will converge to the open-loop robust optimal trajectory.

In order to compare the results in different scenarios more intuitively, the specific values ​​of each item in the closed-loop planning objective function Equation (32) are shown in Table 4.

From the results in Table 4, it is evident that when considering only the correction of state deviation without accounting for terminal dispersion of the trajectory, the difference in control effort consumption compared to the deterministic optimal trajectory is relatively small. As ${\omega }_{\mathrm{D}}$ increases, the robustness of the trajectory improves, accompanied by a corresponding increase in control effort consumption. Therefore, for gliding-guided projectiles, a trade-off needs to be made between robust optimality and energy optimality, fully exploiting their limited control capabilities to maximize comprehensive performance.

3.4 The Influence of the Uncertain factors' Deviation Degree on Closed-Loop Guidance

In the process of closed-loop robust trajectory planning, all uncertain factors are biased to the semi-extreme situation, which is a conservative approach. According to relevant simulation experiences and measured data, the degree of deviation of each random parameter during the actual flight of the gliding projectile is not the same. Therefore, based on the degree of deviation of uncertain factors, this paper distinguishes the following four conditions.

Condition 1: The deviations of uncertain factors are all set to corresponding $+3\sigma$ values ​​to simulate the upper bound of the extreme deviation situation.

Condition 2: The deviations of uncertain factors are alternately set to the corresponding $+\sigma$ and $-\sigma$ values ​​to simulate the actual flight environment.

Condition 3: The deviations of uncertain factors are all set to corresponding $-3\sigma$ values ​​to simulate the lower bound of the extreme deviation situation.

Condition 4: Except for the velocity deviation at the control-start point set to $-{\sigma }_V$, the deviations of other uncertain factors are still set to the corresponding $-3\sigma$ values. Condition 4 serves ​​as an additional control group for Condition 3.

Taking Scenario 3 in Section 3.3 as an example, using the closed-loop guidance law based on the PID controller designed in this paper, the tracking effects of the gliding projectile on the planned trajectory under different conditions are shown in Figure 7 (Ref-Trj in the legend represents the reference trajectory, Act-Trj represents the actual trajectory, and the suffix C is used to distinguish conditions).

The results from Figure 7 indicate that under different conditions, the projectile has successfully tracked the planned trajectory in the lateral plane, while the tracking situation in the longitudinal plane is relatively complex. Specifically, the tracking performance of the projectile trajectory is marginal under Condition 1, good under Condition 2, and under both Condition 3 and Condition 4 the projectile fails to track the planned trajectory. It can be inferred that based on the PID controller designed in this paper, the lateral tracking of the gliding projectile is easier compared to altitude tracking, the degree of deviation of uncertain factors leads to different effects on trajectory tracking. To further compare and analyze the influence of uncertain factor deviations on closed-loop guidance, the tracking performance of the planned trajectory under different conditions in different scenarios is shown in Table 5 (where ✓ indicates good tracking performance, × indicates tracking failure, and △ indicates marginal tracking performance).

The analysis in Section 3.3 reveals that the robustness of the planned trajectory increases from Scenario 1 to Scenario 3, with both the trajectory curvature and projectile flight time increasing sequentially, while the control margin of the projectile decreases accordingly. Therefore, based on the results in Table 5, it can be inferred that within a reasonable deviation range, the closed-loop guidance law designed in this study is most sensitive to changes in the projectile velocity deviation at the control-start point. When the velocity value is too large (Condition 1), for Scenarios 1 and 2 with relatively high control margins, the projectile can accurately track the planned trajectory through proper control. However, for Scenario 3 with insufficient control margin, the tracking performance is marginal. When velocity deviation at the control-start point is reasonable and other parameter deviations are similar to real-world conditions (Condition 2), the projectile in different scenarios can precisely track the planned trajectory. When the value of velocity is too small (Condition 3), the effective range of the projectile decreases, making it impossible to hit the target or track the planned trajectory. When the velocity is reduced by a normal amount but other parameter deviations are significant (Condition 4), the tracking performance of the projectile on the desired trajectory is negatively correlated with the control margin. In summary, to achieve good closed-loop guidance effects, the projectile velocity at the control-start point should not be lower than the design value of the reference trajectory, nor should it deviate significantly from the design value. Therefore, to alleviate the design pressure on the guidance control system, when planning the full trajectory of a gliding projectile, it is a good choice to slightly increase the design value of the projectile velocity at the control-start point.

4. Conclusions

This paper addresses the trajectory planning problem of gliding-guided projectiles with low maneuverability attacking fixed target. In order to reduce the sensitivity of the planned trajectory to uncertainties, enhance the compatibility between the planned trajectory and the guidance control system, and fully exploit its limited control capability, a closed-loop robust trajectory planning method based on the NIPCE method and PID controller is proposed. Numerical simulations are conducted using direct shooting method. The simulation results indicate that:

• When quantifying uncertainty propagation, compared to the traditional MCS method, the NIPCE-based method in this paper significantly enhances computational efficiency while ensuring accuracy. On the MATLAB simulation platform, FOPCE under single-core computation reduces the time by 79.5% compared to MCS under multi-core parallel computation. By removing unnecessary high-order cross terms using the basis truncation strategy, COPCE reduces the problem size by 83.3% and decreases computation time by 84.2%, facilitating robust planning.
• Open-loop robust planning can effectively reduce the sensitivity of gliding projectile trajectories to uncertainties. However, due to the limited control capability of the gliding projectile and the lack of closed-loop feedback from a guidance control system, even with minimized objective function, open-loop robust planning cannot eliminate terminal dispersion. Increasing curvature in the middle section of the trajectory improves the robustness of the planned trajectory but consumes additional control effort. Blindly pursuing robust optimality can lead to projectile control saturation, which is detrimental to the compatibility between the planned trajectory and the guidance control system.
• For gliding-guided projectiles, a trade-off between robust optimality and control effort optimality is necessary. Closed-loop robust planning considers the impact of uncertainties and the coupling between the planned trajectory and the guidance control system. It enhances trajectory robustness while effectively reducing control effort consumption.
• Within a reasonable deviation range, based on the PID controller designed in this paper, the closed-loop guidance law is most sensitive to changes in projectile velocity at the control-start point. To achieve optimal trajectory tracking, the projectile velocity at the control-start point should not be lower than the design value and should not deviate significantly from it.

The closed-loop robust trajectory planning method proposed in this paper considers the compatibility between the planned trajectory and the guidance control system on the basis of existing open-loop planning methods. This method helps to fully exploit the limited control capability of the gliding projectile, achieve maximize comprehensive performance, and can provide a reference for related research and engineering applications.