Reentry Tracking Guidance Method for Aircraft Based on Reinforcement Learning Algorithm

Abstract

The present invention discloses a reentry tracking guidance method for an aircraft based on a reinforcement learning algorithm, which is used to solve the problems of high dependence on the system model, poor adaptability and guidance performance in the prior art. The implementation scheme is as follows: establishing the continuous optimal control problem of aircraft reentry in the semi-speed coordinate system; transforming the continuous optimal control problem of aircraft reentry into a sequential convex optimal control problem; transforming the sequential convex optimal control problem into a sequential second-order cone programming problem; solving the sequential second-order cone programming problem to obtain the optimal bank angle of the aircraft; sampling the optimal bank angle trajectory to obtain a reference trajectory training data set; constructing a neural network and a reward function; using the training data set to perform offline training on the neural network until the cumulative reward converges to the maximum value to obtain a guidance network; obtaining the reentry guidance command of the aircraft online through the guidance network and tracking the target altitude in real time. The present invention has strong adaptability and high guidance accuracy, and can be used for rocket recovery.

Description

Technical Field

[0001] The present invention belongs to the technical field of aircraft guidance and control, and relates to a reentry tracking guidance method for an aircraft, which can be used for rocket recovery. Background Art

[0002] As the core technology of the aircraft overall and control system, trajectory and guidance design is crucial for an aircraft, and the two are gradually showing a trend of integration. The research on trajectory optimization technology early came from top-level design requirements such as coordinating overall parameter design, maximizing resource utilization rate, and reducing mission operation costs. The main task of traditional guidance design is to ensure the completion of the flight mission by offline designing a fixed reference trajectory and online tracking. Since the traditional tracking guidance scheme is based on a simplified motion model, there are problems of strong dependence on the model and poor adaptability. When the modeling deviation and external disturbance are large, the guidance performance will deteriorate severely. For example, the nominal trajectory guidance scheme during the reentry process of the "Apollo" spacecraft. The reentry flight environment of this scheme's aircraft is very complex. It has to go through outer space, the rarefied atmosphere, and the dense atmosphere. During this process, mechanical energy is rapidly converted into heat energy, the aircraft body will heat up violently, and at the same time, it faces the severe test of huge dynamic pressure and overload. Moreover, due to the characteristics of multi-constraints, strong coupling, fast time-variation, strong non-linearity, strong uncertainty, large envelope, and long flight time during the reentry process, it brings great challenges to the design of the guidance and control system, greatly restricting the flight performance that the aircraft should have. Summary of the Invention

[0003] The purpose of the present invention is to propose a reentry tracking guidance method for an aircraft based on a reinforcement learning algorithm in view of the defects existing in the above-mentioned prior art, so as to quickly generate an angle of attack increment command according to the real-time flight state, effectively enhance the task execution ability of the aircraft by tracking the desired altitude curve, reduce the dependence on the aircraft model, and improve the adaptability and guidance performance.

[0004] To achieve the above purpose, the technical solution of the present invention includes the following:

[0005] (1) Describe the trajectory optimization of the aircraft reentry section as a continuous optimal control problem P0 composed of a mathematical model, boundary conditions, admissible control, performance index, and process constraints;

[0006] (2) Perform convexification processing on P0 by means of form replacement, slack variable, softening constraint, and successive linearization methods to obtain a sequence of convex optimal control problems P1, and use the pseudospectral method to perform discrete parameterization processing on this P1 to obtain a sequence of second-order cone programming problems P2;

[0007] (3) Solve the sequence of second-order cone programming problems P2 by the interior point method to obtain the optimal reference trajectory;

[0008] (4) Sample the optimal reference trajectory to obtain a reference trajectory training data set;

[0009] (5) Construct a neural network, the Actor network:

[0010] (5a) Establish an action evaluation sub-network Actor_eval composed of a first input layer, a first hidden layer, and a first output layer connected in sequence. The state variables input to the first input layer are the five state variables of the geocentric distance, longitude, latitude, track angle, and heading angle of the current state agent. The information output by the first output layer is the angle of attack increment command for the current state;

[0011] (5b) Establish an action target sub-network Actor_target composed of a second input layer, a second hidden layer, and a second output layer connected in sequence. The state variables input to the second input layer are the five state variables of the geocentric distance, longitude, latitude, track angle, and heading angle of the next state of the agent. The information output by the second output layer is the target angle of attack increment command for the next state;

[0012] (5c) Parallelly connect the action evaluation sub-network Actor_eval network and the action target sub-network Actor_target network to form the Actor network, which is used to receive the state information of the experience replay pool and output the angle of attack increment command information;

[0013] (6) Construct a neural network, the Critic network:

[0014] (6a) Establish a value evaluation sub-network Critic_eval composed of a third input layer, a third hidden layer, and a third output layer connected in sequence. The variables input to the third input layer are the five state variables of the geocentric distance, longitude, latitude, track angle, and heading angle of the current state agent and the angle of attack increment command output by Actor_eval. The information output by the third output layer is the cumulative reward generated after the agent takes the command in this state;

[0015] (6b) Establish a value target sub-network Critic_target composed of a fourth input layer, a fourth hidden layer, and a fourth output layer connected in sequence. The variables input to the fourth input layer are the five state variables of the geocentric distance, longitude, latitude, track angle, and heading angle of the next state of the agent and the angle of attack increment command output by Actor_target. The information output by the fourth output layer is the target cumulative reward generated after the agent takes the command in the next state;

[0016] (6c) The value evaluation sub-network Critic_eval network and the value target sub-network Critic_target network are connected in parallel to form the Critic network, which is used to receive the state information of the experience replay pool and the angle of attack increment command information output by the Actor network, and output the cumulative reward information generated after the agent takes the command;

[0017] (7) The reward functions R of the Actor network and the Critic network are both designed as:

[0018] R = R 1 + R 2 , where R 1 is the target reward, and R 2 is the constraint reward;

[0019] (8) Use the reference trajectory training data set obtained in (4) to train the Actor network and the Critic network in parallel offline using the Deep Deterministic Policy Gradient algorithm DDPG. When the cumulative reward of the agent converges to a maximum value, the trained guidance network is obtained;

[0020] (9) Online obtain and real-time track the guidance command of the aircraft during reentry:

[0021] During the flight of the aircraft, the on-board computer reads the real-time flight state variables measured by the aircraft navigation system, uses the real-time flight state variables as the input of the guidance network for forward propagation, and obtains the real-time guidance command; the aircraft generates an angle of attack change according to the command to real-time track the predetermined altitude curve.

[0022] Compared with the prior art, the present invention has the following advantages:

[0023] 1) Since the present invention constructs a neural network model and transfers a large amount of optimization calculations to the offline training process, the online calculation amount is reduced, the speed of generating the guidance command is increased, and the aircraft can respond to the environment in real time.

