Gliding aircraft chaotic state excitation method based on Lorenz system

Abstract

The invention discloses a chaotic state excitation method for a gliding aircraft based on a Lorenz system, and belongs to the technical field of aircraft control. Selecting a nominal trajectory, and recording a change sequence of state variables of the nominal trajectory along with time; a complete Lorenz system differential equation is obtained; the expected three-axis speed increment of the aircraft reentry system is obtained through calculation; generating a transverse additional overload and a longitudinal additional overload of the aircraft under the chaotic trajectory reentry system; obtaining a flight state variable of the next step length; and repeating until the chaotic maneuvering flight is finished, and completing the chaotic flight excitation for the aircraft. According to the method, the maneuvering instruction is designed by using the Lorenz system, and the unpredictability of the gliding section trajectory of the aircraft is improved. The generated chaotic trajectory automatically meets constraint conditions such as heat flow, overload, dynamic pressure and the like, and the overall controllability is ensured. The maneuvering instruction is generated based on the real-time state and deviation, dependence on enemy aircraft information is weak, online instruction generation is supported, and tactical flexibility is enhanced.

Description

Technical Field

[0001] The invention relates to a gliding aircraft chaotic state excitation method based on a Lorenz system, belonging to the technical field of aircraft control. Background Art

[0002] Hypersonic glide vehicles, with their high speed and high maneuverability, have broad application prospects in the strategic and tactical fields. However, with the rapid development of modern interception technology, especially the continuous improvement of full-dimensional and multi-level near-space detection systems, the flight trajectories of such vehicles are facing increasingly severe challenges in identification and tracking.

[0003] The flight trajectory of traditional gliders often has a certain regularity and predictability, which allows modern detection systems to quickly lock and track their flight paths and then effectively intercept them. Therefore, the survivability and penetration effectiveness of gliders are greatly threatened.

[0004] In order to meet this challenge, an innovative method is urgently needed to break the predictability of the flight trajectory of the glider and improve its concealment and penetration capability. Based on this, the present invention proposes a method for exciting the chaotic state of a glider based on the Lorenz system. The Lorenz system is a typical chaotic system, and its motion trajectory is highly complex and unpredictable. By introducing the chaotic characteristics of the Lorenz system into the flight control of the glider, its flight trajectory can be given chaotic characteristics, making it difficult to be accurately predicted and tracked by the detection system.

[0005] The invention aims to improve the penetration efficiency of gliders and enhance their survivability in modern warfare environments by stimulating the chaotic state. This method not only has important theoretical value, but also has significant engineering practical significance, and is expected to open up new ideas for the design and application of gliders. Summary of the invention

[0006] In order to solve the problems existing in the background technology, the present invention provides a method for exciting the chaotic state of a gliding aircraft based on a Lorenz system.

[0007] To achieve the above object, the present invention adopts the following technical solution: a method for exciting the chaotic state of a gliding aircraft based on a Lorenz system, the method comprising the following steps:

[0008] S1: Select a reentry glide trajectory as the nominal trajectory for constructing chaotic trajectories, and record the sequence of changes of the state variables of the nominal trajectory over time t;

[0009] S2: According to the differential equation of the Lorenz system, the corresponding relationship between the state variables of the nominal trajectory and the state variables of the reentry dynamics of the aircraft is established, and the parameters of the Lorenz system are set to obtain the complete differential equation of the Lorenz system;

[0010] The complete Lorenz system differential equation described in S2 is as follows:

[0011]

[0012] In formula (1):

[0013] [x, y, z] T is the three-dimensional state variable of the Lorenz system;

[0014] is the differential value corresponding to the three-dimensional state variable;

[0015] Pr>0 is the Prandtl number, which is related to the flow field;

[0016] Ra>0 is the Rayleigh number, which is related to the flow field;

[0017] β>0 is a spatial parameter that reflects the aspect ratio of the volume in the convection problem.

[0018] S3: Based on the differential equations of the complete Lorenz system, the three-axis position deviation of the vehicle in the vehicle reentry system relative to the nominal trajectory is used as input to calculate the desired three-axis velocity increment of the vehicle reentry system;

[0019] The calculation process of the three-axis velocity increment of the desired reentry system of the aircraft described in S3 is as follows:

[0020] S301: Taking the three-axis position deviation of the spacecraft reentry system as input and the desired three-axis velocity increment of the spacecraft reentry system as output, the complete Lorenz system differential equation is rewritten as follows:

[0021]

[0022] In formula (2):

