An online guidance method for rocket substage attitude flip landing

Abstract

The invention provides an online guidance method for attitude flip landing of a rocket sub-stage, comprising: S1, establishing a six-degree-of-freedom landing dynamics model of a rocket's final stage; S2, constructing a six-degree-of-freedom powered soft landing trajectory optimization model of a rocket's final stage that meets constraint conditions; S3, converting non-convex constraints in the trajectory optimization model of S2 into convex constraints by linearization, so as to obtain a convex trajectory optimization model; S4, discretizing the convex trajectory optimization model in S3; S5, adding dynamic relaxation variables to the linearized dynamics equations in S3; S6, designing a trust region constraint to limit the range of variation of a reference trajectory; S7, determining an initial iterative reference trajectory; S8, solving a discrete convex model; S9, repeating S8 to iteratively solve so that the trajectory converges to an optimal trajectory, and completing the trajectory optimization of a sampling point in a guidance cycle; S10, updating the optimal instruction using the trajectory optimization result, and directly using it as a guidance signal, and finally completing the online guidance for attitude flip landing of the rocket sub-stage.

Description

Technical Field

[0001] The invention belongs to the technical field of aircraft guidance, and in particular relates to an online guidance method for rocket sub-stage attitude flip landing. Background Art

[0002] With the development of aerospace technology, reuse has become one of the main development directions of future aerospace transportation systems. At present, the reuse technology of the first stage of a launch vehicle is developing rapidly, but for rocket stages with higher potential energy (such as the last stage of a rocket), the recovery is more difficult, and its flight profile is quite different from that of the first stage of a launch vehicle. The vertical power landing stage of the first stage of a traditional launch vehicle does not require thrust to participate in the attitude control of the aircraft (such as Falcon9), and the attitude adjustment is mainly achieved through RCS (Reaction Control System) and grid rudders. Therefore, most scholars usually use a three-degree-of-freedom dynamic model that does not consider the attitude for the optimization model of the landing trajectory of the first stage of a rocket. However, since the last stage of a reusable rocket descends in a state of aerodynamic deceleration at a large angle of attack at the initial moment of powered landing, it is necessary to use engine thrust to perform attitude flip maneuvers and complete vertical powered landing in the final landing stage, so a new landing guidance method needs to be designed. Summary of the invention

[0003] In order to overcome the deficiencies in the prior art, the inventors have conducted intensive research and provided an online guidance method for rocket substage attitude flip landing based on six-degree-of-freedom sequence convex optimization. On the basis of the traditional three-degree-of-freedom return landing guidance method, attitude variables were introduced, and a six-degree-of-freedom return landing problem was established to solve the attitude flip maneuver problem during landing. At the same time, the position, speed, and attitude meet the terminal constraint requirements, thereby completing the present invention.

[0004] The technical solution provided by the present invention is as follows:

[0005] At present, most of the power landing trajectory optimization methods are based on the 3-DOF center-of-mass dynamics model similar to the rocket stage. The attitude change is not considered during the landing process. This is because the attitude change of the rocket stage in the vertical power landing phase is small, and thrust is usually not required to participate in the attitude control of the aircraft. The attitude adjustment is mainly achieved through RCS and grid rudder. However, since the initial attitude of the final stage landing of the "Starship" rocket is in a state of aerodynamic deceleration at a large angle of attack, it needs to use engine thrust to perform attitude flip maneuvers and complete vertical powered landing in the final landing phase. Therefore, the vertical landing trajectory optimization of this type of aircraft needs to consider the impact of thrust on attitude.

[0006] Aiming at the characteristics of the "Starship" rocket sub-stage powered landing process using thrust to perform attitude flip maneuvers, a six-degree-of-freedom trajectory optimization method based on improved sequential convex optimization is proposed. This method introduces attitude quaternion and angular velocity variables to establish the six-degree-of-freedom landing trajectory optimization problem, and transforms the non-convex trajectory optimization problem into a series of finite-dimensional second-order cone problems through convex processing and discretization. Virtual control and dynamic trust region update strategies are designed to ensure the feasibility and convergence of the second-order cone problem, and then gradually converge to the optimal trajectory by solving a series of second-order cone problems. Finally, based on the six-degree-of-freedom powered landing trajectory optimization method, a rocket powered soft landing rolling time-domain MPC (Model Predictive Control) guidance framework is established, and the effectiveness of the established MPC guidance method based on trajectory sequence convex optimization is verified through simulation experiments.