[0024] 2) Since the present invention designs a reward function for the target altitude and the angle of attack increment, the guidance accuracy is improved, and the aircraft can adapt to complex environments and achieve fine tracking of the altitude. Brief Description of the Drawings

[0025] Figure 1 is the implementation flowchart of the present invention;

[0026] Figure 2 is the original non-convex set used in constructing the optimal trajectory of the aircraft in the present invention;

[0027] Figure 3 is the convex admissible control set obtained in convexifying the optimal trajectory of the aircraft in the present invention;

[0028] Figure 4 It is the reference trajectory of the bank angle in the reentry section obtained by solving the optimal trajectory of the aircraft in the present invention;

[0029] Figure 5 It is the cumulative benefit curve graph of the positive fixed offset of the aerodynamic parameters by using the present invention;

[0030] Figure 6 It is the tracking effect diagram of the positive fixed offset of the aerodynamic parameters by using the present invention;

[0031] Figure 7 It is the tracking effect diagram of the positive random offset of the aerodynamic parameters by using the present invention;

[0032] Figure 8 It is the reward curve graph when the positive random offset of the aerodynamic parameters is carried out by using the present invention;

[0033] Figure 9 It is the cumulative benefit curve graph of the negative fixed offset of the aerodynamic parameters by using the present invention;

[0034] Figure 10 It is the tracking effect diagram of the negative fixed offset of the aerodynamic parameters by using the present invention;

[0035] Figure 11 It is the tracking effect diagram of the negative random offset of the aerodynamic parameters by using the present invention;

[0036] Figure 12 It is the reward curve graph of the negative random offset of the aerodynamic parameters by using the present invention. Specific embodiments

[0037] The following further details the specific embodiments and effects of the present invention in conjunction with the accompanying drawings.

[0038] Refer to Figure 1 , the implementation steps of this example are as follows:

[0039] Step 1, establish the aircraft motion model.

[0040] The motion equation of the aircraft centroid is the basis for studying its motion characteristics. In this example, the reusable launch vehicle RLV without power is taken as the research object, and its reentry centroid motion equation is established in the semi-velocity system. The trajectory optimization problem usually does not consider the action of the control force, and since the additional Coriolis force level is relatively small, its influence is ignored. The RLV reentry process mainly relies on aerodynamic force and earth gravity to change the flight state. Assuming that the earth is an ideal sphere, in the semi-velocity coordinate system, the state differential equation f1 of the RLV is established as follows:

[0041]

[0042] where r is the dimensionless geocentric distance of the RLV, is the dimensionless velocity of the RLV, Ω is the angular velocity of the Earth's rotation, R0 is the radius of the Earth, γ is the flight path angle, ψ is the heading angle, σ is the bank angle, L = 0.5R0ρV 2 SC L / m is the lift acceleration, D = 0.5R0ρV 2 SC D / m is the drag acceleration; e is the dimensionless energy of the RLV, C D and C L are the drag coefficient and the lift coefficient respectively. Both of these coefficients are functions of the angle of attack α and the Mach number Ma. S is the reference area of the RLV, ρ = ρ0exp(-(r - 1) / H) represents the atmospheric density, which is a function of the geocentric distance r. ρ0 is the atmospheric density at sea level, and H is the dimensionless scale height.

[0043] For the state differential equation, in this embodiment, the state variables are selected as s = (r, θ, φ, γ, ψ), and the control variable is the bank angle u = σ.

[0044] Step 2, establish the aircraft constraint conditions and performance indicators.

[0045] 2.1) Establish the aircraft constraint conditions:

[0046] 2.1.1) Establish the process constraint K:

[0047] Since a large amount of heat is generated by the friction between the airframe and the atmosphere when the reusable launch vehicle without power (RLV) flies in the atmosphere, considering the structural strength of the RLV and the safe working conditions of the equipment, for flight safety, constraints are imposed on the dynamic pressure q, the heat flux rate and the normal overload n, denoted as K:

[0048]

[0049] where k Q is the heat flux rate calculation coefficient, g0 is the sea-level gravitational acceleration, q max , and n max are the maximum dynamic pressure, the maximum heat flux rate, and the maximum overload that the RLV can withstand respectively;

[0050] 2.1.2) Establish the initial state constraint L and the terminal state constraint Φ1:

[0051] According to the flight mission of the RLV, the initial state information and the terminal state information of the RLV can be determined, and then constraints are imposed on the initial state and the terminal state of the RLV:

[0052] Initial state constraint L:

[0053]

[0054] Among them, \(x\) represents five state variables of the aircraft, represents the dimensionless initial state of the aircraft, and \(e_0\) represents the initial dimensionless energy of the RLV;

[0055] Terminal state constraint \(\varPhi_1\):

[0056]

[0057] Among them, \(e\) f represents the terminal dimensionless energy of the RLV, respectively represent the set values of the terminal dimensionless geocentric distance, longitude, latitude, flight path angle and heading angle of the RLV;

[0058] 2.1.3) Establish the admissible control \(C\):

[0059] During the flight of the RLV, in order to ensure that the RLV will not perform large-angle flips, the control variable bank angle should also be restricted. Therefore, an admissible control constraint \(C\) is established to restrict the amplitude of the bank angle;

[0060] Admissible control \(C\):

[0061] \(C\): \(\sigma\) min \(\leq|\sigma|\leq\sigma\) max

[0062] Among them, \(\sigma\) min and \(\sigma\) max respectively represent the lower and upper limits of the bank angle.

[0063] 2.2) Establish the performance index \(J\)

[0064] In the trajectory optimization of the aircraft reentry, the usual expected requirements for the aircraft include the shortest reentry time, the maximum terminal velocity, the longest flight range, and the minimum energy consumption. In this example, the performance index \(J\) is the shortest reentry time, that is:

[0065]

[0066] Step 3, construct the optimal control problem \(P_0\) of the RLV reentry trajectory optimization model.

[0067] According to the above state differential equation \(f_1\), process constraint \(K\), initial constraint \(L\), terminal constraint \(\varPhi_1\), admissible control \(C\) and performance index \(J\), the trajectory optimization model of the RLV reentry section is described as an optimal control problem \(P_0\) including these parameters:

[0068]

[0069] Step 4: Transform the reentry continuous optimal control problem P0 of the aircraft into a sequential convex optimal control problem P1.

[0070] 4.1) Form transformation of process constraint K:

[0071] Redescribe the process constraint K as a linear inequality form with respect to the state variable r, meeting the requirements of a convex problem:

[0072]

[0073] The above formula is equivalent to

[0074]

[0075] where l Q (e), l q (e), l n (e) represent the functions of the heat flux rate, dynamic pressure, and overload of the aircraft with respect to the dimensionless geocentric distance r respectively. denotes the maximum value among l Q (e), l q (e), l n (e);.

[0076] 4.2) Relaxation of control constraints:

[0077] 4.2.1) Transform control variables