[0023] [Δx, Δy, Δz] T is the three-axis position deviation of the spacecraft reentry system;

[0024] ΔV=[ΔV x , ΔV y , ΔV z ] T is the three-axis velocity increment of the desired spacecraft reentry system;

[0025] S302: Introducing range control parameters [R x , Ry , R z ] T The new differential equation obtained by combining the chaos construction speed control parameter k is as follows:

[0026]

[0027] S303: Set the parameters of the Lorenz system according to the maneuvering range requirements, maneuvering time requirements and overload constraints

[0028] Pr=10,Ra=28,β=8 / 3,R x =5000,R z =10000,R y =5000, k=1100, the final input-output relationship is as follows:

[0029]

[0030] S304: The three-axis position deviation [Δx, Δy, Δz] of the aircraft reentry system calculated at each simulation step T , that is, to obtain the desired three-axis velocity increment of the reentry system of the spacecraft ΔV = [ΔV x ,ΔV y ,ΔV z ] T

[0031] S4: Combine the state variables of the nominal trajectory at the corresponding moment and the expected velocity increment to calculate the expected three-axis velocity of the chaotic ballistic reentry system. Combine the actual velocity in the current state to generate the lateral additional overload and longitudinal additional overload of the aircraft in the chaotic ballistic reentry system by the PD tracking method.

[0032] The calculation process of the three-axis expected velocity of the chaotic ballistic reentry system described in S4 is as follows:

[0033] S401: Define the current nominal trajectory state variables as [θ, φ, h, V, γ, ψ] T 0 , where θ is longitude, φ is latitude, h is altitude, V is velocity, γ is track angle, and ψ is heading angle. The velocity vector in the North Celestial East coordinate system is:

[0034]

[0035] S402: Convert the velocity vector in the North Celestial East coordinate system to the aircraft reentry system and express it as follows:

[0036]

[0037] In formula (6):

[0038] is the attitude transformation matrix of the spacecraft reentry system relative to the Earth-centered fixed coordinate system;

[0039] is the attitude transformation matrix of the geocentric fixed coordinate system relative to the North Celestial East coordinate system;

[0040] is the velocity component of the nominal trajectory in the reentry system at the current moment;

[0041] S403: Similarly to S401-S402, define the state variables of the current chaotic flight as [θ, φ, h, V, γ, ψ] T , get the velocity component of chaotic flight in the reentry system at the current moment

[0042] S404: The velocity component of the nominal trajectory at the current moment in the reentry system and the desired three-axis velocity increment of the reentry system of the spacecraft ΔV = [ΔV x , ΔV y , ΔV z ] T Obtain the expected velocity component of chaotic flight:

[0043]

[0044] S405: Expected velocity component from chaotic flight and the velocity component of the chaotic flight in the reentry system at the current moment The additional lateral overload Δn is obtained by PD tracking method z and additional longitudinal overload Δn y :

[0045]

[0046] In formula (8):

[0047] K p is the proportionality coefficient;

[0048] K d is the differential coefficient;

[0049] is the y-axis velocity component of the reentry system during the actual flight at the previous simulation step;

[0050] is the expected y-axis velocity component of the reentry system at the previous simulation step;

[0051] is the z-axis velocity component of the reentry system during the actual flight at the previous simulation step;

[0052] is the expected z-axis velocity component of the reentry system at the previous simulation step.

[0053] S5: Combine the nominal ballistic maneuvering overload and additional overload to solve the required maneuvering overload of chaotic maneuvering flight, inversely solve the corresponding angle of attack and roll angle, substitute them into the reentry glide dynamics model, and obtain the next step of flight state variables through the fourth-order Runge-Kutta integration method;

[0054] S501: Calculation of lateral overload n from nominal ballistic state variables z0 and longitudinal overload n y0 , combined with additional lateral overload Δn z and additional longitudinal overload Δn y Calculate the required lateral overload n for chaotic maneuvering flight z and required longitudinal overload n y :

[0055]

[0056] S502: Calculate the required lateral overload n of the glider z and required longitudinal overload n y :

[0057]

[0058] In formula (10):

[0059] is the aircraft dynamic pressure, where: ρ is the atmospheric density, V is the aircraft speed;

[0060] S is the reference aerodynamic area of ​​the aircraft;

[0061] C L is the lift coefficient, which is a function of the angle of attack α;

[0062] C D is the drag coefficient, which is a function of the angle of attack α;

[0063] mg is the gravity on the aircraft;