[0007] The specific plan is as follows

[0008] A rocket substage attitude flip landing online guidance method based on six-degree-of-freedom sequence convex optimization comprises the following steps:

[0009] S1, introduce attitude quaternion and rotational angular velocity to describe the attitude motion of the rocket's final stage, and establish a six-degree-of-freedom landing dynamics model for the rocket's final stage;

[0010] S2, taking the shortest landing time as the objective function, constructs a six-degree-of-freedom powered soft landing trajectory optimization model for the rocket's final stage that meets the constraints;

[0011] S3, converting the non-convex constraints in the trajectory optimization model in step S2 into convex constraints by linearization to obtain a convex trajectory optimization model, wherein the non-convex constraints of the trajectory optimization model include dynamic equations and thrust magnitude constraints;

[0012] S4, using a first-order retainer method to discretize the convex trajectory optimization model in step S3 to obtain a discrete convex model;

[0013] S5, adding dynamic slack variables to the linearized dynamic equation in step S3, and adding the slack variables to the objective function in the form of a penalty function;

[0014] S6, design trust region constraints to limit the range of variation of the reference trajectory;

[0015] S7, obtain the initial iterative reference trajectory through simple linear interpolation of initial and terminal state constraints;

[0016] S8, the initial iteration reference trajectory in step S7 is used as the initial iteration value of the sequence convexification, and the constraints of steps S5 and S6 are considered to solve the discrete convexification model obtained in step S4;

[0017] S9, repeating step S8 to iteratively solve so that the trajectory converges to the optimal trajectory, completing the trajectory optimization of a guidance cycle sampling point;

[0018] S10, using the trajectory optimization result of a guidance cycle sampling point in step S9 as the guidance instruction generator, updating and executing steps S2-S9 according to the current state, updating the optimal instruction with the trajectory optimization result, and directly using it as the guidance signal, and finally completing the online guidance of the rocket sub-stage attitude flip landing.

[0019] The rocket substage attitude flip landing online guidance method based on six-degree-of-freedom sequence convex optimization provided by the present invention has the following beneficial effects:

[0020] (1) In view of the characteristics of the rocket's final stage attitude flip landing, the present invention introduces attitude variables under the traditional three-degree-of-freedom powered landing model to establish a six-degree-of-freedom powered landing trajectory optimization problem, thereby improving mission adaptability.

[0021] (2) The present invention improves the traditional fixed trust region sequence convexification method and designs a dynamic trust region strategy based on linearization error measurement, thereby improving the solution efficiency of the trajectory optimization method.

[0022] (3) The present invention establishes a model prediction guidance framework based on the six-degree-of-freedom powered landing trajectory optimization method, successfully applies the online trajectory optimization method to the guidance and control system, and realizes the online guidance of the rocket substage attitude flip landing throughout the whole process. BRIEF DESCRIPTION OF THE DRAWINGS

[0023] Figure 1 This is a schematic diagram of the power landing of the present invention;

[0024] Figure 2 It is a framework structure diagram of the guidance method of the present invention;

[0025] Figure 3 It is a simulation result diagram of the guidance method of the present invention, wherein Figure 3 a is a 6DOF three-dimensional trajectory diagram; Figure 3 b is the 6DOF longitudinal plane trajectory projection; Figure 3 c is the position change curve diagram; Figure 3 d is the velocity component variation curve; Figure 3 e is the curve diagram of thrust size variation; Figure 3 f is the mass change curve; Figure 3 g is the velocity change curve; Figure 3 h is the curve diagram of the change of angular velocity component. DETAILED DESCRIPTION

[0026] The following detailed description of the present invention will make the features and advantages of the present invention more clear and explicit.

[0027] The word "exemplary" is used exclusively herein to mean "serving as an example, embodiment, or illustration." Any embodiment described herein as "exemplary" is not necessarily to be construed as preferred or advantageous over other embodiments. Although various aspects of the embodiments are shown in the drawings, the drawings are not necessarily drawn to scale unless otherwise noted.

[0028] The present invention provides a rocket substage attitude flip landing online guidance method based on six-degree-of-freedom sequence convex optimization, comprising the following steps:

[0029] S1, introduce attitude quaternion and angular velocity to describe the attitude motion of the rocket's last stage, and establish a six-degree-of-freedom landing dynamics model for the rocket's last stage.