[0078] Since choosing the bank angle σ as the control variable is prone to produce bad oscillations during the successive linearization of the differential equation system, it is necessary to replace the control variable to eliminate such bad oscillations, which is described in detail as follows:

[0079] First, according to the non-convex property of the centroid motion equation f1 with respect to both the state variable x and the control variable u = σ when choosing the bank angle σ as the control variable, use the successive linearization method to convexify the equation f1. The solution at the k-th step during the successive solution process includes the state variable solution x (k) (e) and the control variable solution u (k) (e), denoted as {x (k) (e); u (k) (e)};

[0080] Second, according to the solution obtained in the k-th iteration, linearize the equation in the (k + 1)-th iteration as:

[0081] x′ = F x (x (k) , u (k) , e)x + F u (x (k) , u (k) , e)u + b(x (k) , u(k) , e)

[0082] where x′ represents the derivative of the five state variables with respect to the dimensionless energy e, represents the Jacobian matrix of the differential equation system f with respect to x, represents the Jacobian matrix of the differential equation system f with respect to u, and b(x (k) , u (k) , e) = f(x (k) , u (k) , e) - F x (x (k) , u (k) , e)x (k) - F u (x (k) , u (k) , e)u (k) represents the co-vector.

[0083] Next, based on the above coefficient matrices F x , F u and the co-vector b related to u (k) , in the iterative process, the control variable u (k) is prone to oscillation, that is, in the iterative process, the oscillation of the control variable u (k) is used to control the control variable u (k+1) obtained in the (k + 1)-th iteration. As the iteration progresses, the oscillation may be further amplified. These oscillations are not inherent in the solution of the problem and are manifestations of the instability of the numerical solution process. To eliminate the relationship between the linearized dynamic equation and u (k) , the following control variable is introduced to replace the original control variable: where u1 represents the cosine value of the tilt angle and u2 represents the sine value of the tilt angle. The replaced control variable is separated from the centroid motion equation, and the differential equation system is linearized based on the new control variable. Its Jacobian coefficient matrix is approximately a 0 matrix, and the adverse oscillation transmission can be eliminated.

[0084] 4.2.2) Establish new control variable constraints

[0085] After the control variable transformation, the control constraint is expressed as: cosσ max ≤ u1 ≤ cosσ min . Since the control variables u1 and u2 are not independent of each other, the following equation also needs to be satisfied:

[0086]

[0087] This equality constraint is a quadratic constraint, and the original admissible control set determined by the constraint is as Figure 2 shown. This original admissible control set is on the unit circle at u1 = ω l = cosσmax and \(u_1 = \omega\) h \(= \cos\sigma\) min The two arcs between them form a non - convex set, and convexification is needed, that is, the constraints are relaxed, and the original constraints are relaxed into the following second - order cone inequality constraints:

[0088]

[0089] The set determined by relaxing the constraints expands the original admissible control set, making the original non - convex set become a convex set, as shown in the shaded area of Figure 3.

[0090] 4.3) Successive linearization of the performance index:

[0091] 4.3.1) Substitute the drag acceleration expression into the performance index and express the performance index as a function of the following state variables:

[0092]

[0093] where \(D = 0.5R_0\rho V\) 2 \(S_C\) D / m is the drag acceleration, and the coefficient does not contain any state variables, is the non - dimensional velocity of the RLV, m is the mass of the aircraft, \(\Omega\) is the angular velocity of the Earth's rotation, \(R_0\) is the radius of the Earth, e is the non - dimensional energy of the RLV, \(C\) D is the drag coefficient, which is a function of the angle of attack \(\alpha\) and the Mach number Ma, S is the reference area of the RLV, \(\rho=\rho_0\exp(-(r - 1) / H)\) represents the atmospheric density, which is a function of the geocentric distance r, \(\rho_0\) is the atmospheric density at sea level, and H is the non - dimensional scale height.

[0094] 4.3.2) Linearize the performance index

[0095] The atmospheric density \(\rho\) in the performance index formula is the only quantity that can be optimized, and the related state component is the geocentric distance r. The atmospheric density approximately adopts an exponential model. The integrand of this performance index is a convex function of r, but to finally obtain a second - order cone programming problem, it needs to be transformed into a linear function of r. Therefore, near the solution \(r\) (k) (e) of the previous iteration, perform a first - order Taylor expansion on the integrand of this performance index:

[0096]

[0097] where \(\rho\) (k) represents the atmospheric density calculated using \(r\) (k) (e), represents the value of the derivative of the atmospheric density \(\rho\) with respect to r at \(r\) (k) (e). respectively represent the solution r in the previous iteration (k) (e) The first-order term and the remainder of the performance index near the dimensionless geocentric distance r;

[0098] 4.3.3) Introduce a regularization term to ensure the constraints remain unchanged

[0099] After linearizing the original performance index, a regularization term needs to be introduced to ensure that the relaxation of the constraints in step 4.2) has a positive effect, that is, introduce it into the linearized performance index

[0100] where ε ψ ≠0 is a constant. To ensure that the obtained solution is approximate to the solution of the original problem, the magnitude of ε ψ needs to be set small enough.

[0101] 4.4) Relaxation of terminal constraints:

[0102] Although the terminal constraints are all linear equality constraints, in the initial optimization process, it may be difficult to satisfy the two hard constraints of longitude and latitude. Therefore, the following penalty terms are introduced into the performance index to replace these two terminal constraints:

[0103]

[0104] where d(θ(e f ), φ(e f )) represents the penalty term, c θ >0 and c φ >0 are given constants,

[0105] Since this penalty term is non-linear, two slack variables and need to be introduced to transform it into a linear function and incorporate it into the performance index: and subject to the linear inequality constraints:

[0106] Update the performance index to:

[0107]

[0108] where c1 and c0 respectively represent the first-order term and the remainder of the performance index with respect to the dimensionless geocentric distance r during successive linearization, respectively represent the set values of the dimensionless geocentric distance, longitude, latitude, flight path angle, and heading angle at the RLV terminal.

[0109] 4.5) Successive linearization of the centroid motion equation:

[0110] According to the replaced control quantity, the centroid motion equation can be expressed as:

[0111]

[0112] Among them, f1(x, e) is the coefficient matrix of the terms without the control quantity in the original differential equation system, and B(x, e) is the coefficient matrix of the terms with the control quantity in the original differential equation system;

[0113] The successive linearization method is used to convexify the centroid motion equation. In the successive solution process, the solution at the k-th time includes the state variable solution x (k) (e) and the control quantity solution u (k) (e), denoted as {x (k) (e); u (k) (e)};

[0114] According to the solution obtained in the k-th iteration, the equation is linearized for the (k + 1)-th iteration as:

[0115] x′ = F x (x (k) , u (k) , e)x + B(x (k) , e)u + b(x (k) , u (k) , e)

[0116] Among them, represents the Jacobian matrix;