[0064] S503: According to formula (9), the required lateral overload n is known at this time z and required longitudinal overload n y , inversely solve the corresponding angle of attack α and roll angle σ through Newton iteration method;

[0065] S504: Combine the state variables of the current chaotic flight [θ, φ, h, V, γ, ψ] T , angle of attack α and roll angle σ, and calculate the differential values ​​of each state variable by the reentry glide dynamics equation:

[0066]

[0067] In formula (11):

[0068] r=R 0 +h is the distance from the center of the earth to the aircraft, where: R 0 is the radius of the Earth;

[0069] g is the acceleration due to gravity of the Earth;

[0070] m is the weight of the aircraft;

[0071] S505: Based on the differential value of each state variable The next state variable is obtained through the fourth-order Runge-Kutta integration method, and then the aircraft state update is completed.

[0072] S6: Repeat S2-S5 until the chaotic maneuvering flight ends, completing the chaotic flight excitation for the aircraft.

[0073] Compared with the prior art, the present invention has the following beneficial effects:

[0074] The present invention uses the classic Lorenz system to design maneuvering instructions for the reentry glide process of a hypersonic glide vehicle, actively stimulates the chaotic state of the vehicle, and the trajectory designed based on this method has state sensitivity, disorder, and global boundedness. Under the condition of overall controllability, the deviation generated during the flight process is used to greatly improve the unpredictability of the trajectory of the glide segment of the vehicle. At the same time, the generated chaotic trajectory is above the nominal trajectory and automatically meets the constraints of heat flow, overload, dynamic pressure, etc., ensuring overall controllability. In addition, this method generates maneuvering instructions based on real-time status and deviations, has weak dependence on enemy aircraft information, supports online generation of instructions, and enhances tactical flexibility. BRIEF DESCRIPTION OF THE DRAWINGS

[0075] Figure 1 is a flow chart of the present invention;

[0076] Figure 2 It is a schematic diagram of the three-dimensional trajectory of the complete Lorenz system;

[0077] Figure 3 It is a schematic diagram of the longitude and latitude height curves of the nominal trajectory and chaotic maneuvering trajectory;

[0078] Figure 4 It is a schematic diagram of the position curve of the reentry system of the nominal ballistic trajectory and chaotic maneuvering trajectory. DETAILED DESCRIPTION

[0079] The technical solution of the present invention will be clearly and completely described below in conjunction with the drawings in the embodiments of the present invention. Obviously, the described embodiments are only part of the embodiments of the invention, rather than all of the embodiments. Based on the embodiments of the present invention, all other embodiments obtained by ordinary technicians in this field without making creative work are within the scope of protection of the present invention.

[0080] A method for exciting a chaotic state of a gliding aircraft based on a Lorenz system, the method comprising the following steps:

[0081] S1: A reentry glide trajectory is selected as the nominal trajectory for constructing chaotic trajectories. The nominal trajectory includes three stages: active stage, glide stage and pressure stage. Its range is about 12000km and the flight time is about 1100s. The change sequence of the state variables (longitude θ, latitude φ, altitude h, velocity V, track angle γ and heading angle ψ) of the nominal trajectory with time t is recorded;

[0082] S2: According to the differential equation of the Lorenz system, the corresponding relationship between the state variables of the nominal trajectory and the state variables of the reentry dynamics of the aircraft is established. At the same time, the parameters of the Lorenz system are set according to the range, speed and other characteristics of the excited chaos to obtain the complete differential equation of the Lorenz system; according to the basic conditions of the nominal trajectory, the altitude maneuvering range of the chaotic trajectory is set to 5km, the lateral maneuvering range is set to 20km, the continuous maneuvering time is set to 500s, and the maximum maneuvering overload is set to no more than 10.

[0083] The complete Lorenz system differential equation described in S2 is as follows:

[0084]

[0085] In formula (1):

[0086] [x, y, z] T is the three-dimensional state variable of the Lorenz system;

[0087] is the differential value corresponding to the three-dimensional state variable;

[0088] Pr>0 is the Prandtl number, which is related to the flow field;

[0089] Ra>0 is the Rayleigh number, which is related to the flow field;

[0090] β>0 is a spatial parameter that reflects the aspect ratio of the volume in the convection problem.

[0091] Under this equation, the three-dimensional trajectory is formed as Figure 2 shown.

[0092] S3: Based on the differential equations of the complete Lorenz system, the three-axis position deviation of the vehicle in the vehicle reentry system relative to the nominal trajectory is used as input to calculate the desired three-axis velocity increment of the vehicle reentry system;

