Time-coordinated robust guidance method for multiple hypersonic vehicles reentry

Abstract

This invention discloses a time-coordinated robust guidance method for hypersonic vehicle reentry. Based on the three-degree-of-freedom equations of motion of energy, a model of the coordinated reentry guidance problem of multiple hypersonic vehicles under time constraints is established. A piecewise function is used to transform the coordinated time constraint into a normalized error constraint, which is then introduced into the roll angle amplitude update process. Simultaneously, the time constraint is transformed into a heading angle error corridor width constraint. An unscented Kalman filter algorithm is used to correct the roll angle command based on the identification results of aerodynamic parameter uncertainties during flight. This invention features high robustness and strong controllability.

Description

Technical Field

[0001] This invention relates to the field of aircraft control, and more specifically to a cooperative reentry guidance method. Background Technology

[0002] Hypersonic reentry vehicles are aircraft that achieve long-distance non-ballistic maneuvering reentry through unpowered gliding motion in near space. They can reach speeds exceeding Mach 20 and possess rapid response, resistance to interception, and precise terminal strike capabilities, making them extremely promising for applications. Hypersonic reentry vehicle cooperative operations refer to multiple aircraft working together in a network to replace traditional individual vehicles, leveraging collective advantages to more effectively accomplish reconnaissance, penetration, and strike missions.

[0003] Hypersonic reentry vehicles involve multiple technologies, including overall design, structural materials, aerodynamics, and flight control, with guidance technology being the most critical. Currently, cooperative guidance technology for hypersonic reentry vehicles is still in the conceptual research and technological exploration stage, with related cooperative guidance law research focusing on its practical application in cooperative combat missions. Hypersonic vehicles fly at high speeds in near-space for extended periods, requiring them to meet set mission requirements while also being constrained by numerous complex conditions that must be considered during guidance design. Based on the characteristics of these constraints, they can be categorized into two main types: process constraints and endpoint constraints. The numerous constraints and complex mission requirements significantly increase the difficulty of designing guidance systems, and the challenge of multi-vehicle cooperative guidance missions is even greater than that of single-vehicle missions. During the glide phase maneuvering guidance of a hypersonic reentry vehicle, the flight altitude and speed change drastically, and the aerodynamic and aerothermal characteristics are complex. Its guidance and control system is inevitably affected by various disturbances, such as atmospheric density, aerodynamic parameters, vehicle mass, and initial reentry state, all of which are subject to varying degrees of disturbance. Because traditional flight corridor and reference flight profile designs typically do not consider the impact of parameter disturbances, or the disturbances considered differ significantly from the actual situation, and aircraft are inevitably affected by various disturbance factors during actual flight.

[0004] Furthermore, these disturbances and perturbations cause fluctuations in the guidance system parameters, resulting in guidance information that is often not a completely stationary signal. Under these circumstances, the stability and robustness of the guidance system decrease, directly affecting the guidance effect and potentially leading to instability during the gliding maneuver guidance process of hypersonic vehicles.

[0005] Therefore, a new technical solution is needed to solve the above problems. Summary of the Invention

[0006] To address the problems arising from existing technologies, this invention provides a robust and controllable time-coordinated reentry robust guidance method for multiple hypersonic vehicles.

[0007] To achieve the above objectives, the present invention provides a time-coordinated robust reentry guidance method for multiple hypersonic vehicles, employing the following technical solution:

[0008] A time-coordinated reentry robust guidance method for multiple hypersonic vehicles includes the following steps:

[0009] S1. Establish the dimensionless three-degree-of-freedom motion equations based on energy for multiple hypersonic reentry vehicles, considering endpoint constraints, process constraints, and time constraints, and establish a dynamic model for the reentry cooperative guidance problem of multiple hypersonic vehicles.

[0010] S2. In each guidance cycle, the design error coefficient weight p i and normalized error δ i Trajectory integral correction normalization error δ i Simultaneously, the pan angle parameterization model is updated and corrected to adjust the pan angle amplitude |σ| under the current flight condition. i (e i )|, until the normalized error δ i =0;