[0117] b(x (k ), u (k) , e) = f1(x (k) , e) - F x (x (k) , u (k) , e)x (k) represents the co-vector;

[0118] represents the partial derivative of B(x, e)u with respect to x, represents the partial derivative of f1(x, e) with respect to x;

[0119] L / D = C L / C D , both of which can be regarded as single-variable functions of the energy e. Therefore, the only element in the matrix B(x, e) related to the state variable x is cosγ, resulting in having a non-zero term. Since during the reentry process of the aircraft, the flight path angle |γ| is very small, that is, this non-zero term is very small. Therefore, the linearized equation can be simplified as:

[0120] x′ = F x (x (k) , e)x + B(x(k) , e) u + b (x (k) , e);

[0121] Furthermore, the terms related to the Earth's rotation with very small magnitudes are separated from the centroid motion equation, and the centroid motion equation is updated to the following concise form:

[0122] x′ = f1(x, e) + B(x, e)u = f0(x, e) + f Ω (x, e) + B(x, e)u

[0123] where, represents the terms in the differential equation without the control quantity, and f Ω (x, e) is the term related to the Earth's angular velocity;

[0124] Furthermore, the updated centroid motion equation is linearly approximated successively, and the dimensionless energy e in the symbols is omitted. According to the solution obtained from the k-th iteration, the equation for the (k + 1)-th iteration is linearized as:

[0125] x′ = A(x (k) ) x + B(x (k) ) u + b(x (k) )

[0126] where, b(x (k) ) = f0(x (k) ) - A(x (k) ) x (k) + f Ω (x (k) ) represents the co-vector, A(x (k) ) and B(x (k) ) represent the Jacobian matrices of f0(x (k) ) with respect to x and u respectively, f0(x (k) ) represents the state differential equation obtained from the k-th iteration of f0(x, e), and f Ω (x (k) ) represents the state differential equation obtained from the k-th iteration of f Ω (x, e);

[0127] To ensure the local effectiveness of successive linearization, a trust region constraint is imposed: |x - x (k) | ≤ δ, where δ is a given five-dimensional constant vector.

[0128] 4.6) According to the convexification processing of the constraints and performance indices in steps 4.1) to 4.5), the trajectory optimization model for the reentry phase of the RLV is re-described as a sequential convex optimal control problem P1:

[0129]

[0130] Step 5: Convert the sequential convex optimal control problem P1 into a sequential second-order cone programming problem P2.

[0131] According to the principle of the Gauss pseudospectral method, the sequential convex optimal control problem P1 is discretized into a sequential convex programming problem P2, and the implementation steps are as follows:

[0132] 5.1) Set collocation points and discrete points:

[0133] Convert the energy independent variable e to the pseudo-energy independent variable Ε ∈ [-1, 1]:

[0134]

[0135] where e0 and e f respectively represent the initial dimensionless energy and the terminal dimensionless energy of the aircraft;

[0136] Set N LG points (Ε1, Ε2, …, Ε m , …, Ε N ) in Ε ∈ (-1, 1) as collocation points, and the coordinates of each collocation point are the roots of each Nth-order Legendre polynomial , where m = 1, 2, …, N represents each LG point;

[0137] Set discrete points, including the two endpoints Ε0 and Ε f of the pseudo-energy independent variable and all LG points inside the intervals, that is, a total of N + 2 discrete points Ε0, Ε1, …, Ε n , …, Ε N+1 , where n = 0, 1, …, N + 1 represents each discrete point, and Ε N+1 = Ε f represents the right endpoint of the pseudo-energy independent variable.

[0138] 5.2) Discretization of the differential equation of the center-of-mass motion:

[0139] 5.2.1) Use the Lagrange interpolation polynomial to represent the constraint of the differential equation of the center-of-mass motion for the pseudo-energy independent variable Ε ∈ [-1, 1):

[0140] Use N + 1 Lagrange interpolation polynomials as basis functions to approximate the state variables at each LG point where i = 0, 1, …, N represents each interpolation polynomial;

[0141] Use the approximate state variable x(Ε m ) to differentiate the pseudo-energy independent variable Ε , that is where, represents the basis function L i(Ε m ) The derivative at the LG point;

[0142] Using to replace x′ in the differential equation of the centroid motion, the constraint of the differential equation of the centroid motion is transformed into an algebraic equation constraint at the LG collocation points:

[0143]

[0144] where represents the solution of the state variable at the k - th iteration at the m - th LG point, x m represents the solution of the state variable at the m - th LG point, u m represents the solution of the control quantity at the m - th LG point;

[0145] 5.2.2) Using Gauss integration to estimate the constraint of the differential equation of the centroid motion at the right - hand endpoint of the pseudo - energy independent variable:

[0146] Since the pseudo - energy interval corresponding to the approximate expression of the state variable is [-1, 1), which does not include the terminal state variable x N+1 , so Gauss integration is used to estimate x N+1 , that is, the constraint of the differential equation of the centroid motion that the terminal state x N+1 needs to satisfy:

[0147]

[0148] where w m represents the LG integration weight.

[0149] 5.3) Set other discrete constraint conditions:

[0150] 5.3.1) Trust - region constraint:

[0151]

[0152] where represents the solution of the state variable at the k - th iteration at the n - th discrete point;

[0153] 5.3.2) Control constraint:

[0154]

[0155] where, (u1) m and (u2) m respectively represent the values of the control quantities u1 and u2 at the m - th LG point, ω l and ω h respectively represent the minimum and maximum values of the cosine of the bank angle;

[0156] 5.3.3) Process constraint:

[0157]

[0158] Among them, r m represents the value of the non - dimensional geocentric distance of the aircraft at the m - th LG point, represents the maximum values of the heat flux rate, dynamic pressure, and overload of the aircraft at the m - th LG point;

[0159] 5.3.4) Boundary conditions:

[0160]

[0161]

[0162]

[0163]

[0164]

[0165]

[0166] Among them, represents the initial value set for the state variables of the aircraft, r N+1 represents the terminal state of the non - dimensional geocentric distance of the aircraft, θ N+1 represents the terminal state of the longitude of the aircraft, φ N+1 represents the terminal state of the latitude of the aircraft, γ N+1 represents the terminal state of the flight path angle of the aircraft, ψ N+1 represents the terminal state of the heading angle of the aircraft, respectively represent the terminal non - dimensional geocentric distance, longitude, latitude, flight path angle, and heading angle of the aircraft, and respectively represent the relaxation variables for longitude and latitude;

[0167] 5.3.5) Performance index:

[0168] The performance index function is approximated by Gauss integration as:

[0169]

[0170] Among them, c θ 、 respectively represent the penalty coefficients for longitude and latitude, respectively represent the first - order term and the remainder term with respect to the non - dimensional geocentric distance r in the successive linearization process of the performance index, ε ψ represents the coefficient of the heading angle regularization term, ψ m represents the value of the heading angle at the m - th LG point.