[0093] The calculation process of the three-axis velocity increment of the desired reentry system of the aircraft described in S3 is as follows:

[0094] S301: Taking the three-axis position deviation of the spacecraft reentry system as input and the desired three-axis velocity increment of the spacecraft reentry system as output, the complete Lorenz system differential equation is rewritten as follows:

[0095]

[0096] In formula (2):

[0097] [Δx, Δy, Δz] T is the three-axis position deviation of the spacecraft reentry system;

[0098] ΔV=[ΔV x , ΔV y , ΔV z ] T is the three-axis velocity increment of the desired spacecraft reentry system;

[0099] S302: Introducing range control parameters [R x , R y , R z ] T The new differential equation obtained by combining the chaos construction speed control parameter k is as follows:

[0100]

[0101] S303: Set the parameters of the Lorenz system according to the maneuvering range requirements, maneuvering time requirements and overload constraints

[0102] Pr=10,Ra=28,β=8 / 3,R x =5000,R z =10000,R y =5000, k=100, the final input-output relationship is as follows:

[0103]

[0104] S304: The three-axis position deviation [Δx, Δy, Δz] of the aircraft reentry system calculated at each simulation step T , that is, to obtain the desired three-axis velocity increment of the spacecraft reentry system

[0105] ΔV=[ΔV x, ΔV y ,ΔV z ] T .

[0106] S4: Combine the state variables of the nominal trajectory at the corresponding moment and the expected velocity increment to calculate the expected three-axis velocity of the chaotic ballistic reentry system. Combine the actual velocity in the current state to generate the lateral additional overload and longitudinal additional overload of the aircraft in the chaotic ballistic reentry system by the PD tracking method.

[0107] The calculation process of the three-axis expected velocity of the chaotic ballistic reentry system described in S4 is as follows:

[0108] S401: Define the current nominal trajectory state variables as [θ, φ, h, V, γ, ψ] T 0 , where θ is longitude, φ is latitude, h is altitude, V is velocity, γ is track angle, and ψ is heading angle. The velocity vector in the North Celestial East coordinate system is:

[0109]

[0110] S402: Convert the velocity vector in the North Celestial East coordinate system to the aircraft reentry system and express it as follows:

[0111]

[0112] In formula (6):

[0113] is the attitude transformation matrix of the spacecraft reentry system relative to the Earth-centered fixed coordinate system;

[0114] is the attitude transformation matrix of the geocentric fixed coordinate system relative to the North Celestial East coordinate system;

[0115] is the velocity component of the nominal trajectory in the reentry system at the current moment;

[0116] S403: Similarly to S401-S402, define the state variables of the current chaotic flight as [θ, φ, h, V, γ, ψ] T , get the velocity component of chaotic flight in the reentry system at the current moment

[0117] S404: The velocity component of the nominal trajectory at the current moment in the reentry system and the desired three-axis velocity increment of the reentry system of the spacecraft ΔV = [ΔV x , ΔV y , ΔV z ] T Obtain the expected velocity component of chaotic flight:

[0118]

[0119] S405: Expected velocity component from chaotic flight and the velocity component of the chaotic flight in the reentry system at the current moment The additional lateral overload Δn is obtained by PD tracking method z and additional longitudinal overload Δn y :

[0120]

[0121] In formula (8):

[0122] K p is the proportionality coefficient;

[0123] K d is the differential coefficient;

[0124] is the y-axis velocity component of the reentry system during the actual flight at the previous simulation step;

[0125] is the expected y-axis velocity component of the reentry system at the previous simulation step;

[0126] is the z-axis velocity component of the reentry system during the actual flight at the previous simulation step;

[0127] is the expected z-axis velocity component of the reentry system at the previous simulation step.

[0128] S5: Combine the nominal ballistic maneuvering overload and additional overload to solve the required maneuvering overload of chaotic maneuvering flight, inversely solve the corresponding angle of attack and roll angle, substitute them into the reentry glide dynamics model, and obtain the next step of flight state variables through the fourth-order Runge-Kutta integration method;

[0129] S501: Calculation of lateral overload n from nominal ballistic state variables z0 and longitudinal overload n y0 , combined with additional lateral overload Δn z and additional longitudinal overload Δn y Calculate the required lateral overload n for chaotic maneuvering flight z and required longitudinal overload n y :

[0130]

[0131] S502: Calculate the required lateral overload n of the glider z and required longitudinal overload n y :