[0030] The target point ground coordinate system F L and the body coordinate system F B is the reference coordinate system, describing the relative position and attitude relationship between the aircraft and the target point, such as Figure 1 Considering that the rocket's last stage needs to perform an attitude flip maneuver through engine thrust in the final powered landing phase to achieve vertical powered landing, the attitude quaternion (q B / L ) and the angular velocity (ω B ) is used to describe the attitude motion of the rocket body, and the specific form of the six-degree-of-freedom landing dynamics model is as follows:

[0031]

[0032] Where m is the mass of the aircraft; α m =1 / (I sp g0), I sp is the specific impulse, g0 is the standard earth gravity acceleration; t is the running time, r L 、v L F is the ground coordinate system of the aircraft at the target point L Position and velocity vector under D L =-0.5ρSC D ||v L ||v L is the ground coordinate system F at the target point L aerodynamic drag under, ρ atmospheric density, S is the reference area, C D is the resistance coefficient; T B is the coordinate system F B The thrust under B is the coordinate system F B Lower body coordinate system F B Relative to the target point ground coordinate system F L Angular velocity vector; qB / L F is the ground coordinate system of the target point L Relative body coordinate system F B The unit attitude quaternion of g L is the gravity vector; r T,B J is the arm from the thrust point to the center of mass moment; B is the inertia matrix; C T / B (t) is the body coordinate system F B To the target point ground coordinate system F L The rotation matrix of is expanded as follows:

[0033]

[0034] Among them, q0, q1, q2, q3 are the corresponding attitude quaternions;

[0035] Ω(ω B ) is the multiplication expansion form of the quaternion:

[0036]

[0037] Among them, ω B =[ω Bx ω By ω Bz ]

[0038] [ζ × ] is the cross product expansion form:

[0039]

[0040] Where ζ=[ζ x ζ x ζ x ] is a variable.

[0041] S2, taking the shortest landing time as the objective function, constructs a six-degree-of-freedom powered soft landing trajectory optimization model for the rocket's final stage that meets the constraints.

[0042] For the powered landing trajectory optimization problem, the minimum landing time (MLT) is used as the objective function of the trajectory optimization problem, namely:

[0043]

[0044] t f for landing time;

[0045] At the same time, the following constraints must be met during the powered soft landing process:

[0046] a) Fuel constraints

[0047] mdry ≤m(t)

[0048] Among them, m dry is the dry weight of the aircraft.

[0049] b) Thrust restraint

[0050] T min ≤||T B (t)||2≤T max

[0051] cosδ max ||T B (t)||2≤e1·T B (t)

[0052] T min and T max are the lower and upper limits of the thrust, δ max is the maximum thrust deflection angle, and e1 is the unit matrix.

[0053] c) Slope angle constraint

[0054] Similarly, in order to prevent the rocket's final stage from colliding with ground slopes, bumps and other terrain during landing, the following slope constraints are set:

[0055]

[0056] H 23 =[e2;e3]

[0057] Among them, e i represents a unit vector whose i-th element is 1, γ gs is the maximum slope angle.

[0058] d) Tilt angle constraint

[0059] In order to avoid the optimized trajectory having too large an inclination angle, it is necessary to limit the inclination angle.

[0060]

[0061] Among them, θ max is the maximum tilt angle. In order to meet the requirements of attitude flip maneuver and prevent the attitude from swinging back and forth, the maximum tilt angle θ max The constraint is set to shrink dynamically over time, and cosθ max Set it to a cubic polynomial related to time t, that is:

[0062]

[0063] Among them, k0, k1, k2, and k3 are polynomial coefficients, which are selected according to the control characteristics of the aircraft.

[0064] e) Angular velocity constraint

[0065] ||ω B (t)||2≤ω max

[0066] Among them, ω max is the maximum angular velocity;

[0067] f) Boundary Constraints

[0068] The initial and final state constraints of rocket powered soft landing are:

[0069]

[0070] Among them, m wet is the initial total mass of the rocket including fuel, r L,i , v L,i are the initial position and velocity vector, q B / L,i ,ω B,i are the initial attitude quaternion and angular velocity vector, v L,f is the terminal velocity, q B / L,f is the terminal attitude quaternion.

[0071] S3, converting the non-convex constraints in the trajectory optimization model in step S2 into convex constraints by linearization to obtain a convex trajectory optimization model, wherein the non-convex constraints in the trajectory optimization model include dynamic equations and thrust size constraints.