[0011] S3, based on flight time error △t togo,i Design heading angle error corridor ψ side,i (s i Perform trajectory integration and correct the heading angle error corridor ψ side,i (s i The sign of the roll angle is reversed based on the correction value until the time-of-flight error Δt is reached. togo,i =0;

[0012] S4. Design a gas parameter perturbation model, based on the unscented Kalman filter algorithm for the unknown perturbation parameter C of the drag coefficient. Dd,i The unknown perturbation parameter C of the lift coefficient Ld,i Accurately identify and correct the aerodynamic parameters of the aircraft during the trajectory prediction integration process in steps S2 and S3; simultaneously, address the heading angle error corridor ψ. side,i (s i Compensation will be provided.

[0013] S5. Repeat steps S2-S4 until the hypersonic vehicle reaches the predetermined target.

[0014] Furthermore, the dimensionless three-degree-of-freedom motion equations described in step S1, with energy as the independent variable, are as follows:

[0015]

[0016] Among them, R i V represents the dimensionless geocentric distance of the aircraft. i θ represents the dimensionless velocity of the aircraft. i It is the longitude of the aircraft, φ i It is the latitude of the aircraft, γ i ψ represents the inclination angle of the aircraft's trajectory. i σ represents the azimuth angle of the flight path. i D represents the yaw angle. i L represents the dimensionless drag of an aircraft. i s represents the dimensionless lift of an aircraft. i Indicates the aircraft's ready-to-fly range, e i e is an energy variable i The specific expression is:

[0017]

[0018] Furthermore, the process constraints described in step S1 are as follows:

[0019]

[0020] q i =0.5ρV i 2 g0r0≤q max

[0021]

[0022] in, q represents the current heat flux density of each aircraft. i Indicates the current dynamic pressure of each aircraft, n i Indicates the current overload of each aircraft, k Q The coefficients for the heat flux density model are determined by the aircraft structure; r0 represents the Earth's radius, g0 represents the gravitational acceleration at sea level, ρ represents the atmospheric density, and Q... max n represents the maximum allowable heat flux. max For the maximum permissible overload, q max Maximum permissible dynamic pressure;

[0023] The time constraints are as follows:

[0024] t f,i =t f

[0025] Among them, t f,i The total flight time of all aircraft is represented by t, and the coordinated arrival time of multiple aircraft is t.f , t f ∈[max{t min,1 ,…,t min·i}, min{t max,1 ,…,t max,i The shortest gliding time from the starting point to the terminal target point of the i-th spacecraft is t. min,i The longest time is t max,i ;

[0026] The endpoint constraints are as follows:

[0027]

[0028] Where, r f For the target flight altitude, v f Let θ be the target flight speed. f Indicates the target longitude. Represents the target latitude, r(t) f,i ) indicates that each aircraft is at t f,i The height at time, v(t) f,i ) indicates that each aircraft is at t f,i The velocity at time t, θ(t) f,i ) indicates that each aircraft is at t f,i Longitude of time Indicates that each aircraft is at t f,i The latitude of time.

[0029] Furthermore, in step S2, the normalization error δ i Regarding the time error Δt togo,i and landing point error Δs i As shown in the following formula:

[0030] δ i =(1-p i )Δs i +p i Δt togo,i

[0031] in,

[0032]

[0033]

[0034] In the formula, p i For the design of normalized error coefficient weights, e 0,i e represents the current energy of the spacecraft. f,i This indicates the terminal energy of the aircraft.

[0035] Furthermore, the error coefficient weight designed in step S2 is p.i The expression is as follows:

[0036]

[0037] In the formula, s i This represents the aircraft's ready-to-fly range, specifically the length of the great circle between the aircraft's real-time nadir position and the target position, expressed in seconds (s). set To pre-design the switching parameters;

[0038] The tilt angle parameterization model is designed for a specified energy level as follows;

[0039] |σ i (e i )|=A i (e i -B i ) 2 +C i

[0040] In the formula, C i For the parameter to be corrected, A i and B i By C i Uniquely determined, |σ i (e i )| indicates the specified energy e i The tilt angle amplitude at that location.

[0041] Furthermore, in step S2, the specific process for correcting the normalization error is as follows:

[0042] Combining the three-degree-of-freedom motion equations and the tilt angle parameterization model of the aircraft, from the current flight state of the aircraft to the point where the energy reaches the terminal energy e... f,i By performing trajectory prediction integration, the normalized error δ is obtained. i It is the yaw angle parameter C i Functions:

[0043] δ i =δ i (C i )

[0044] Find C using the iterative method i Size, and update the corresponding tilt angle magnitude;

[0045]

[0046] In the formula, C j,i C represents the tilt angle parameter value in the j-th iteration. j+1,i C represents the tilt angle parameter value in the (j+1)th iteration. j-1,i This represents the tilt angle parameter value at the (j-1)th iteration.

[0047] Furthermore, the heading angle error Δψ in step S3... i Let ψ be the line-of-sight azimuth angle of the target from the current position of the aircraft. LOS,i Angle with the current direction of velocity in the horizontal plane:

[0048] △ψ i =ψ i -ψ LOS,i

[0049] In the formula, ψ LOS,i This indicates the line-of-sight azimuth angle between the aircraft's current position and the target position;

[0050] The aforementioned heading angle error corridor ψ side,i (s i The flight range s of the aircraft i The changes are as follows:

[0051]

[0052] Among them, K i ψ1 is the corridor width coefficient, s1 is the design of the near-distance boundary condition of the heading angle error corridor, s2 is the design of the far-distance boundary condition of the heading angle error corridor, ψ1 is the width of the near-distance heading angle error corridor, and ψ2 is the width of the far-distance heading angle error corridor.

[0053] Furthermore, step S3 corrects the heading angle error corridor ψ. side,i (s i The details are as follows:

[0054] Combining the equation of motion and the tilt angle amplitude coefficient A i B i C i The trajectory prediction integration for the reentry flight path is completed; the trajectory prediction integration process terminates when the spacecraft enters the terminal region, and the flight time error Δt is output. togo,i Then the time error Δt togo,i It is the corridor boundary width coefficient K i The function, that is:

[0055] △t togo,i =△t togo,i (K i )

[0056] Based on the difference between the predicted and expected terminal times, the boundary width coefficient K of the corridor is adjusted using an iterative method. i And correct the corresponding heading angle error corridor ψ side,i (s i ):

[0057]

[0058] In the formula, K j,i K represents the width coefficient for the j-th iteration. j+1,i K represents the width coefficient of the (j+1)th iteration. j-1,i This represents the width coefficient of the (j-1)th iteration.

[0059] During reentry, if the heading angle error exceeds the corridor boundary, the aircraft will perform a roll angle reversal maneuver. The reversal logic is as follows:

[0060]

[0061] In the formula, sgn(σ 0,i ) represents the sign of the aircraft's current tilt angle, presgn(σ) 0,i The symbol () represents the tilt angle of the aircraft at a given moment.

[0062] Furthermore, in step S4, for the trajectory prediction integration process of S2 and S3, based on the uncertainty of aerodynamic parameters, a perturbation model of the reentry vehicle's aerodynamic parameters is designed:

[0063] C Df,i =C Do,i +C Dd,i

[0064] C Lf,i =C Lo,i +C Ld,i

[0065] In the formula, C Df,i The drag coefficient after disturbance, C Lf,i C is the lift coefficient after disturbance. Do,i For standard drag coefficient, C Lo,i C is the standard lift coefficient. Dd,i For the unknown disturbance parameter of the drag coefficient, C Ld,i Unknown perturbation parameters of the lift coefficient;

[0066] Introduce the unknown aerodynamic parameter disturbance into the state variable, and denote the augmented state variable x. i for:

[0067] x i =[R i θ i φ i γ i ψ i V i C Dd,i C Ld,i λ D,i λ L,i ] T

[0068] In the formula, the state variable is R. i V i θ i φ i γ i ψ i , λ D,i , λ L,i Relevant variables for unknown disturbance parameters in the design;

[0069] Therefore, the state equation of the augmented system is:

[0070]

[0071] In the formula, ω is the system noise vector, and f i (x i ,σ i ) represents the dynamic model of the system in S1.

[0072] Reentry vehicle measures axial acceleration D caused by aerodynamic forces b and normal acceleration L b The measurement equation is:

[0073] y i =h i (x i )+v i +η i

[0074] In the formula, y i For accelerometer measurement data, h i (x i ) represents the accelerometer measurement function of the aircraft, determined by the aircraft model, v i For the constant drift of the accelerometer, η i The random drift of the accelerometer follows a Gaussian distribution;

[0075] Based on the system state equation and measurement equation, the unscented Kalman filter algorithm is used to analyze the unknown disturbance parameter C. Dd,i and C Ld,i Accurate identification is used to estimate aerodynamic parameters and correct the aerodynamic parameter perturbation model of the aircraft.

[0076] Furthermore, in step S4, the compensation coefficient K is designed. L / D For the trajectory direction angle error corridor ψ side,i (s i Compensation will be provided.

[0077] ψ side,i (s i ) L / D =ψ side,i (si )*K L / D

[0078]

[0079] In the formula, ψ side,i (s i ) L / D This indicates the revised re-entry corridor.

[0080] The present invention has the following beneficial effects:

[0081] 1. By using a dimensionless energy equation, guidance accuracy is improved and computational pressure is reduced. By using a piecewise function to transform the cooperative time constraint into a normalized error constraint and introducing it into the roll angle amplitude update process, the flight time of the aircraft is controlled without affecting the guidance accuracy. At the same time, the time constraint is transformed into a heading angle error corridor width constraint, which further improves the range and accuracy of flight time control during the reentry process, resulting in strong controllability.

[0082] 2. By introducing unknown disturbances into the hypersonic reentry vehicle model and augmenting the system state, the system uncertainty is estimated and compensated using the unscented Kalman filter method, which improves the robustness and accuracy of the cooperative guidance system. Attached Figure Description

[0083] Figure 1 This is a flowchart of the time-coordinated hypersonic vehicle reentry robust guidance method of the present invention;

[0084] Figure 2 A diagram showing the landing point location of a cooperative guidance mission using the guidance method of this invention under disturbed conditions;

[0085] Figure 3 Flight trajectory diagram of a cooperative guidance mission using the guidance method of this invention under disturbed environment;

[0086] Figure 4 A diagram showing the tilt angle deflection of a cooperative guidance mission using the guidance method of this invention under disturbed conditions;

[0087] Figure 5 Heat flux density-time diagram for a cooperative guidance mission using the guidance method of this invention under disturbed environment;

[0088] Figure 6 The overload-time diagram is shown for a cooperative guidance mission using the guidance method of this invention under a disturbed environment. Detailed Implementation

[0089] The present invention will be further illustrated below with reference to the accompanying drawings and specific embodiments. It should be understood that the following specific embodiments are for illustrative purposes only and are not intended to limit the scope of the invention. After reading the present invention, any modifications of the present invention in various equivalent forms by those skilled in the art will fall within the scope defined by the appended claims.

[0090] Please see Figure 1 As shown, this invention discloses a time-coordinated reentry robust guidance method for multiple hypersonic vehicles, comprising the following steps:

[0091] S1. Establish the dimensionless three-degree-of-freedom motion equations based on energy for multiple hypersonic reentry vehicles, and consider endpoint constraints, process constraints, and time constraints to establish a model for the cooperative guidance problem of multiple hypersonic vehicle reentry.

[0092] The modeling process for the cooperative guidance problem of multiple hypersonic vehicles reentry is as follows:

[0093] 1) The dimensionless three-degree-of-freedom equation of motion with energy as the independent variable is as follows: The dimensionless three-degree-of-freedom equation of motion for the i-th member of a hypersonic reentry vehicle is:

[0094]

[0095] In the formula, R i V represents the dimensionless geocentric distance of the aircraft. i θ represents the dimensionless velocity of the aircraft. i It is the longitude of the aircraft, φ i It is the latitude of the aircraft, γ i ψ represents the inclination angle of the aircraft's trajectory. i σ represents the azimuth angle of the flight path. i D represents the yaw angle. i L represents the dimensionless drag of an aircraft. i s represents the dimensionless lift of an aircraft. i Indicates the aircraft's ready-to-fly range, e i e is an energy variable i The specific expression is:

[0096]

[0097] 2) The process constraints are as follows:

[0098]

[0099] In the formula, q represents the current heat flux density of each aircraft. i Indicates the current dynamic pressure of each aircraft, n i Indicates the current overload of each aircraft, k QThe coefficients for the heat flux density model are determined by the aircraft structure; r0 represents the Earth's radius, g0 represents the gravitational acceleration at sea level, ρ represents the atmospheric density, and Q represents the thermal flux density. max n represents the maximum allowable heat flux. max For the maximum permissible overload, q max Maximum permissible dynamic pressure;

[0100] 3) The time constraints are as follows:

[0101] t f,i =t f (3)

[0102] In the formula, t f,i The total flight time of all aircraft is represented by t, and the coordinated arrival time of multiple aircraft is t. f , t f ∈[max{t min,1 ,…,t min·i}, min{t max,1 ,…,t max,i The shortest gliding time from the starting point to the terminal target point of the i-th spacecraft is t. min,i The longest time is t max,i ;

[0103] 4) The endpoint constraints are as follows:

[0104]

[0105] Where, r f For the target flight altitude, v f Let θ be the target flight speed. f Indicates the target longitude. Represents the target latitude, r(t) f,i ) indicates that each aircraft is at t f,i The height at time, v(t) f,i ) indicates that each aircraft is at t f,i The velocity at time t, θ(t) f,i ) indicates that each aircraft is at t f,i Longitude of time Indicates that each aircraft is at t f,i The latitude of time.