[0171] 5.4) According to the discretization of each constraint and performance index in steps 5.1) to 5.3), the trajectory optimization model of the RLV reentry section is re-described as a sequential convex optimal control problem P2:

[0172]

[0173] where n = 0, 1, 2... N + 1 represents each discrete point, and m = 1, 2, 3... N represents each LG point.

[0174] Step 6, perform offline solution for the sequential second-order cone programming problem P2 of the unpowered reusable launch vehicle (RLV) reentry section. The process is as follows:

[0175] 6.1) Let the iteration count k = 0, and given the initial state variables at the discrete points represents the solution of the state variables at the nth discrete point in the kth iteration, n = 0, 1, 2... N + 1;

[0176] 6.2) Increment the iteration count k by one, and use the interior point method to solve problem P2 to obtain the solution where z (k+1) represents the solution vector of the (k + 1)th iteration,

[0177] represents the solution of the state variables of the (k + 1)th iteration, represents the solution of the state variables at the nth discrete point in the (k + 1)th iteration, n = 0, 1, 2... N + 1,

[0178] represents the solution of the control variables of the (k + 1)th iteration, represents the solution of the control variables at the mth discrete point in the (k + 1)th iteration, m = 1, 2, 3... N,

[0179] and represents the solution of the slack variables of the (k + 1)th iteration;

[0180] 6.3) Let a five-dimensional column vector ε represent the convergence judgment condition, and judge whether the solution of the state variables at each discrete point satisfies the convergence condition:

[0181] If not satisfied, let x (k) = x (k+1) , and return to step 6.2;

[0182] Otherwise, output the optimal solution The calculation ends.

[0183] Finally, according to the u in the optimal solution (k+1) obtain the reference trajectory of the bank angle of the reentry section control quantity, as Figure 4 shown.

[0184] Step 7: Sample the optimal reference trajectory obtained in Step 6 to obtain a reference trajectory training data set.

[0185] Step 8: Construct a neural network Actor and a neural network Critic.

[0186] 8.1) Construct the Actor network:

[0187] 8.1.1) Establish an action evaluation sub-network Actor_eval composed of a first input layer, a first hidden layer, and a first output layer connected in sequence. The state variables input to the first input layer are the five state variables of the geocentric distance, longitude, latitude, track angle, and heading angle of the current state agent. The information output by the first output layer is the angle of attack increment command of the current state;

[0188] 8.1.2) Establish an action target sub-network Actor_target composed of a second input layer, a second hidden layer, and a second output layer connected in sequence. The state variables input to the second input layer are the five state variables of the geocentric distance, longitude, latitude, track angle, and heading angle of the next state of the agent. The information output by the second output layer is the target angle of attack increment command of the next state;

[0189] 8.1.3) Connect the action evaluation sub-network Actor_eval network and the action target sub-network Actor_target network in parallel to form the Actor network, which is used to receive the state information of the experience replay pool and output the angle of attack increment command information;

[0190] 8.2) Construct the Critic network:

[0191] 8.2.1) Establish a value evaluation sub-network Critic_eval composed of a third input layer, a third hidden layer, and a third output layer connected in sequence. The variables input to the third input layer are the five state variables of the geocentric distance, longitude, latitude, track angle, and heading angle of the current state agent and the angle of attack increment command output by the action evaluation sub-network Actor_eval. The information output by the third output layer is the cumulative reward generated after the agent takes the command in this state;

[0192] 8.2.2) Establish a value target sub-network Critic_target composed of a fourth input layer, a fourth hidden layer, and a fourth output layer connected in sequence. The variables input to the fourth input layer are the five state variables of the geocentric distance, longitude, latitude, track angle, and heading angle of the next state agent, as well as the angle of attack increment command output by the action target sub-network Actor_target. The information output by the fourth output layer is the target cumulative reward generated after the agent takes the command in the next state;

[0193] 8.2.3) Connect the value evaluation sub-network Critic_eval network and the value target sub-network Critic_target network in parallel to form the Critic network, which is used to receive the state information of the experience replay pool and the angle of attack increment command information output by the Actor network, and output the cumulative reward information generated after the agent takes the command;

[0194] 8.3) Set all hidden layers in the Actor network and the Critic network to 100 neurons, use the relu function for the activation functions of the hidden layers, and use the tanh function for the output layers.

[0195] Step 9, set the reward functions of the neural network Actor and the neural network Critic.

[0196] The reward functions include target rewards and constraint rewards, and are set as follows:

[0197] 9.1) According to the task objective, set the target reward R 1 :

[0198]

[0199] Among them, let h next represent the height at the next time step after executing the action Δα, Δα represents the angle of attack increment command, and h0 represents the nominal height corresponding to the next time step;

[0200] 9.2) According to the magnitude of the angle of attack increment, set the constraint reward R 2 :

[0201] R 2 =-10|Δα|

[0202] 9.3) Add the target reward R 1 and the constraint reward R 2 to obtain the reward functions of the neural network Actor and the neural network Critic as: R = R 1 +R 2 .

[0203] Step 10: Synchronously train the neural network Actor and the neural network Critic offline to obtain a trained guidance network.

[0204] 10.1) Parameter initialization:

[0205] Randomly initialize the parameters θ of the Actor_eval network μ and the parameters θ of the Actor_target network μ′ ;

[0206] Randomly initialize the parameters θ of the Critic_eval network Q and the parameters θ of the Critic_target network Q′ ;

[0207] Initialize the relevant hyperparameters of the training process: that is, set the parameter for replacing the update frequency of the target network as τ, the number of episodes for the target network update period as T, the size of the experience replay pool as M, the total number of transition processes sampled from the experience replay pool at each time step as batch, and the reward discount rate as λ;

[0208] 10.2) Use the Markov decision process MDP model to describe the single-step state transition process of the agent:

[0209] 10.2.1) Describe the single-step transition process in combination with the reentry motion model:

[0210] According to the aircraft reentry motion model constructed in Step 1, let the environmental state of the agent be the trajectory state variable S = [r, θ, φ, γ, ψ], and the action of the agent be the angle of attack increment Δα;

[0211] In each state transition process, the agent first receives a current environmental state information S t = [r t , θ t , φ t , γ t , ψ t , and then makes a corresponding action A t = Δα t . After the model executes this action, according to the reference trajectory training dataset obtained in Step 7 and S t and A t , the current state S t is transferred to the next state S t+1 through Runge-Kutta integration. At the same time, the agent will receive a reward R t+1 , where t represents the current state of the agent, and t + 1 represents the next state the agent transfers to after executing the action;

[0212] 10.2.2) Add random noise to the agent's action and extract a new action to replace the original action:

[0213] To enable the agent to explore more actions and have stronger environmental adaptability, add random noise to the action A t = Δα t selected by the agent;