[0072] Because the above-mentioned powered soft landing trajectory optimization problem contains many non-convex constraints, mainly existing in the dynamic equation and thrust size constraints, this step mainly converts the non-convex constraints of these two parts into convex constraints. The present invention adopts the method of sequential convex optimization to convert the non-convex trajectory optimization problem into a convex second-order cone problem (SOCP), and finally converges the optimization result to the optimal trajectory by solving the second-order cone problem through multiple iterations.

[0073] a) Linearization of the kinetic model

[0074] For the non-convex form of the dynamic equation, the linearization method is used to convexify the dynamic equation constraints. For the convenience of representation, the state variable x(t)∈R is defined 14 and control variable u(t)∈R 3 .

[0075]

[0076] u(t)=T B (t)

[0077] The kinetic equation can be simplified as:

[0078]

[0079] Normalizing the time variable t to τ∈[0,1], the above formula can be expressed as:

[0080]

[0081] Let τ be the scaling factor of t. Taking τ as the independent variable, the above kinetic equation can be expressed as:

[0082]

[0083] In order to transform the entire nonlinear dynamics equation into a convex constraint form, the dynamics equation can be transformed into a convex constraint form on the reference trajectory (given by and The first-order Taylor expansion of the equation is performed to obtain the approximate linearized dynamic equation as follows:

[0084]

[0085] b) Control Constraint Linearization

[0086] From the thrust constraint form, we can see that the minimum thrust constraint is a non-convex constraint. We can make the minimum thrust constraint convex by linearization and define the function g:R 3 →R.

[0087] g(u(τ))=T min -||u(τ)||2≤0

[0088] Similarly, the above formula is applied to the reference trajectory (given by and A first-order Taylor expansion is performed at the position of the equation ( ) to obtain an approximate linearized minimum thrust constraint.

[0089]

[0090] in, is the reference trajectory control variable.

[0091] S4, using the first-order retainer method to discretize the convex trajectory optimization model in step S3 to obtain a discrete convex model.

[0092] In order to transform the continuous-time optimization control problem into an optimization problem with optimization parameters in finite dimensions, the entire problem needs to be discretized. This step uses the discretization method of the first-order retainer to discretize the standardized time into K equally spaced discrete points.

[0093]

[0094] The control variable is interpolated for each step by using a first-order holder. k ,τ k+1 ],make

[0095] u(τ)=α k (τ)u k +β k (τ)u k+1 ,τ∈[τ k ,τ k+1 ],k=0,1,…,K-2

[0096]

[0097]

[0098] Let Φ(τ k+1 ,τ k ) indicates the state x(τ when there is no input k ) to x(τ k+1 ) state transfer matrix. The state transfer matrix Φ(τ,τ k ) can be expressed in differential form as:

[0099]

[0100] Through the above discretization method, the dynamic equation can be discretized as follows:

[0101]

[0102] ξ is the integration variable.

[0103] S5, adding dynamic relaxation variables to the linearized dynamic equation in step S3, and adding the relaxation variables to the objective function in the form of a penalty function.

[0104] Since the linearization process of non-convex dynamic equations easily causes sub-problems to be artificially infeasible, it is necessary to introduce dynamic relaxation variables into the linear dynamic equations. Used to compensate for the linearization approximation error, that is

[0105]

[0106] Where superscript i represents the iteration number, and k represents the position of the discrete point. In order to make the dynamic linearization error gradually converge, the dynamic relaxation is added to the objective function in the form of a penalty function, and the objective function is modified to the following form:

[0107]

[0108] Among them, w v is the weight coefficient.

[0109] S6, design trust region constraints to limit the range of variation of the reference trajectory.

[0110] Due to the local validity of the linear expansion, that is, the optimization variables can only guarantee a good linear approximation effect when they take values ​​near the reference point, so the range of variation of the reference trajectory is limited by designing trust region constraints.

[0111] The reference trajectory change is defined as:

[0112]

[0113]

[0114]

[0115] Since the landing dynamics equations are dimensionless, the following trust region constraints can be added to the subproblem:

[0116]