[0106] S2. In each guidance cycle, the design error coefficient weight p i and normalized error δ i Trajectory integral correction normalization error δ i Simultaneously, the pan angle parameterization model is updated and corrected to adjust the pan angle amplitude |σ| under the current flight condition. i (e i )|, until the normalized error δ i =0;

[0107] in,

[0108] 1) Normalization error δ i As shown in the following formula:

[0109] δ i =(1-p i )Δs i +p i Δt togo,i (5)

[0110] In the formula, the error coefficient weight is p i The time error is Δt togo,i The landing point error is Δs i

[0111] Among them, the designed error coefficient weight p i Time error Δt togo,i Landing point error Δs i The calculation process is as follows:

[0112] Design error coefficient weight p i The expression is as follows:

[0113]

[0114] In the formula, p i For the design of normalized error coefficient weights, s i This represents the aircraft's ready-to-fly range, specifically the length of the great circle between the aircraft's real-time nadir position and the target position, expressed in seconds (s). set To pre-design the switching parameters;

[0115] If we disregard the effect of Earth's rotation, the derivative of e with respect to dimensionless time is:

[0116]

[0117] The derivative with respect to energy is:

[0118]

[0119] The flight time error Δt can be obtained by integration. togo,i :

[0120]

[0121] In the formula, e 0,i e represents the current energy of the spacecraft. f,i This indicates the terminal energy of the aircraft.