[0214] Set the parameter v ar as the standard deviation, and use the original action A t selected by the reinforcement learning algorithm as the mean value to construct a normal distribution function, and then randomly extract a new action a t from the normal distribution function to replace A t .

[0215] 10.3) Store the single-step transfer process in step 10.2 in the form of (S t , a t , R t+1 , S t+1 ) into the experience replay pool, where S t represents the current state, a t represents the angle of attack increment action taken at the current moment, R t+1 represents the reward value obtained after executing the action, and S t+1 represents the next state;

[0216] 10.4) Randomly extract some samples from the experience replay pool and input them into the Actor network and the Critic network. Update the parameters of the Actor_eval sub-network in the Actor network by maximizing the cumulative reward through the policy gradient of the performance metric , and update the parameters of the Critic_eval sub-network in the Critic network by minimizing the loss function , that is, the mean square error between the current cumulative reward and the target cumulative reward. Update the parameters of the Actor_target sub-network in the Actor network and the parameters of the Critic_target sub-network in the Critic network by weighting according to τ and the original network parameters every T rounds ,

[0217] where N represents the number of samples sampled from the experience replay pool, s i is the current state quantity, a i is the angle of attack increment command selected by the neural network according to the current state quantity, y i = R i+1 + λQ′(s i+1 , μ′(s i+1 |w μ′ )w Q′), representing the target cumulative reward, Q(s i , a i | w Q ) represents the output value of the Critic_eval network, w represents the weight and threshold parameters of the network, and the superscript of w indicates the network to which this parameter belongs. The superscript Q represents the Critic_eval network, the superscript Q′ represents the Critic_target network, the superscript μ represents the Actor_eval network, and the superscript μ′ represents the Actor_target network. μ′(s i+1 | w μ′ ) represents the output value of the Actor_target network; Q′(s i+1 , μ′(s i+1 | w μ′ ) | w Q′ ) represents the output value of the Critic_target network, and the output value of the Actor_eval network is μ(s i | w μ ).

[0218] 10.5) Repeat step 10.4). When the cumulative reward of the agent converges to a maximum value, the trained Actor network and Critic network are obtained, and the action evaluation sub-network Actor_eval network in the Actor network is used as the guidance network.

[0219] Step 11, obtain the reentry guidance command of the aircraft online.

[0220] During the flight of the aircraft, the on-board computer reads the real-time flight state variables measured by the aircraft navigation system, uses the real-time flight state variables as the input of the guidance network for forward propagation to obtain the real-time guidance command, and according to the command, the aircraft will generate an angle of attack change to track the predetermined altitude curve in real time.

[0221] The technical effects of this embodiment are further described below in combination with simulation experiments:

[0222] I. Simulation conditions

[0223] The initial conditions of the aircraft are:

[0224] Initial geocentric distance Initial longitude Initial latitude Initial track angle Initial heading angle Initial speed

[0225] The terminal conditions of the aircraft are:

[0226] Terminal geocentric distance Terminal longitude Terminal latitude Terminal track angle Terminal heading angle Terminal speed

[0227] All numerical simulations were carried out on a PC with an I7 - 8700F (3.20GHz) CPU, and the software environment was MATLAB and PYCHARM.

[0228] II. Simulation content

[0229] Simulation experiment 1:

[0230] Under the above - mentioned simulation conditions, the lift coefficient and drag coefficient at each discrete time step in the present invention were fixed at a 20% positive deviation, and the convergence of the present invention was tested. The results are as Figure 5 shown. It can be seen from Figure 5 that the cumulative reward curve of the present invention converges to a maximum value within 200 rounds, ensuring the convergence and data accuracy under the positive deviation condition;

[0231] Simulation experiment 2:

[0232] Under the above - mentioned simulation conditions, the lift coefficient and drag coefficient at each discrete time step in the present invention were fixed at a 20% positive deviation, and the tracking effect of the guidance network of the present invention on altitude under the positive deviation condition was tested. The results are as Figure 6 shown, where the dashed line is the optimal altitude trajectory calculated under the nominal aerodynamic parameter conditions, and the solid line is the altitude curve generated by the guidance network of the present invention. According to Figure 6 calculations, the terminal altitude guidance accuracy of the guidance network of the present invention under the positive deviation condition is 0.52%.

[0233] Simulation experiment 3:

[0234] Due to the complex and variable flight environment of the RLV, to simulate a more realistic flight environment, under the above - mentioned simulation conditions, the aerodynamic coefficients at each discrete time step in the present invention were randomly deviated within the range of [0, 20%], and the deviation amounts of the lift coefficient and drag coefficient were different. The tracking effect of the guidance network of the present invention on altitude under the complex environment of positive random deviation was tested. The results are as Figure 7 shown, where the dashed line is the optimal altitude trajectory calculated under the nominal aerodynamic parameter conditions, and the solid line is the altitude curve generated by the guidance network of the present invention. According to Figure 7 calculations, the terminal altitude guidance accuracy of the guidance network of the present invention under the complex environment of positive random deviation is 0.1%.

[0235] Simulation experiment 4:

[0236] Under the above simulation conditions, the aerodynamic coefficients at each discrete time step in the present invention are randomly deviated within the range of [0, 20%], and the deviation amounts of the lift coefficient and the drag coefficient are different. The difficulty of the aircraft in the present invention to track altitude in a complex environment with positive random deviation is tested. The results are as Figure 8 shown. As can be seen from Figure 8 , the reward curve of the present invention fluctuates between -1000 and -1150, and the reward at the end of the trajectory drops rapidly, indicating that it is more difficult for the aircraft in the present invention to track at the end of the trajectory in a complex situation with positive random deviation.

[0237] Simulation Experiment 5:

[0238] Under the above simulation conditions, the lift coefficient and the drag coefficient at each discrete time step in the present invention are fixed at a reverse deviation of 20%. The convergence of the present invention is tested. The results are as Figure 9 shown. As can be seen from Figure 9 , the cumulative reward curve of the present invention converges to a maximum value within 200 rounds, ensuring the convergence in the case of reverse deviation and the accuracy of the data.

[0239] Simulation Experiment 6:

[0240] Under the above simulation conditions, the lift coefficient and the drag coefficient at each discrete time step in the present invention are fixed at a reverse deviation of 20%. The tracking effect of the guidance network of the present invention on altitude in the case of reverse deviation is tested. The results are as Figure 10 shown, where the dashed line is the optimal altitude trajectory solved and calculated under the condition of nominal aerodynamic parameters, and the solid line is the altitude curve generated by the guidance network of the present invention. According to Figure 10 , it can be calculated that the terminal altitude guidance accuracy of the guidance network of the present invention in the case of reverse deviation is 2.9%.

[0241] Simulation Experiment 7:

[0242] Under the above simulation conditions, the aerodynamic coefficients at each discrete time step in the present invention are randomly deviated within the range of [-20%, 0], and the deviation amounts of the lift coefficient and the drag coefficient are different. The tracking effect of the guidance network of the present invention on altitude in a complex environment with reverse random deviation is tested. The results are as Figure 11 shown, where the dashed line is the optimal altitude trajectory solved and calculated under the condition of nominal aerodynamic parameters, and the solid line is the altitude curve generated by the guidance network of the present invention. According to Figure 11 , it can be calculated that the terminal altitude guidance accuracy of the guidance network of the present invention in a complex situation with reverse random deviation is 2.4%.

[0243] Simulation Experiment 8:

[0244] Under the above simulation conditions, the aerodynamic coefficients at each discrete time step in the present invention are randomly offset within the range of [-20%, 0], and the offset amounts of the lift coefficient and the drag coefficient are different. The difficulty of the aircraft in the present invention to track the altitude in a complex environment with reverse random offset is tested. The results are as Figure 12 shown. It can be seen from Figure 12 that the reward curve of the present invention fluctuates between -1000 and -1025, and the reward at the end of the trajectory slightly increases, indicating that the tracking difficulty of the aircraft at the end of the trajectory is relatively small in the complex situation of reverse random offset in the present invention.

[0245] It can be seen from the simulation experiment under the above random offset conditions that the solution proposed in the present invention has good tracking effect and high accuracy for altitude when the aircraft is in a complex environment, and the trained guidance network can better adapt to the deviation of aerodynamic parameters.

Claims

1. A reentry tracking guidance method for an aircraft based on a reinforcement learning algorithm, characterized in that, It includes the following steps: (1) Optimize the reentry trajectory of the aircraft and describe it as a continuous optimal control problem P0 composed of a mathematical model, boundary conditions, admissible control, performance indicators, and process constraints; (2) Perform convexification processing on P0 by means of form replacement, slack variables, softened constraints, and successive linearization methods to obtain a sequence of convex optimal control problems P1. Use the pseudospectral method to perform discrete parameterization processing on P1 to obtain a sequence of second-order cone programming problems P2; (3) Use the interior point method to solve the sequence of second-order cone programming problems P2 to obtain the optimal reference trajectory; (4) Sample the optimal reference trajectory to obtain a reference trajectory training data set; (5) Construct a neural network Actor network: (5a) Establish an action evaluation sub-network Actor_eval composed of a first input layer, a first hidden layer, and a first output layer connected in sequence. The state variables input by the first input layer are the five state variables of the geocentric distance, longitude, latitude, track angle, and heading angle of the current state agent. The information output by the first output layer is the angle of attack increment command for the current state; (5b) Establish an action target sub-network Actor_target composed of a second input layer, a second hidden layer, and a second output layer connected in sequence. The state variables input by the second input layer are the five state variables of the geocentric distance, longitude, latitude, track angle, and heading angle of the next state of the agent. The information output by the second output layer is the target angle of attack increment command for the next state; (5c) Connect the action evaluation sub-network Actor_eval network and the action target sub-network Actor_target network in parallel to form the Actor network, which is used to receive the state information of the experience replay pool and output the angle of attack increment command information; (6) Construct a neural network Critic network: (6a) Establish a value evaluation sub-network Critic_eval composed of a third input layer, a third hidden layer, and a third output layer connected in sequence. The variables input by the third input layer are the five state variables of the geocentric distance, longitude, latitude, track angle, and heading angle of the current state agent and the angle of attack increment command output by Actor_eval. The information output by the third output layer is the cumulative reward generated after the agent takes the command in this state; (6b) Establish a value target sub-network Critic_target composed of a fourth input layer, a fourth hidden layer, and a fourth output layer connected in sequence. The variables input by the fourth input layer are the five state variables of the geocentric distance, longitude, latitude, track angle, and heading angle of the next state of the agent and the angle of attack increment command output by Actor_target. The information output by the fourth output layer is the target cumulative reward generated after the agent takes the command in the next state; (6c) Connect the value evaluation sub-network Critic_eval network and the value target sub-network Critic_target network in parallel to form the Critic network, which is used to receive the state information of the experience replay pool and the angle of attack increment command information output by the Actor network, and output the cumulative reward information generated after the agent takes the command; (7) The reward functions R for both the Actor network and the Critic network are as follows: R = R 1 + R 2 , where R 1 is the target reward and R 2 is the constraint reward; (8) Use the reference trajectory training dataset obtained in (4) to train the Actor network and the Critic network in parallel offline using the Deep Deterministic Policy Gradient algorithm DDPG. When the cumulative reward of the agent converges to a maximum value, the trained guidance network is obtained. (9) Online obtain the reentry guidance command of the aircraft and track it in real time: During the flight of the aircraft, the on-board computer reads the real-time flight state variables measured by the aircraft navigation system, uses the real-time flight state variables as the input of the guidance network for forward propagation, and obtains the real-time guidance command. The aircraft will generate an angle of attack change according to the command and track the predetermined altitude curve in real time.

2. The method according to claim 1, wherein: The continuous optimal control problem P0 in step (1) is expressed as: where f1 is the state differential equation, expressed as: r, V, Ω, L, and D are dimensionless geocentric distance, velocity, angular velocity of the Earth's rotation, and lift and drag accelerations respectively. J represents the performance index of the continuous optimal control problem. x represents the five state variables of the aircraft. and represent the initial and final dimensionless states of the aircraft respectively. e represents the dimensionless energy of the aircraft. e0 and e f represent the initial and final dimensionless energies of the aircraft respectively. θ and φ are longitude and latitude respectively. R0 is the radius of the Earth. γ is the flight path angle. ψ is the heading angle. σ is the bank angle. σ min and σ max are the lower and upper limits of the bank angle respectively. q and n are the heat flux rate, dynamic pressure, and overload of the aircraft respectively. q max and n max are the maximum values of the heat flux rate, dynamic pressure, and overload that the aircraft can withstand respectively. k Q and ρ represent the heat flux rate calculation coefficient and atmospheric density respectively.