[0117] in, is the trust region constraint, |δx i |1,|δu i |1,|δσ i |1 is the change of state variables, control variables and time variables. Although the fixed trust region constraint can obtain better optimization results, the fixed trust region constraint is a "rigid hard constraint". When the quality of the initial reference trajectory is low, the optimization calculation of the sub-problem will be infeasible. Therefore, the fixed trust region constraint is set too large, which is easy to increase the number of iterative calculations of the sequential convex optimization. On the other hand, the perturbation allowable range changes continuously during the iteration of the sequential convex optimization algorithm, which is closely related to the approximation degree of the original nonlinear function and the corresponding linearized function. Therefore, the trust region change can be adjusted by measuring the linearization error, improving the convergence speed and expanding the optimization range.

[0118] Design a simple and direct linear error measurement method. Define the following two functions and Indicates that the optimization variable value is taken as When , the "degree of satisfaction" of the original nonlinear equality and inequality constraints; and the function Indicates that the reference trajectory is taken as The disturbance is taken as The linearization approximates the “degree of satisfaction” of the equality and inequality constraints when .

[0119] Based on this, the linearization error metric coefficient ρ for the kth iteration result is defined as k as follows:

[0120]

[0121] in, It means: for the original problem, the result of the kth iteration is Relative to the k-1th iteration result Improvements in constraint satisfaction; It means: for the linearized approximate subproblem, the improvement in constraint satisfaction of the k-th iteration result relative to the k-1-th iteration result. k When it is close to 1, it is considered that the current linear approximation accuracy is high; and when ρ k If it is much smaller than 1, the current linearization accuracy is considered to be too low. Therefore, the following dynamic trust region strategy can be designed based on the linearization error measurement coefficient: Dynamic adjustment is performed according to the following piecewise function:

[0122]

[0123] Among them, σ1, σ2 are trust region penalty parameter update proportional coefficients, satisfying 0<σ2<1, σ1>1. r0, r1, r2 are trust region update criterion parameters, satisfying r0<r1<r2<1.

[0124] S7, obtain the initial iterative reference trajectory through simple linear interpolation of the initial and terminal state constraints, use the initial iterative reference trajectory as the initial value of the sequence convexification iteration, consider the constraints of steps S5 and S6, and solve the discrete convexification model obtained in step S4; repeat the iterative solution to make the trajectory converge to the optimal trajectory, and complete the trajectory optimization of a guidance cycle sampling point.

[0125] The above discretization and convexification methods are used to transform the continuous optimal control problem into a sequential convex problem. The inner loop is a second-order cone programming problem, and the outer loop is an iterative loop process. The interior point method solver is used to solve a series of SOCP sub-problems. The design variables include the state variables, control variables and landing time of the optimal control problem. The constraints are dynamic differentials, path constraints and boundary conditions.

[0126] The initialization trajectory of the sequential convexification algorithm can be obtained by simple linear interpolation of the initial and terminal state constraints, and the sequential convexification algorithm is used to gradually iterate and converge to the optimal trajectory.

[0127] The sequential convexification algorithm solution process of the rocket's final stage powered soft landing problem is as follows:

[0128]

[0129]

[0130] S8, uses the trajectory optimization result of a guidance cycle sampling point as the guidance instruction generator, updates the execution steps S2 to S7 according to the current state, uses the trajectory optimization result to update the optimal instruction, and directly uses it as a guidance signal to finally complete the online guidance of the rocket sub-stage attitude flip landing.

[0131] With the above six-degree-of-freedom trajectory optimization method as the core, the optimization method is embedded in the model predictive control (MPC) framework. MPC guidance uses the trajectory optimization algorithm as a guidance instruction generator, continuously updates the optimal instructions according to the current state by performing optimization calculations in the online rolling time domain, and directly uses the latest instructions as guidance signals. MPC guidance updates control instructions by performing optimization calculations through real-time sampling of state information, so it has strong adaptability to the uncertainty of environmental parameters and disturbance deviations. For the MPC guidance method, its stability, accuracy and robustness depend on the rapidity of solving the corresponding optimization problem, that is, the optimization calculation update frequency is an important evaluation indicator of the MPC algorithm capability.

[0132] Consider the following nonlinear autonomous system model:

[0133]

[0134] Among them, x is the state variable and u is the control variable. Assume that at the sampling time t i The actual state of the aircraft is x R (t i ), based on x R (t i ) According to the previous trajectory optimization method, the optimal control u is obtained * (t,x R (t i )), t∈[t i ,t fi ]. According to x R (t i ) to obtain the optimal control u * (t,x R (t i )) is used as the guidance control command in this guidance cycle [t i ,t i+1 ], the next guidance cycle is based on the sampling time t i+1 The actual system state is x R (t i+1 ) as the initial state of the optimization algorithm, that is, the actual state x of the system at the initial moment of each guidance cycle R As the initial state of each optimization update calculation. It can be seen that the optimal guidance algorithm of the MPC framework can be regarded as updating the initial state in each guidance cycle to solve the trajectory optimization problem and thus update the guidance instructions.