[0122] Similarly, the landing point error is also known as the pre-flight range error Δs. i Obtained by integration

[0123]

[0124] 2) The parameterized model for the tilt angle is as follows;

[0125] |σ i (e i )|=A i (e i -B i ) 2 +C i (9)

[0126] In the formula, C i For the yaw angle parameter to be corrected, A i and B i By C i The only certainty.

[0127]

[0128] Where, σ 0,i The current energy e to be corrected 0,i The amplitude of the tilt angle at that point, σ f,i Represents the terminal energy e f,i The tilt angle amplitude.

[0129] 3) The specific process for correcting the normalization error is as follows:

[0130] Combining the three-degree-of-freedom motion equations (1) and the tilt angle parameterization model (9) of the aircraft, from the current flight state of the aircraft to the energy reaching the terminal energy e f,i By performing trajectory prediction integration, the normalized error δ is obtained. i It is the yaw angle parameter C i Functions:

[0131] δ i =δ i (C i )

[0132] Find C using the iterative method i size;

[0133]

[0134] In the formula, C j,i C represents the tilt angle parameter value in the j-th iteration. j+1,i C represents the tilt angle parameter value in the (j+1)th iteration. j-1,i This represents the tilt angle parameter value at the (j-1)th iteration.

[0135] 4) Then, the C obtained by the iterative method i Calculate Ai and B i And substitute it into (9) to obtain the tilt angle amplitude |σ i (e i )|.