3. The method according to claim 1, wherein: The sequential convex optimal control problem P1 obtained in step (2) is expressed as: s.t. x′ = A(x (k) )x + B(x (k) )u + b(x (k) ) |x(e) - x (k) (e)| ≤ δ Among them, \(u = [u_1, u_2]=[\cos\sigma,\sin\sigma]\) represents the control variable, \(r\) and \(V\) are the dimensionless geocentric distance and velocity respectively, \(J\) represents the performance index of the continuous optimal control problem, and \(x\) represents the five state variables of the aircraft. represents the dimensionless initial state of the aircraft, \(e\) represents the dimensionless energy of the aircraft, \(e_0\) and \(e\) f represent the initial and final dimensionless energies of the aircraft respectively, \(\theta\) and \(\varphi\) are the longitude and latitude respectively, \(\gamma\) is the flight path angle, \(\psi\) is the heading angle, \(\sigma\) is the bank angle, \(\sigma\) min and \(\sigma\) max are the lower and upper limits of the bank angle respectively. represent the slack variables of the longitude and latitude respectively, the superscript \((k)\) indicates that the variable is in the \(k\)-th iteration, \(k\) represents the iteration count, \(A(x\) (k) ) and \(B(x\) (k) ) represent the coefficient matrices of the state variable and the control variable respectively, \(b(x\) (k) ) represents the co-vector matrix, \(c\) θ 、 represent the penalty coefficients of the longitude and latitude respectively, \(c_1\) and \(c_0\) represent the first-order term and the remainder term with respect to the dimensionless geocentric distance \(r\) in the successive linearization process of the performance index, \(\varepsilon\) ψ represents the coefficient of the heading angle regularization term, \(\delta\) represents the trust region constraint vector, represents the convexified process constraint, \(\omega\) l 、\(\omega\) h are the cosine values corresponding to the upper and lower limits of the bank angle respectively. represent the dimensionless geocentric distance, longitude, latitude, flight path angle and heading angle of the final state of the aircraft respectively.

4. The method according to claim 1, characterized in that: In step (2), the P1 is discretized and parameterized using the pseudospectral method, and the resulting sequential second-order cone programming problem P2 is expressed as: Among them, \(u = [u_1, u_2]=[\cos\sigma,\sin\sigma]\) represents the control variable, \(r\) and \(V\) are the dimensionless geocentric distance and velocity respectively, \(J\) represents the performance index of the continuous optimal control problem, \(x\) represents the five state variables of the aircraft, represents the dimensionless initial state of the aircraft, \(e\) represents the dimensionless energy of the aircraft, \(e_0\) and \(e\) f represent the initial and final dimensionless energies of the aircraft respectively, \(\theta\) and \(\varphi\) are the longitude and latitude respectively, \(\gamma\) is the track angle, \(\psi\) is the course angle, \(\sigma\) is the bank angle, \(\sigma\) min and \(\sigma\) max are the lower and upper limits of the bank angle respectively, represent the relaxation variables of longitude and latitude respectively, the superscript \((k)\) indicates that the variable is in the \(k\)-th iteration, \(k\) represents the iteration count, \(A(x\) (k) ) and \(B(x\) (k) ) represent the coefficient matrices of the state variable and the control variable respectively, \(b(x\) (k) ) represents the co-vector matrix, \(c\) θ 、 represent the penalty coefficients of longitude and latitude respectively, \(c_1\) and \(c_0\) represent the first-order term and the remainder term of the dimensionless geocentric distance \(r\) in the successive linearization process of the performance index, \(\varepsilon\) ψ represents the coefficient of the course angle regularization term, \(\delta\) represents the trust region constraint vector, represents the convexified process constraint, \(\omega\) l 、\(\omega\) h represent the cosine values corresponding to the upper and lower limits of the bank angle respectively, represent the dimensionless geocentric distance, longitude, latitude, track angle and course angle of the final state of the aircraft respectively, \(n = 0,1,2...N + 1\) represents each discrete point, \(m = 1,2,3…N\) represents each LG point, \(N\) represents the number of selected LG points, \(w\) represents the integration weight corresponding to each LG point, \(D\) i represents the derivative of the basis function at the LG point.

5. The method according to claim 1, characterized in that: Step (3) uses the interior point method to solve the sequential second-order cone programming problem P2, as follows: (3a) Let the iteration count k = 0, and the initial state quantity at the given discrete points represents the solution of the state variable at the k-th iteration at the n-th discrete point, where n = 0, 1, 2…N+1; (3b)Increment the iteration count k by one, and solve problem P2 using the interior point method to obtain the solution where z (k+1) represents the solution vector for the (k + 1)-th iteration. denotes the solution of the state variable at the (k + 1)-th iteration denotes the solution of the state variable at the (k + 1)-th iteration at the n-th discrete point, where n = 0, 1, 2... N + 1 Denote the control quantity solution of the (k + 1)-th iteration, Denote the control quantity solution of the (k + 1)-th iteration at the m-th discrete point, where m = 1, 2, 3…N, and represents the solution of the slack variable at the (k + 1)-th iteration; (3c) Let a five-dimensional column vector ε represent the convergence judgment condition, and determine whether the state variable solutions at each discrete point satisfy the convergence condition: If not satisfied, let x (k) = x (k+1) , and return to (3b); Otherwise, output the optimal solution The calculation is completed.

6. The method according to claim 1, wherein: Both the hidden layers in the Actor network and the Critic network contain 100 neurons, and the activation function of the hidden layers both uses the relu function, and the output layers both use the tanh function.

7. The method according to claim 1, wherein: In step (7), the target reward and the constraint reward are set as follows: Set the target reward R according to the task objective 1 : where h next represents the height at the next time step after executing the action Δα, Δα represents the angle of attack increment command of the current state, and h0 represents the nominal height corresponding to the next time step; Set the constraint reward R according to the magnitude of the angle of attack increment 2 : R 2 = -10|Δα|。 8. The method according to claim 1, wherein: In step (8), the DDPG algorithm is used to train the Actor network and the Critic network in parallel offline, as follows: (8a) Parameter initialization: Randomly initialize the parameters θ of the Actor_eval network μ and the parameters θ of the Actor_target network μ′ ; Randomly initialize the parameters θ of the Critic_eval network Q and the parameters θ of the Critic_target network Q′ ; Initialize the relevant hyperparameters of the training process: Set the target network replacement update frequency parameter to τ, the target network update cycle number of episodes to T, the size of the experience replay pool to M, the total number of transition processes sampled from the experience replay pool at each time step to batch, and the reward discount rate to λ. (8b) Use the Markov decision process (MDP) model to describe the single-step state transition process of the agent, that is, the agent receives the current environmental state information S t , and based on S t executes the angle of attack increment α t action and then transfers to the new environmental state S t+1 , and at the same time, the agent receives a scalar reward R t+1 ; (8c) Store the single-step transfer process in (8b) into the experience replay pool in the form of (S t , a t , R t+1 , S t+1 ); Among them, S t represents the state at the current moment, a t represents the angle of attack increment action taken at the current moment, R t+1 represents the reward value obtained after executing the action, S t+1 represents the state at the next moment; (8d) During each training process, randomly draw some samples from the experience replay pool and input them into the Actor network and the Critic network. Update the parameters of the Actor_eval sub-network in the Actor network by maximizing the cumulative reward through the policy gradient of the performance index, and update the parameters of the Critic_eval sub-network in the Critic network by minimizing the mean square error between the current cumulative reward and the target cumulative reward. Update the parameters of the Actor_target sub-network in the Actor network and the parameters of the Critic_target sub-network in the Critic network according to τ every T episodes. (8e) Repeat step (8d). When the cumulative reward of the agent converges to a maximum value, the trained guidance network is obtained.

Images