[0135] To facilitate the discussion of the following methods, the trajectory optimization problem of each guidance cycle in the MPC framework can be abstracted into a general numerical optimization problem (MPC0):

[0136]

[0137]

[0138]

[0139]

[0140] in, is the state variable of the trajectory optimization problem, t i is the update sampling time of the current problem, is the terminal flight time of the current problem, represents the process constraint, u(t)∈U is the control constraint, represents the initial state constraint; Represents a terminal state constraint.

[0141] In the optimization update of MPC guidance instructions, there may be situations where the optimization calculation is not feasible, especially when approaching the landing point, because full thrust control is usually used in the late stage of landing, which reduces the adjustment ability at the end of landing. R (t i ) relative to the (i-1)th calculation to obtain x i-1 (t i ) results have a large deviation, which may lead to the inability to meet the various constraints of the trajectory optimization problem, making the current optimization calculation infeasible.

[0142] In order to avoid the above situation and prevent the trajectory optimization algorithm from being infeasible at the end of landing, the terminal state constraints of the trajectory optimization problem (MPC0) can be relaxed and penalized in the performance index to establish a new trajectory optimization problem RMPC:

[0143]

[0144]

[0145]

[0146]

[0147] in, is the relaxed terminal state constraint, J0 is the original performance index, and ω is the penalty parameter for the deviation between the terminal state variable and the desired terminal constraint. At the end of landing, the feasibility of the optimization problem is ensured by converting the strict terminal state into a soft constraint. At the same time, the terminal state error is introduced into the objective function to minimize the terminal state error.

[0148] On the other hand, when the MPC guidance instructions are updated, the flight time is less than the last calculation, so the number of discrete points currently used can be reduced in the next guidance cycle without affecting the discrete accuracy. Therefore, in order to improve the computational efficiency of the trajectory optimization problem, a dynamic adjustment strategy for the number of discrete points is designed:

[0149]

[0150] Among them, floor(·) is the floor rounding function, N i is the number of discrete points for the i-th optimization calculation, is the flight time of the i-th optimization result, Δt is the guidance cycle time, is an estimate of the flight time of the (i+1)th optimization calculation. κ is the discrete point number adjustment factor, and its value is less than 1, so that the rate of decrease in the number of discrete points is less than the rate of decrease in the flight time. Therefore, as the flight time continues to decrease, the number of discrete points decreases, making the problem scale gradually smaller and the calculation speed gradually increase.

[0151] In summary, the entire MPC guidance algorithm is as follows Figure 2 As shown, the specific process is as follows:

[0152]

[0153]

[0154] The simulation results of the entire guidance method are as follows: Figure 3 shown.

[0155] Figure 3 a and Figure 3 b shows the trajectory and attitude changes of the aircraft during the landing process. The dotted line represents the landing trajectory, the thick solid line represents the axis of the rocket body, the thin solid line represents the thrust vector, and the length of the thin solid line represents the thrust size. Figure 3 e shows the thrust change curve, and the thrust change form is similar to the "Bang-Bang" control form. Figure 3 c. Figure 3 d. Figure 3 f- Figure 3h respectively gives the changes of the position, velocity component, mass, velocity magnitude and angular velocity state variables of the rocket's final stage landing process. The entire MPC calculation guidance simulation flight time is 14.707s. The simulation results show that: under the set initial position and attitude, it can land at the designated location through thrust control, and the position, velocity and attitude meet the requirements. The position deviation in the horizontal direction of the terminal is 0.028m; the vertical velocity Vz of the landing point is -0.563m / s, the horizontal velocity Vx is 0.171m / s, and Vy is 0.228m / s; the terminal angular velocity is 0.775° / s, and the inclination angle is 0.423°. The final mass of the rocket's final stage when landing is 149293.629kg, and the fuel consumption is 10706.371kg, which meets the fuel consumption constraint. It can be seen from the guidance simulation results that the guidance method of the present invention can achieve high guidance accuracy by implementing online update guidance instructions, and the fuel consumption is low, and the overall guidance effect is good.