[0132]

[0133] In formula (10):

[0134] is the aircraft dynamic pressure, where: ρ is the atmospheric density, V is the aircraft speed;

[0135] S is the reference aerodynamic area of ​​the aircraft;

[0136] C L is the lift coefficient, which is a function of the angle of attack α;

[0137] C D is the drag coefficient, which is a function of the angle of attack α;

[0138] mg is the gravity on the aircraft;

[0139] S503: According to formula (9), the required lateral overload n is known at this time z and required longitudinal overload n y , inversely solve the corresponding angle of attack α and roll angle σ through Newton iteration method;

[0140] S504: Combine the state variables of the current chaotic flight [θ, φ, h, V, γ, ψ] T , angle of attack α and roll angle σ, and calculate the differential values ​​of each state variable by the reentry glide dynamics equation:

[0141]

[0142] In formula (11):

[0143] r=R 0 +h is the distance from the center of the earth to the aircraft, where: R 0 is the radius of the Earth;

[0144] g is the acceleration due to gravity of the Earth;

[0145] m is the weight of the aircraft;

[0146] S505: Based on the differential value of each state variable The next state variable is obtained through the fourth-order Runge-Kutta integration method, and then the aircraft state update is completed.

[0147] S6: Repeat S2-S5 until the chaotic maneuvering flight ends, completing the chaotic flight excitation for the aircraft.

[0148] It will be apparent to those skilled in the art that the invention is not limited to the details of the exemplary embodiments described above and that the invention can be implemented in other forms without departing from the spirit or essential features of the invention. Therefore, the embodiments should be considered in all respects as exemplary and non-limiting, and the scope of the invention is defined by the appended claims rather than the foregoing description, and it is intended that all variations within the meaning and range of equivalents of the claims be included in the invention. Any reference numeral in a claim should not be considered as limiting the claim to which it relates.

[0149] In addition, it should be understood that although the present specification is described according to implementation modes, not every implementation mode contains only one independent technical solution. This description of the specification is only for the sake of clarity. Those skilled in the art should regard the specification as a whole. The technical solutions in each embodiment may also be appropriately combined to form other implementation modes that can be understood by those skilled in the art.

Claims

1. A method for exciting the chaotic state of a gliding aircraft based on a Lorenz system, characterized in that: The method comprises the following steps: S1: Select a reentry glide trajectory as the nominal trajectory for constructing chaotic trajectories, and record the sequence of changes of the state variables of the nominal trajectory over time t; S2: According to the differential equation of the Lorenz system, the corresponding relationship between the state variables of the nominal trajectory and the state variables of the reentry dynamics of the aircraft is established, and the parameters of the Lorenz system are set to obtain the complete differential equation of the Lorenz system; S3: Based on the differential equations of the complete Lorenz system, the three-axis position deviation of the vehicle in the vehicle reentry system relative to the nominal trajectory is used as input to calculate the desired three-axis velocity increment of the vehicle reentry system; S4: Combine the state variables of the nominal trajectory at the corresponding moment and the expected velocity increment to calculate the expected three-axis velocity of the chaotic ballistic reentry system. Combine the actual velocity in the current state to generate the lateral additional overload and longitudinal additional overload of the aircraft in the chaotic ballistic reentry system by the PD tracking method. S5: Combine the nominal ballistic maneuvering overload and additional overload to solve the required maneuvering overload of chaotic maneuvering flight, inversely solve the corresponding angle of attack and roll angle, substitute them into the reentry glide dynamics model, and obtain the next step of flight state variables through the fourth-order Runge-Kutta integration method; S6: Repeat S2-S5 until the chaotic maneuvering flight ends, completing the chaotic flight excitation for the aircraft.

2. The method for exciting the chaotic state of a gliding aircraft based on a Lorenz system according to claim 1, characterized in that: The complete Lorenz system differential equation described in S2 is as follows: In formula (1): [x, y, z] T is the three-dimensional state variable of the Lorenz system; is the differential value corresponding to the three-dimensional state variable; Pr>0 is the Prandtl number, which is related to the flow field; Ra>0 is the Rayleigh number, which is related to the flow field; β>0 is a spatial parameter that reflects the aspect ratio of the volume in the convection problem.