[0136] S3, based on flight time error △t togo,i Design heading angle error corridor ψ side,i (s i Perform trajectory integration and correct the heading angle error corridor ψ side,i (s i The sign of the roll angle is reversed based on the correction value until the time-of-flight error Δt is reached. togo,i =0;

[0137] in,

[0138] 1) Heading angle error △ψ i The angle between the line-of-sight azimuth of the target from the aircraft's current position and the velocity direction in the current horizontal plane:

[0139] △ψ i =ψ i -ψ LOS,i

[0140] In the formula, ψ LOS,i The azimuth angle between the aircraft's current position and the target position is expressed by the following formula:

[0141]

[0142] Design the heading angle error corridor ψ that varies with the aircraft's waiting flight range using piecewise functions. side,i (s i )as follows:

[0143]

[0144] Among them, K i Let ψ1 be the corridor width coefficient, s1 be the near-range boundary condition of the designed heading angle error corridor, s2 be the far-range boundary condition of the designed heading angle error corridor, ψ1 be the width of the near-range heading angle error corridor, and ψ2 be the width of the far-range heading angle error corridor. In this embodiment, the heading angle error corridor widths ψ1 and ψ2 are 10° and 20° respectively, the boundary s1 is 10°, and the boundary s2 is 40°.

[0145] 2) Correction of heading angle error corridor ψ side,i (s i The specific process is as follows:

[0146] Combining the equation of motion (1) and the tilt angle amplitude coefficient A i B i C iThe trajectory prediction integration for the reentry flight path is completed; the trajectory prediction integration process terminates when the spacecraft enters the terminal region, and the flight time Δt is output. togo,i Then the flight time error Δt togo,i It is the corridor boundary width coefficient K i The function, that is:

[0147] △t togo,i =△t togo,i (K i )

[0148] Based on the difference between the predicted and expected terminal times, the boundary width coefficient K of the corridor is adjusted using an iterative method. i And update the heading angle error corridor ψ accordingly. side,i (s i ):

[0149]

[0150] In the formula, K j,i K represents the width coefficient for the j-th iteration. j+1,i K represents the width coefficient of the (j+1)th iteration. j-1,i This represents the width coefficient of the (j-1)th iteration.

[0151] Next, the calculated boundary width coefficient K of the corridor will be... i Substitute (10) to update the heading angle error corridor ψ side,i (s i ).

[0152] 3) During reentry, when the heading angle error Δψ i If the aircraft exceeds the corridor boundary, it will perform a roll angle reversal maneuver. The reversal logic is as follows, which facilitates the aircraft in executing signed roll angle commands:

[0153]

[0154] In the formula, sgn(σ 0,i ) represents the sign of the aircraft's current tilt angle, presgn(σ) 0,i The symbol () represents the tilt angle of the aircraft at a given moment.

[0155] S4. Design a gas parameter perturbation model, based on the unscented Kalman filter algorithm for the unknown perturbation parameter C of the drag coefficient. Dd,i The unknown perturbation parameter C of the lift coefficient Ld,i Accurately identify and correct the aerodynamic parameters of the aircraft during the trajectory prediction integration process in steps S2 and S3; simultaneously, address the heading angle error corridor ψ. side,i (s i Compensation will be provided.

[0156] 1) For the trajectory prediction integration process of S2 and S3, based on the uncertainty of aerodynamic parameters, design an aerodynamic parameter disturbance model for the reentry vehicle:

[0157]

[0158] In the formula, C Df,i The drag coefficient after disturbance, C Lf,i C is the lift coefficient after disturbance. Do,i For standard drag coefficient, C Lo,i C is the standard lift coefficient. Dd,i For the unknown disturbance parameter of the drag coefficient, C Ld,i Unknown perturbation parameters of the lift coefficient.

[0159] 2) Place C Dd C Ld The motions considered as standard data are described using a first-order Gaussian Markov process, specifically in the following form:

[0160]

[0161] In the formula, It is a Gaussian white noise sequence.

[0162] Introduce the unknown aerodynamic parameter disturbance into the state variable, and denote the augmented state variable x. i for:

[0163] x i =[R i θ i φ i γ i ψ i V i C Dd,i C Ld,i λ D,i λ L,i ] T

[0164] In the formula, the state variables include R i V i θ i φ i γ i ψ i , λ D,i and λ L,i To design the relevant variables for unknown disturbance parameters.

[0165] Therefore, the state equation of the augmented system is:

[0166]

[0167] In the formula, ω is the system noise vector, and f i (x i ,σ i ) represents the dynamic model of the system.

[0168] Reentry vehicle measures axial acceleration D caused by aerodynamic forces b and normal acceleration L b The measurement equation is:

[0169] y i =h i (x i )+v i +η i

[0170] In the formula, y i For accelerometer measurement data, h i (x i ) represents the accelerometer measurement function of the aircraft, determined by the aircraft model, v i For the constant drift of the accelerometer, η i The accelerometer's random drift follows a Gaussian distribution.

[0171] In this embodiment, the covariance matrix P is taken as...

[0172]

[0173] Based on the system state equation (12) and measurement equation (13), the unknown disturbance parameter C of the drag coefficient is obtained using the unscented Kalman filter algorithm. Dd,i The unknown perturbation parameter C of the lift coefficient Ld,i Accurate identification is used to estimate aerodynamic parameters, and the aerodynamic parameter disturbance model in the aircraft trajectory prediction integration process is corrected (11).