[0156] The present invention has been described in detail above in conjunction with specific implementations and exemplary examples, but these descriptions cannot be understood as limiting the present invention. Those skilled in the art understand that, without departing from the spirit and scope of the present invention, a variety of equivalent substitutions, modifications or improvements may be made to the technical solution of the present invention and its implementation methods, all of which fall within the scope of the present invention. The scope of protection of the present invention shall be subject to the attached claims.

[0157] The contents not described in detail in the specification of the present invention belong to the common knowledge of those skilled in the art.

Claims

1. A rocket substage attitude flip landing online guidance method based on six-degree-of-freedom sequence convex optimization, characterized in that: The steps include: S1, introduce attitude quaternion and rotational angular velocity to describe the attitude motion of the rocket's final stage, and establish a six-degree-of-freedom landing dynamics model for the rocket's final stage; S2, taking the shortest landing time as the objective function, constructs a six-degree-of-freedom powered soft landing trajectory optimization model for the rocket's final stage that meets the constraints; S3, converting the non-convex constraints in the trajectory optimization model in step S2 into convex constraints by linearization to obtain a convex trajectory optimization model, wherein the non-convex constraints in the trajectory optimization model include dynamic equations and thrust magnitude constraints; S4, using a first-order retainer method to discretize the convex trajectory optimization model in step S3 to obtain a discrete convex model; S5, adding dynamic slack variables to the linearized dynamic equation in step S3, and adding the slack variables to the objective function in the form of a penalty function; S6, design trust region constraints to limit the range of variation of the reference trajectory; S7, obtain the initial iterative reference trajectory through simple linear interpolation of initial and terminal state constraints; S8, the initial iteration reference trajectory in step S7 is used as the initial iteration value of the sequence convexification, and the constraints of steps S5 and S6 are considered to solve the discrete convexification model obtained in step S4; S9, repeating step S8 to iteratively solve so that the trajectory converges to the optimal trajectory, completing the trajectory optimization of a guidance cycle sampling point; S10, using the trajectory optimization result of a guidance cycle sampling point in step S9 as the guidance instruction generator, updating and executing steps S2-S9 according to the current state, using the trajectory optimization result to update the optimal instruction, and directly using it as a guidance signal, finally completing the online guidance of the rocket sub-stage attitude flip landing.

2. The rocket substage attitude flip landing online guidance method based on six-degree-of-freedom sequence convex optimization according to claim 1 is characterized in that: In step S1, the six-degree-of-freedom landing dynamics model is expressed as follows: Where m is the mass of the aircraft; I sp is the specific impulse, g0 is the standard earth gravity acceleration; t is the running time, r L 、v L F is the ground coordinate system of the aircraft at the target point L Position and velocity vector under D L =-0.5ρSC D ||v L ||v L is the ground coordinate system F at the target point L aerodynamic drag under, ρ atmospheric density, S is the reference area, C D is the resistance coefficient; T B is the coordinate system F B The thrust under B is the coordinate system F B Lower body coordinate system F B Relative to the target point ground coordinate system F L Angular velocity vector; q B / L F is the ground coordinate system of the target point L Relative body coordinate system F B The unit attitude quaternion of g L is the gravity vector; r T,B J is the arm from the thrust point to the center of mass moment; B is the inertia matrix; C T / B (t) is the body coordinate system F B To the target point ground coordinate system F L The rotation matrix of is expanded as follows: Among them, q0, q1, q2, q3 are the corresponding attitude quaternions; Ω(ω B ) is the multiplication expansion form of the quaternion: Among them, oh B =[ω Bx oh By oh Bz ] [ζ × ] is the cross product expansion form: Where, ζ=[ζ x g x g x ] for variables.

3. The method for online guidance of rocket substage attitude flip landing based on six-degree-of-freedom sequence convex optimization according to claim 1, characterized in that: In step S2, the constraints include: fuel constraint, thrust constraint, slope angle constraint, inclination angle constraint, angular velocity constraint, and boundary constraint; the shortest landing time is used as the objective function of the trajectory optimization model, that is: t f For landing time.