3. The method for exciting the chaotic state of a gliding aircraft based on a Lorenz system according to claim 2, characterized in that: The calculation process of the three-axis velocity increment of the desired reentry system of the aircraft described in S3 is as follows: S301: Taking the three-axis position deviation of the spacecraft reentry system as input and the desired three-axis velocity increment of the spacecraft reentry system as output, the complete Lorenz system differential equation is rewritten as follows: In formula (2): [Δx, Δy, Δz] T is the three-axis position deviation of the spacecraft reentry system; ΔV=[ΔV x , ΔV y , ΔV z ] T is the three-axis velocity increment of the desired spacecraft reentry system; S302: Introducing range control parameters [R x , R y , R z ] T The new differential equation obtained by combining the chaos construction speed control parameter k is as follows: S303: According to the maneuvering range requirement, maneuvering time requirement and overload constraint, the parameters of the Lorenz system are set as P=10, Ra=28, β=8 / 3, Rx=5000, Rz=10000, Ry=5000, k=100, and the final input-output relationship is as follows: S304: The three-axis position deviation [Δx, Δy, Δz] of the aircraft reentry system calculated at each simulation step T , that is, to obtain the desired three-axis velocity increment of the reentry system of the spacecraft ΔV = [ΔV x , ΔV y , ΔV z ] T .

4. The method for exciting the chaotic state of a gliding aircraft based on a Lorenz system according to claim 3, characterized in that: The calculation process of the three-axis expected velocity of the chaotic ballistic reentry system described in S4 is as follows: S401: Define the current nominal trajectory state variables as [θ, φ, h, V, γ, ψ] T 0, where θ is longitude, φ is latitude, h is altitude, V is velocity, γ is track angle, and ψ is heading angle. The velocity vector in the North Celestial East coordinate system is: S402: Convert the velocity vector in the North Celestial East coordinate system to the aircraft reentry system and express it as follows: In formula (6): is the attitude transformation matrix of the spacecraft reentry system relative to the Earth-centered fixed coordinate system; is the attitude transformation matrix of the geocentric fixed coordinate system relative to the North Celestial East coordinate system; is the velocity component of the nominal trajectory in the reentry system at the current moment; S403: Similarly to S401-S402, define the state variables of the current chaotic flight as [θ, φ, h, V, γψ] T , get the velocity component of chaotic flight in the reentry system at the current moment S404: The velocity component of the nominal trajectory at the current moment in the reentry system and the desired three-axis velocity increment of the reentry system of the spacecraft ΔV = [ΔV x , ΔV y , ΔV z ] T Obtain the expected velocity component of chaotic flight: S405: Expected velocity component from chaotic flight and the velocity component of the chaotic flight in the reentry system at the current moment The additional lateral overload Δn is obtained by PD tracking method z and additional longitudinal overload Δn y : In formula (8): K p is the proportionality coefficient; K d is the differential coefficient; is the y-axis velocity component of the reentry system during the actual flight at the previous simulation step; is the expected y-axis velocity component of the reentry system at the previous simulation step; is the z-axis velocity component of the reentry system during the actual flight at the previous simulation step; is the expected z-axis velocity component of the reentry system at the previous simulation step.

5. The method for exciting the chaotic state of a gliding aircraft based on a Lorenz system according to claim 4, characterized in that: The S5 comprises the following steps: S501: Calculation of lateral overload n from nominal ballistic state variables z0 and longitudinal overload n y0 , combined with additional lateral overload Δn z and additional longitudinal overload Δn y Calculate the required lateral overload n for chaotic maneuvering flight z and required longitudinal overload n y : S502: Calculate the required lateral overload n of the glider z and required longitudinal overload n y : In formula (10): is the aircraft dynamic pressure, where: ρ is the atmospheric density, V is the aircraft speed; S is the reference aerodynamic area of ​​the aircraft; C L is the lift coefficient, which is a function of the angle of attack α; C D is the drag coefficient, which is a function of the angle of attack α; mg is the gravity on the aircraft; S503: According to formula (9), the required lateral overload n is known at this time z and required longitudinal overload n y , inversely solve the corresponding angle of attack α and roll angle σ through Newton iteration method; S504: Combine the state variables of the current chaotic flight [θ, φ, h, V, γ, ψ] T , angle of attack α and roll angle σ, and calculate the differential values ​​of each state variable by the reentry glide dynamics equation: In formula (11): r = R0 + h is the distance from the center of the earth to the aircraft, where: R0 is the radius of the earth; g is the acceleration due to gravity of the Earth; m is the weight of the aircraft; S505: Based on the differential value of each state variable The next state variable is obtained through the fourth-order Runge-Kutta integration method, and then the aircraft state update is completed.

Images