[0174] 4) Design compensation coefficient K L / D For the heading angle error corridor ψ side,i (s i Compensation will be provided.

[0175] ψ side,i (s i ) L / D =ψ side,i (s i )*K L / D

[0176]

[0177] In the formula, ψ side,i (s i ) L / D This indicates the revised re-entry corridor.

[0178] S5. Repeat steps S2-S4 until the hypersonic vehicle reaches the predetermined target.

[0179] This embodiment uses a simulation method based on the aforementioned guidance method to verify the feasibility of the guidance method. In this embodiment, three CAV-H spacecraft with different initial positions and a cooperative reentry mission are used. The heat flux density constraint of the CAV-H spacecraft is 6000 Kw / m^2, the dynamic pressure constraint is 40 kPa, and the overload constraint is 4. The cooperative flight time t is set. f The timeframe is 1435s. The target coordinates for the gliding phase are set to (105°, 5°, 20km), and the terminal velocity is greater than 1500m / s. The initial states for the reentry phase are selected as shown in the table below.

[0180]

[0181]

[0182] To compare the robustness of the proposed guidance method, a total of 100 simulations were performed. The boundary conditions under the Gaussian distribution of the variable perturbation in the Monte Carlo simulation results are shown in the table below.

[0183]

[0184] In both tables, the unit for aircraft altitude is kilometers, the unit for speed is meters per second, and the unit for other angles is degrees. Specifically, the state variables H represent flight altitude, θ represents the longitude of the nadir point, Φ represents latitude, V represents the aircraft speed, Γ represents the track inclination angle, Ψ represents the track azimuth angle, m represents mass in kilograms, cl represents the lift coefficient, and cd represents the drag coefficient.

[0185] Please see Figures 2 to 6 As shown, under this guidance method, the ground projection curves of the reentry terminal landing point and reentry trajectory, the vehicle's roll angle curve, and the heat flux density constraint and overload constraint diagrams are obtained. The vehicle's roll angle curve shows that there are not many reversals, meeting engineering requirements. In the heat flux density constraint and overload constraint diagrams, the average time error calculated from the flight time is -0.66s, meeting the terminal time constraint requirements of the cooperative attack mission and also satisfying the vehicle performance constraints. These results demonstrate that the method designed in this invention can effectively control the reentry terminal time for cooperative missions of hypersonic vehicles.

[0186] In summary, this invention provides a time-coordinated reentry robust guidance method for multiple hypersonic vehicles, which features high robustness and strong controllability.

Claims

1. A time-coordinated multi-hypersonic vehicle re-entry robust guidance method, characterized in that, Includes the following steps: S1. Establish the dimensionless three-degree-of-freedom motion equations based on energy for multiple hypersonic reentry vehicles, considering endpoint constraints, process constraints, and time constraints, and establish a dynamic model for the reentry cooperative guidance problem of multiple hypersonic vehicles. S2. In each guidance cycle, the design error coefficient weighting and normalized error Trajectory integral correction normalization error Meanwhile, the roll angle parameterization model is updated and corrected to adjust the roll angle amplitude under the current flight condition. until the normalization error The normalization error is 0. Regarding time error and landing point error As shown in the following formula: in, In the formula, For the design of normalized error coefficient weights, Indicates the current energy of the aircraft. Indicates the terminal energy of the aircraft. This represents the dimensionless drag of an aircraft. This represents the dimensionless velocity of the aircraft. For energy variables, Indicates the tilt angle of the aircraft's trajectory. The dimensionless geocentric distance of the aircraft; The error coefficient weights of the design are The expression is as follows: In the formula, This represents the aircraft's ready-to-fly range, which is the length of the great circle between the aircraft's real-time nadir position and the target position. To pre-design the switching parameters; Indicates the coordinated arrival time of multiple aircraft. Indicates time-of-flight error; The tilt angle parameterization model is designed for a specified energy level as follows; In the formula, , and This is the tilt angle amplitude coefficient. The parameters to be corrected are... and Depend on The only certainty, Indicates specified energy The tilt angle amplitude at that location; S3, based on flight time error Design heading angle error corridor Perform trajectory integration and correct heading angle error corridor The sign of the roll angle is reversed based on the correction value until the time-of-flight error is reached. =0; S4. Design a gas parameter perturbation model, based on the unscented Kalman filter algorithm for the unknown perturbation parameters of the drag coefficient. Unknown perturbation parameters of lift coefficient Accurately identify and correct the aerodynamic parameters of the aircraft during the trajectory prediction integration process in steps S2 and S3; simultaneously, address the heading angle error corridor. Provide compensation; S5. Repeat steps S2-S4 until the hypersonic vehicle reaches the predetermined target.

2. The robust guidance method according to claim 1, characterized in that, The dimensionless three-degree-of-freedom equation of motion described in step S1 uses energy as the independent variable, and the third degree of motion in a hypersonic reentry vehicle... The dimensionless three-degree-of-freedom equations of motion for each member are: in, This represents the dimensionless geocentric distance of the aircraft. This represents the dimensionless velocity of the aircraft. It is the longitude of the aircraft. It is the latitude of the aircraft. Indicates the tilt angle of the aircraft's trajectory. Indicates the azimuth angle of the flight path. Indicates the tilt angle. This represents the dimensionless drag of an aircraft. This represents the dimensionless lift of an aircraft. Indicates the aircraft's ready-to-fly range. For energy variables, The specific expression is: 。 3. The robust guidance method according to claim 2, characterized in that, The process constraints described in step S1 are as follows: in, This indicates the current heat flux density of each aircraft. This indicates the current dynamic pressure of each aircraft. This indicates that each aircraft is currently overloaded. These are the coefficients for the heat flux density model, determined by the aircraft structure. Represents the radius of the Earth. This represents the gravitational acceleration at sea level. Indicates atmospheric density. Indicates the maximum allowable heat flux. For the maximum permissible overload, Maximum permissible dynamic pressure; The time constraints are as follows: in, This represents the total flight time of all aircraft, and the coordinated arrival time of multiple aircraft is... , , No. The shortest gliding time from the starting point to the terminal target point for each aircraft is: The longest time is ; The endpoint constraints are as follows: in, For the target flight altitude, For the target flight speed, Indicates the target longitude. Indicates the target latitude. Indicates that each aircraft is in The height of time, Indicates that each aircraft is in The speed of time, Indicates that each aircraft is in Longitude of time Indicates that each aircraft is in The latitude of time.

4. The robust guidance method according to claim 3, characterized in that, In step S2, the specific process for correcting the normalization error is as follows: Combining the three-degree-of-freedom motion equations and the parameterized model of the tilt angle of the aircraft, from the current flight state of the aircraft to the point where the energy reaches the terminal energy... Perform trajectory prediction integration to obtain the normalized error. It is the yaw angle parameter. Functions: Using iterative methods to find Size, and update the corresponding tilt angle magnitude; In the formula, Indicates the first The tilt angle parameter value is determined in the next iteration. Indicates the middle The tilt angle parameter value of the next iteration Indicates the first The tilt angle parameter value at the next iteration.

5. The robust guidance method according to claim 4, characterized in that, The heading angle error in step S3 The azimuth angle of the target from the current position of the aircraft. Angle with the current direction of velocity in the horizontal plane: In the formula, This indicates the line-of-sight azimuth angle between the aircraft's current position and the target position; The aforementioned heading angle error corridor The flight range of the aircraft The changes are as follows: in, This is the corridor width coefficient. For the design of the heading angle error corridor, near-range boundary conditions, For the long-distance boundary conditions of the designed heading angle error corridor, The width of the near-range heading angle error corridor. The width of the long-distance heading angle error corridor.

6. The robust guidance method according to claim 5, characterized in that, The step S3 described above corrects the heading angle error corridor. Specifically as follows: Combining the equations of motion and the tilt angle magnitude coefficient , , The trajectory prediction integration for the reentry flight path is completed; the trajectory prediction integration process terminates when the spacecraft enters the terminal region, and the flight time error is output. Then time error It is the corridor boundary width coefficient. The function, that is: Based on the difference between the predicted and expected terminal times, an iterative method is used to adjust the boundary width coefficient of the corridor. And correct the corresponding heading angle error corridor : In the formula, The width coefficient represents the number of iterations in the j-th iteration. Indicates the first Width coefficient of the next iteration Indicates the first The width coefficient of the next iteration; During reentry, if the heading angle error exceeds the corridor boundary, the aircraft will perform a roll angle reversal maneuver. The reversal logic is as follows: In the formula, The symbol representing the aircraft's current tilt angle. The symbol representing the tilt angle of the aircraft at a given moment.

7. The robust guidance method according to claim 6, characterized in that, In step S4, for the trajectory prediction integration process of S2 and S3, based on the uncertainty of aerodynamic parameters, a perturbation model of the reentry vehicle's aerodynamic parameters is designed: In the formula, For the drag coefficient after disturbance, The lift coefficient after disturbance. For standard drag coefficient, The standard lift coefficient, For unknown disturbance parameters of drag coefficient, The unknown perturbation parameter is the lift coefficient; Introducing the unknown aerodynamic parameter disturbance into the state variables, we denote the augmented state variables. for: In the formula, the state variable is , , , , , , , Relevant variables for unknown disturbance parameters in the design; Therefore, the state equation of the augmented system is: In the formula, Let be the system noise vector. This represents the dynamic model of the system in S1; Reentry vehicle measures axial acceleration caused by aerodynamic forces and normal acceleration The measurement equation is: In the formula, For accelerometer measurement data, The measurement function for the aircraft accelerometer is determined by the aircraft model. For constant drift of the accelerometer, The random drift of the accelerometer follows a Gaussian distribution; Based on the system state equation and measurement equation, the unscented Kalman filter algorithm is used to analyze the unknown disturbance parameters. and Accurate identification is used to estimate aerodynamic parameters and correct the aerodynamic parameter perturbation model of the aircraft.

8. The robust guidance method according to claim 7, characterized in that, In step S4, Design compensation coefficient For the trajectory direction angle error corridor Provide compensation: In the formula, This indicates the revised re-entry corridor.

Images