4. The method for online guidance of rocket substage attitude flip landing based on six-degree-of-freedom sequence convex optimization according to claim 2 is characterized in that: In step S3, the process of converting the dynamic equation constraint into a convex constraint is as follows: define the state variable x(t)∈R 14 and control variable u(t)∈R 3 ; u(t)=T B (t) The kinetic equation is abbreviated as: The dynamic equation with the time variable t normalized to τ∈[0,1] can be expressed as: The dynamic equation is expanded by first-order Taylor on the reference trajectory to obtain the approximate linearized dynamic equation; The process of converting the thrust magnitude constraint into a convex constraint is as follows: The non-convex minimum thrust constraint is: g(u(τ)) = T min -||u(τ)||2≤0, where T min is the lower limit of the thrust size; the non-convex minimum thrust constraint is expanded by a first-order Taylor on the reference trajectory to obtain the convex minimum thrust constraint.

5. The method for online guidance of rocket substage attitude flip landing based on six-degree-of-freedom sequence convex optimization according to claim 4 is characterized in that: In step S4, the convex trajectory optimization model in step S3 is discretized by using a first-order retainer method, which is specifically implemented in the following manner: using a first-order retainer discretization method to discretize the standardized time into K equally spaced discrete points; The control variable is interpolated for each step by using a first-order holder. k ,τ k+1 ],make u(τ)=a k (t)u k +b k (t)u k+1 ,τ∈[τ k ,t k+1 ],k=0,1,…,K-2 Let Φ(τ k+1 ,τ k ) indicates the state x(τ when there is no input k ) to x(τ k+1 ) state transfer matrix, state transfer matrix Φ(τ,τ k ) is expressed in differential form as: Through the above discretization method, the dynamic equation is discretized into: Where ξ is the integral variable.

6. The method for online guidance of rocket substage attitude flip landing based on six-degree-of-freedom sequence convex optimization according to claim 5, characterized in that: In step S5, the adding of dynamic slack variables to the linearized dynamic equation in step S3, and adding the slack variables to the objective function in the form of a penalty function, is specifically implemented in the following manner: Introducing dynamic relaxation variables into linear dynamic equations Used to compensate for the linearization approximation error, namely: Among them, the superscript i represents the iteration number, and k represents the position of the discrete point. In order to make the dynamic linearization error gradually converge, the dynamic relaxation is added to the objective function in the form of a penalty function. The modified objective function is: Among them, w v is the weight coefficient.

7. The method for online guidance of rocket substage attitude flip landing based on six-degree-of-freedom sequence convex optimization according to claim 6, characterized in that: In step S6, the trust region change is adjusted by measuring the linearization error. The trust region constraints are as follows: in, is the trust region constraint; |δx i |1,|δu i |1,|δσ i |1 is the change of state variable, control variable and time variable; According to the linearization error metric coefficient ρ k Design the following dynamic trust region strategy to constrain the trust region Dynamic adjustment is performed according to the following piecewise function: Among them, σ1, σ2 are trust region penalty parameter update proportional coefficients, satisfying 0<σ2<1, σ1>1; r0, r1, r2 are trust region update criterion parameters, satisfying r0<r1<r2<1.

8. The method for online guidance of rocket substage attitude flip landing based on six-degree-of-freedom sequence convex optimization according to claim 7, characterized in that: The linearization error metric coefficient ρ k for: in, Indicates that the optimization variable value is taken as When the original nonlinear equality and inequality constraints are satisfied, the function Indicates that the reference trajectory is taken as The disturbance is taken as The linearization approximates the "degree of satisfaction" of the equality and inequality constraints when .

9. The method for online guidance of rocket substage attitude flip landing based on six-degree-of-freedom sequence convex optimization according to claim 8, characterized in that: Step S10 also includes relaxing the terminal state constraints of the trajectory optimization model at the landing end and penalizing it in the performance index to establish a new trajectory optimization model: in, is the state variable of the trajectory optimization model, t i is the update sampling time of the current model, is the terminal flight time of the current model, represents the process constraint, u(t)∈U is the control constraint, is the relaxed terminal state constraint, J0 is the original performance index, and ω is the penalty parameter for the deviation of the terminal state variable from the expected terminal constraint.

10. The method for online guidance of rocket substage attitude flip landing based on six-degree-of-freedom sequence convex optimization according to claim 9, characterized in that: In step S10, when the guidance instruction is updated, a dynamic adjustment strategy for the number of discrete points is designed: Among them, floor(·) is the floor rounding function, N i is the number of discrete points for the i-th optimization calculation, is the flight time of the i-th optimization result, Δt is the guidance cycle time, is an estimate of the flight time of the (i+1)th optimization calculation; κ is an adjustment factor for the number of discrete points, and its value is less than 1.

Images