CN115686059A - Hypersonic aircraft flight-forbidden region avoidance guidance method based on pseudo-spectrum method - Google Patents
Hypersonic aircraft flight-forbidden region avoidance guidance method based on pseudo-spectrum method Download PDFInfo
- Publication number
- CN115686059A CN115686059A CN202211192338.2A CN202211192338A CN115686059A CN 115686059 A CN115686059 A CN 115686059A CN 202211192338 A CN202211192338 A CN 202211192338A CN 115686059 A CN115686059 A CN 115686059A
- Authority
- CN
- China
- Prior art keywords
- aircraft
- pseudo
- hypersonic aircraft
- reentry
- hypersonic
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Pending
Links
Images
Classifications
-
- Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
- Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
- Y02T—CLIMATE CHANGE MITIGATION TECHNOLOGIES RELATED TO TRANSPORTATION
- Y02T90/00—Enabling technologies or technologies with a potential or indirect contribution to GHG emissions mitigation
Landscapes
- Feedback Control In General (AREA)
Abstract
The invention discloses a hypersonic aircraft flight control forbidden region avoidance guidance method based on a pseudo-spectral method, which comprises the following steps: establishing an unpowered reentry guidance kinetic equation of the hypersonic aircraft; establishing an environment atmospheric density and gravity acceleration model; analyzing a hypersonic aircraft pneumatic model; establishing process constraint models of a hypersonic aircraft, such as heat flow, dynamic pressure, overload, no-fly zone and the like; discretizing a reentry process state variable and a control variable based on a Radau pseudo-spectrum method principle according to the model and the constraint; for the well-constructed optimization problem, the initial state variable and the estimated value of the control variable are solved by adopting a fourth-order Ronge-Kutta method; based on the discrete state and the control variable, the minimum heat flow and the track oscillation are taken as optimization targets, a nonlinear programming problem is established, and the discrete state is subjected to self-adaptive correction according to the precision requirement in the solving process. The method and the device can enhance the success rate of avoiding the multi-forbidden flight areas in the complex environment in the re-entry guidance process.
Description
Technical Field
The invention relates to the technical field of aerospace, in particular to a hypersonic aircraft flight control area avoidance guidance method based on a pseudo-spectrum method.
Background
The hypersonic aerocraft is an important guarantee of the national defense industrial defense system in China, is a national emphasis device and has an unimpacting status. Generally, a hypersonic aerocraft refers to an aerocraft with a flying speed of more than 5Ma, adopts a rotating body with a larger lift-drag ratio or a wave-riding body type, is thrown aloft by a space combat aerocraft or boosted to a preset height by a rocket booster, and then realizes long-distance non-ballistic reentry flight by a gliding warhead depending on aerodynamic lift force, so that the hypersonic aerocraft has long-distance quick arrival capability, and further realizes long-distance accurate strike. The flying span of the aircraft realizes the long-distance high-speed flight in the atmosphere and across the atmosphere, namely below 20km and 20-100 km high altitude, has a plurality of advantages of high flying speed, good concealment, strong penetration resistance and the like, can greatly expand battlefield space, and has an important deterrent effect on military affairs. The future war is embodied by the complexity, informatization and intellectualization of a battlefield environment, and the air strike can be won by depending on the speed and the height. The hypersonic aircraft with ultrahigh speed and ultrahigh speed has the following advantages:
(1) The range of battle is wide, and the battle can reach the world in a short time. The hypersonic flight vehicle can reach any point on the earth within 2h, so that various military targets beyond thousands or tens of thousands of kilometers can be hit quickly, and the battlefield space is greatly expanded.
(2) Good concealment and high detection difficulty. The hypersonic aircraft can rapidly pass through an airspace for fighting against a defense system of the other party by virtue of extremely high flying speed and ultra-strong maneuvering capacity, and the accumulation quantity of the echo is small, so that the detection difficulty of a ground air defense weapon system is greatly increased, and the hypersonic aircraft is convenient for information investigation and battlefield monitoring.
(3) The prevention capacity and the survival capacity are strong. On the basis of good concealment, hypersonic flight causes the opposite side to intercept the weapon without enough reaction time and the flight speed for intercepting, so the burst protection capability is extremely strong, and the unexpected attack effect can be achieved. In addition, the high flight height, the extremely fast flight speed and the strong lateral maneuvering capability of the aircraft make any air defense weapon unable to hit the aircraft at present. If the stealth technology is applied to the hypersonic aircraft, the stealth performance of the hypersonic aircraft is better, and the defense burst and the survival ability are also stronger.
(4) The range is relatively far, and the combat efficiency is high. At present, hypersonic missiles which are researched abroad obtain huge kinetic energy at ranges of hundreds of kilometers or even thousands of kilometers by means of high-speed flight, and the target is killed and killed by direct collision, so that the combat efficiency is extremely high. The warhead power would be greater if a warhead of comparable mass to a subsonic vehicle were designed.
(5) The strategic deterrence effect is obvious. The significant advantage of hypersonic aircrafts in military affairs is to make the hypersonic aircrafts become a new generation of aerospace 'killer mace' with powerful strategic deterrence.
However, the flight environment of the hypersonic aircraft is extremely complex and is limited by various constraints, and particularly, in the process of executing tasks, the flight forbidden region is formed under the influence of factors such as geopolitical or military monitoring and the like, so that the complexity of the reentry trajectory optimization process is further increased. The reentry process of the hypersonic aircraft is limited by various process constraints, and due to the characteristic of unpowered gliding, the energy of the hypersonic aircraft needs to be carefully managed to ensure the subsequent attack destruction capability, while the occurrence of the no-fly zone can further limit the reentry trajectory planning, reduce the feasible flight area of the aircraft and provide further challenges for the reentry process of the aircraft.
Therefore, how to design a reentry state and control quantity sequence which can generate a reentry state and control quantity sequence under multiple constraint conditions, can satisfy process constraint and terminal constraint, can reach a specified area, can penetrate through a complex battlefield environment, and can satisfy smaller heat flow indexes and trajectory oscillation conditions becomes a problem to be solved urgently at present.
Disclosure of Invention
The invention aims to provide a hypersonic aircraft flight forbidding area avoidance guidance method based on a pseudo-spectrum method, so as to solve the problems.
The technical scheme adopted by the invention for solving the technical problems is as follows:
a hypersonic aircraft flight-forbidden zone avoidance guidance method based on a pseudo-spectrum method comprises the following steps:
s1, establishing an unpowered reentry guidance kinetic equation of the hypersonic aircraft;
s2, establishing an environment atmospheric density and gravity acceleration model;
s3, analyzing the hypersonic aircraft pneumatic model based on the S1 and the S2;
s4, establishing process constraint models of the hypersonic aircraft, such as heat flow, dynamic pressure, overload, no-fly zones and the like;
s5, discretizing reentry process state variables and control variables based on a Radau pseudo-spectrum method principle according to the model and the constraint;
s6, for the constructed optimization problem, integrating and solving the initial state variable and the control variable estimation value by adopting a four-order Ronge-Kutta method;
and S7, based on the discrete state and the control variable, establishing a nonlinear programming problem by taking the minimum heat flow and the track oscillation as optimization targets, and performing self-adaptive correction on the discrete state according to the precision requirement in the solving process.
Further, the equation established in step S1 is:
wherein r is the geocentric distance, λ is the longitude, φ is the latitude, S e Representing a range angle, wherein V is a flight speed, sigma is an inclination angle, theta is a ballistic inclination angle psi and a heading deflection angle, is an included angle between the horizontal projection of a velocity vector on the local and the true north direction, and is rotated clockwise to be positive; l, D represent aerodynamic lift and drag, m 0 And g represents the earth gravity acceleration of the mass of the aircraft and the height of the aircraft:
wherein R is 0 Representing the radius of the earth, S ref M, ρ represents the aircraft reference area, mass and atmospheric density, C L ,C D Representing the lift coefficient and the drag coefficient.
Further, the equation established in step S2 is:
ρ＝ρ 0 e -h/7200 (10)
where ρ is 0 ＝1.225kg/m 3 Is standard atmospheric density, g 0 ＝9.81m/s 2 Is standard earth gravity acceleration, h represents the flying height of the aircraft, R 0 And r represents the earth radius and the aircraft-to-geocentric distance.
Further, in step S3, the method for analyzing the hypersonic aircraft pneumatic model based on S1 and S2 is as follows:
assuming that the moment of the reentry gliding missile is in a balanced state and has no lateral force influence in the whole process, the aerodynamic force and the aerodynamic moment of the aircraft in the reentry process are mainly determined by the following parameters:
f and M respectively represent aerodynamic force and aerodynamic moment borne by the aircraft, and Re represents the measurement of the friction force by the Reynolds number; ma is the mach number, alpha, beta,
representing the airflow angle and the airflow angular velocity of the current gliding missile;
representing the rotation angular velocity of the single body of the gliding missile;
representing the current rudder deflection angle of the gliding missile; reference area S of wing ref Obtaining a function expression of pneumatic parameters after normalization processing on the current dynamic pressure q and the reference characteristic length l; shape, scale and power respectively represent the appearance characteristics, size and engine characteristics of the aircraft;
lift coefficient of current gliding missile:
current glide missile drag coefficient:
the above two physical quantities are determined by the following parameters:
C L |C D ＝f{M,Re,α,β,shape,poweron/off} (15)
after reasonable simplification, the final form is:
C L |C D ＝f{M,α,β} (16)
under the condition of hypersonic flight, the following aerodynamic coefficient fitting results are finally shown by combining aerodynamic data:
C L ＝C L00 +C L10 α+C L01 M+C L20 α 2 +C L11 αM+C L02 M 2 +C L21 α 2 M+C L12 αM 2 +C L03 M 3 +C L22 α 2 M 2 +C L13 αM 3 +C L04 M 4
C D ＝C D00 +C D10 α+C D01 M+C D20 α 2 +C D11 αM+C D02 M 2 +C D21 α 2 M+C D12 αM 2 +C D03 M 3 +C D22 α 2 M 2 +C D13 αM 3 +C D04 M 4 (17)
in the formula C Lij Representing the lift polynomial fitting coefficient, C Dij Representing the resistance polynomial fit coefficients.
Further, the method for establishing the process constraint models of the hypersonic aircraft such as heat flow, dynamic pressure, overload, no-fly zone and the like in the step S4 comprises the following steps:
the reentry process of the hypersonic aerocraft is mainly limited by heat flow, dynamic pressure and overload
In the formula (I), the compound is shown in the specification,
is the heat flow density at the stagnation point;
maximum allowable heat flux density; rho, rho 0 Respectively representing the corresponding atmospheric density of actual altitude and the standard atmospheric density at sea level, C1 and R d Respectively, coefficients relating to the aircraft;
respectively representing the current flight speed and the normalized speed; q. q.s max N is the total overload constraint for the maximum allowable dynamic pressure; n is max Maximum total overload allowed;
the no-fly zone is regarded as an infinite-height cylinder, and the track of the aircraft can only be avoided from the left side and the right side without considering the situation of avoiding from the upper side or the lower side; let λ be m ,φ m ,λ z ,φ z ,R z ,R e Respectively representing the current longitude and latitude of the aircraft, the longitude and latitude of the center of the no-fly zone, the radius of the no-fly zone and the radius of the earth, and then satisfying the path constraint shown in the formula:
further, the method of step S5 is:
dispersing a control variable and a state variable on a series of Legendre-Gauss-Rafau points through a Radau pseudo-spectrum method, approximating the control variable and the state variable through an interpolation polynomial, approximating a kinetic differential equation through derivation of the state variable, converting a continuous optimal control problem into a series of NLP problems constrained by algebra, and solving the problems by adopting a corresponding numerical method; the method mainly comprises the following steps:
(1) Time domain transformation;
(2) Approximating the state variable and the control variable polynomial based on time domain transformation;
(3) Based on time domain transformation, carrying out constraint transformation on a differential equation;
(4) Based on time domain transformation, constraint transformation is carried out on the terminal;
(5) Based on the time domain transformation, a performance index function approximation is obtained.
Further, the method of step S6 is:
the Ronge-Kutta integral equation is as follows:
in the formula y n ,y n+1 Respectively representing the current time and the next time of the dependent variable of the functional relationship f, x n ,h d Respectively representing the value of the independent variable at the current moment and the updating step length;
and obtaining a state variable and a control variable estimated value capable of roughly completing a reentry task by presetting a roll angle control quantity, a preset attack angle profile and course angle deviation corridor guidance so as to provide initial value guess for subsequent nonlinear programming problem solving.
Further, the method of step S7 is:
obtaining an NLP problem after dispersing state variables, control variables and constraint conditions, setting an objective function as minimum heat flow and track oscillation, starting to solve an optimization problem, and considering self-adaptive correction of a pseudo-spectrum orthogonal node in a solving process, wherein the method mainly comprises the following steps:
(1) Calculating approximation precision;
(2) Designing a subdivision grid or increasing a judgment standard of distribution points according to the approximation precision calculation;
(3) Determining the number of the distribution points according to the judgment result of the judgment standard;
(4) And according to the judgment result of the judgment standard, determining the number and the position of the divided time regions.
The invention has the beneficial effects that:
the invention discloses a hypersonic aircraft reentry no-fly region evasion guidance method based on an adaptive pseudo-spectrum method, aiming at the high uncertainty, high nonlinearity, high complexity and high constraint environment of the reentry process of a hypersonic aircraft and the reentry minimum heat flow requirement and no-fly region evasion requirement, the original complex track planning problem is converted into a solved nonlinear planning problem on the basis of a discrete state variable and a control variable of the adaptive pseudo-spectrum method, so that the multi-constraint requirement and no-fly region evasion requirement of reentry track planning are realized, and a smaller pneumatic heating effect can be finally obtained.
Drawings
FIG. 1 is a flow chart of the method of the present invention.
Fig. 2 is a flow chart of trajectory planning.
Fig. 3 is a coordinate system transformation diagram.
FIG. 4 is a graph of lift coefficient fit effects.
Fig. 5 is a graph of the effect of the drag coefficient fit.
FIG. 6 is a graph showing the effect of lift-drag ratio fitting.
Fig. 7 is a graph of height versus speed.
Fig. 8 is a graph showing the change of the ballistic inclination angle and the ballistic declination angle.
FIG. 9 is a graph of the angle of attack versus roll angle.
FIG. 10 is a schematic diagram of a constraint satisfaction scenario.
Fig. 11 is a schematic diagram of a no-fly zone avoidance 1.
Fig. 12 is a schematic diagram of a no-fly zone avoidance 2.
Detailed Description
In order to make the objects, technical solutions and advantages of the embodiments of the present invention clearer, the technical solutions in the embodiments of the present invention will be clearly and completely described below with reference to the drawings in the embodiments of the present invention, and it is obvious that the described embodiments are some, but not all, embodiments of the present invention. All other embodiments, which can be derived by a person skilled in the art from the embodiments given herein without making any creative effort, shall fall within the protection scope of the present invention.
Referring to fig. 1, the invention discloses a hypersonic aircraft flight control area avoidance guidance method based on a pseudo-spectrum method, which comprises the following steps:
s1, establishing an unpowered reentry guidance kinetic equation of the hypersonic aircraft;
s2, establishing an environment atmospheric density and gravity acceleration model;
s3, analyzing the hypersonic aircraft pneumatic model based on the S1 and the S2;
s4, establishing process constraint models of a hypersonic aircraft, such as heat flow, dynamic pressure, overload, no-fly zone and the like;
s5, discretizing reentry process state variables and control variables based on a Radau pseudo-spectrum method principle according to the model and the constraint;
s6, for the constructed optimization problem, integrating and solving the initial state variable and the control variable estimation value by adopting a four-order Ronge-Kutta method;
and S7, based on the discrete state and the control variable, establishing a nonlinear programming problem by taking the minimum heat flow and the track oscillation as optimization targets, and performing self-adaptive correction on the discrete state according to the precision requirement in the solving process.
Fig. 2 shows a whole reentry trajectory planning process, which includes the steps of initial value generation, pseudo-spectrum discrete variable generation, solution, and the like.
The method of step S1 is as follows:
according to the relevant angles defined by the coordinate system conversion relation shown in FIG. 3, the equation of motion of the reentry process of the hypersonic flight vehicle is established as follows:
wherein r is geocentric distance, λ is longitude, φ is latitude, S e Representing a range angle, wherein V is a flight speed, sigma is an inclination angle, theta is a ballistic inclination angle psi and a heading deflection angle, is an included angle between the horizontal projection of a velocity vector on the local and the true north direction, and is rotated clockwise to be positive; l, D represent aerodynamic lift and drag, respectively, m 0 G represents the earth gravity acceleration of the aircraft mass and the height of the aircraft:
wherein R is 0 Representing the radius of the earth, S ref M, p denote the aircraft reference planeVolume, mass and atmospheric density, C L ,C D Representing the lift coefficient and the drag coefficient.
The equation established in step S2 is:
ρ＝ρ 0 e -h/7200 (10)
where ρ is 0 ＝1.225kg/m 3 Is standard atmospheric density, g 0 ＝9.81m/s 2 Is standard earth gravity acceleration, h represents the flying height of the aircraft, R 0 And r represents the earth radius and the aircraft-to-geocentric distance.
The method for analyzing the hypersonic aircraft pneumatic model based on the S1 and the S2 in the step S3 comprises the following steps:
assuming that the moment of the reentry gliding missile is in a balanced state and has no lateral force influence in the whole process, the aerodynamic force and the aerodynamic moment of the aircraft in the reentry process are mainly determined by the following parameters:
f and M respectively represent aerodynamic force and aerodynamic moment borne by the aircraft, and Re represents the measurement of the friction force by the Reynolds number; ma is the mach number, alpha, beta,
representing the airflow angle and the airflow angular velocity of the current gliding missile;
representing the rotation angular velocity of the single body of the gliding missile;
representing the current rudder deflection angle of the gliding missile; reference area S of wing ref The function expression of the pneumatic parameters is obtained after normalization processing of the current dynamic pressure q and the reference characteristic length lFormula (I); shape, scale and power respectively represent the appearance characteristics, size and engine characteristics of the aircraft;
lift coefficient of current gliding missile:
current gliding missile resistance coefficient:
the above two physical quantities are determined by the following parameters:
C L |C D ＝f{M,Re,α,β,shape,poweron/off} (15)
after reasonable simplification, the final form is:
C L |C D ＝f{M,α,β} (16)
under the condition of hypersonic flight, the following aerodynamic coefficient fitting results are finally shown by combining aerodynamic data:
C L ＝C L00 +C L10 α+C L01 M+C L20 α 2 +C L11 αM+C L02 M 2 +C L21 α 2 M+C L12 αM 2 +C L03 M 3 +C L22 α 2 M 2 +C L13 αM 3 +C L04 M 4
C D ＝C D00 +C D10 α+C D01 M+C D20 α 2 +C D11 αM+C D02 M 2 +C D21 α 2 M+C D12 αM 2 +C D03 M 3 +C D22 α 2 M 2 +C D13 αM 3 +C D04 M 4 (17)
in the formula C Lij Representing the lift polynomial fitting coefficient, C Dij Representing the resistance polynomial fit coefficients.
As shown in fig. 4, 5 and 6, fitting coefficients related to the aircraft are substituted, and fitting results of lift force, resistance and lift-drag ratio corresponding to the aircraft are obtained through simulation.
The method for establishing the process constraint models of the hypersonic aircraft in the step S4, such as the heat flow, the dynamic pressure, the overload, the no-fly zone and the like, comprises the following steps:
the reentry process of the hypersonic aerocraft is mainly limited by heat flow, dynamic pressure and overload
In the formula (I), the compound is shown in the specification,
is the heat flow density at the stagnation point;
maximum allowable heat flux density; ρ, ρ 0 Respectively representing the corresponding atmospheric density of actual altitude and the standard atmospheric density at sea level, C1 and R d Respectively, coefficients relating to the aircraft;
respectively representing the current flight speed and the normalized speed; q. q of max N is the total overload constraint for the maximum allowable dynamic pressure; n is max Maximum total overload allowed;
the no-fly zone is an area formed by the influence of radar detection, electromagnetic interference, interception, terrain, geopolitical factors and other factors on the reentry aircraft, and the reentry track should avoid the no-fly zone as much as possible. For the convenience of design analysis, the no-fly zone is regarded as an infinite height cylinderThe aircraft track can only be avoided from the left side and the right side, and the situation of avoiding from the upper side or the lower side is not considered; suppose λ m ,φ m ,λ z ,φ z ,R z ,R e Respectively representing the current longitude and latitude of the aircraft, the longitude and latitude of the center of the no-fly zone, the radius of the no-fly zone and the radius of the earth, and then satisfying the path constraint shown in the formula:
the method of step S5 is:
dispersing the control variable and the state variable on a series of Legendre-Gauss-Radau points through a Radau pseudo-spectrum method, approximating the control variable and the state variable through an interpolation polynomial, approximating a dynamic differential equation through derivation of the state variable, converting a continuous optimal control problem into a series of NLP problems subjected to algebraic constraints, and solving the problems by adopting a corresponding numerical method; the method mainly comprises the following steps:
(1) Time domain transformation;
(2) Approximating the state variable and the control variable polynomial based on time domain transformation;
(3) Based on time domain transformation, carrying out constraint transformation on a differential equation;
(4) Based on time domain transformation, constraint transformation is carried out on the terminal;
(5) Based on the time domain transformation, a performance index function approximation is obtained.
The method of step S5 specifically includes:
for the traditional optimal control problem, a state variable x (·) epsilon R is defined n Control variable u (·) epsilon R m The static parameter p is equal to R q Initial test time t 0 Terminal time t f (free or fixed), minimization of Bolza-type performance metrics in general:
satisfying the differential equation constraint:
and (3) path constraint:
C(x(t),u(t),t；p)≤0 (25)
and (3) boundary constraint:
wherein p = [ p ] 1 ,p 2 ,...,p q ] T ，x＝[x 1 ,x 2 ,...,x n ] T ，u＝[u 1 ,u 2 ,...,u m ] T 。
The idea of solving the optimal control problem by RPM is to disperse a control variable and a state variable on a series of Legendre-Gauss-Radau (LGR) points, approximate the control variable and the state variable by an interpolation polynomial, approximate a kinetic differential equation by derivation of the state variable, convert a continuous optimal control problem into a series of NLP problems constrained by algebra, and solve the problem by adopting a corresponding numerical method, and the steps are as follows:
(1) Time domain transformation
The time interval of the track optimization problem needs to be changed from interval [ t ] by adopting a Radau pseudo-spectrum method 0 ,t f ]Conversion to [ -1,1]The transformation is as follows:
where τ is the current time after transformation, t 0 ,t f Representing the current time, the initial time and the terminal time;
(2) Polynomial approximation of state variables to control variables
The matching point of the Radau pseudo-spectrum method is an LGR point, and the LGR point is in the interval (-1, 1)]And (6) taking the value. LGR points of order k are polynomials P k (ζ)+P k-1 (ζ) root of which P k (ζ) is a Legendre polynomial of order k, defined as follows:
in order to make the matching points cover to [ -1,1]Selecting nodes of Radau pseudo-spectrum method as nodes and initial time nodes zeta 0 And (4) = -1. That is, if the number of discrete nodes in the optimization problem is N, LGR nodes with K = N-1 are adopted, and the state variables are discretized as follows by combining the initial nodes:
τ i the discrete points are N, the number of discrete points is X, and X respectively represent the actual state quantity and the discrete approximate state quantity.
Among them are Lagrange's interpolation basis functions
And make the approximate state at the node equal to the actual state, there is X i ＝X(τ i )＝x(τ i )。
The derivative to the control variable is not used in the Radau pseudospectral method, so u (τ) is satisfied k )＝U k N-1), again using LGR nodes, using Lagrange interpolation basis functions to approximate the controlled variables to, U representing the actual controlled variable versus the approximate controlled variable:
(3) Differential equation constrained transformation
Differentiating the formula by
Wherein the differential matrix satisfies:
thereby converting the kinetic differential equation constraints into algebraic constraints:
n-1, where k = 1.
(4) Terminal constraints
The terminal state should satisfy the kinetic equation constraints as follows:
the terminal constraint condition can be discretized and a Gauss integral is adopted to discretize the approximation, so that:
in the formula:
τ f representing the terminal time after time domain conversion;
(5) Performance indicator function approximation
Approximating a Bolza type performance index function integral term through Gauss integration to obtain an approximate performance index function in a Radau pseudo-spectrum method as follows:
in the equation, (+), g () represents a terminal index function and an integral index function.
Through the mathematical transformation, the original continuous optimal control problem is converted into a general nonlinear programming problem.
Solving for state X on discrete nodes i (i =0,1, N-1) and a control variable U k (k =1, 2.., N-1), and an initial terminal time t 0 ,t f And minimizing the performance index formula and meeting the boundary condition formula of the original optimal control problem.
The method of step S6 is:
the Ronge-Kutta integral equation is as follows:
in the formula y n ,y n+1 Respectively representing the current time and the next time of the dependent variable of the functional relationship f, x n ,h d Respectively representing the value of the independent variable at the current moment and the updating step length;
and obtaining a state variable and a control variable estimated value capable of roughly completing a reentry task by presetting a roll angle control quantity, a preset attack angle profile and course angle deviation corridor guidance so as to provide initial value guess for subsequent nonlinear programming problem solving.
Step S7 specifically includes:
obtaining an NLP problem after dispersing state variables, control variables and constraint conditions by adopting a Radau pseudo-spectrum method, setting an objective function as minimum heat flow and trajectory oscillation, and starting to solve an optimization problem. The longitudinal track of the aircraft during reentry is often a leapfrog forward 'Qian schsen' trajectory, and the aircraft has strong maneuvering and defense-breaking capabilities, which also causes the reentry track to vibrate. The Radau pseudo spectrum method described above can be divided into a global pseudo spectrum method (p-RPM) and a local pseudo spectrum method (h-RPM), the global RPM has a fast convergence rate for smooth trajectory optimization, but if the state is oscillating, it is difficult to ensure the accuracy of the optimization process; the local RPM effectively improves the precision by dividing the time domain and fitting each sub-time domain by a low-order polynomial, but the convergence speed is low. The hp self-adaptive pseudo-spectrum method developed in recent years can effectively combine the advantages of the hp self-adaptive pseudo-spectrum method and the hpadaptive pseudo-spectrum method, reasonably divides the whole time domain into smooth sub-time domains, and is suitable for the reentry problem of an aircraft.
The main steps of the hp self-adaptive pseudo-spectrum are as follows:
(1) And (3) approximation precision calculation:
calculating the first order differential residual error and the process constraint residual error of the state variable at each collocation point in each time subinterval:
in the formula, i and j respectively represent collocation and self-interval marks in an interval, f and C are respectively represented by a dynamic equation and process constraints;
from the above, the maximum residual error in this interval is:
setting a maximum allowable deviation e if e j If the value is less than epsilon, the state and the controlled variable at the interval j meet the requirements of differential equation constraint and process constraint, stopping iterative calculation of the interval, and otherwise, further improving the optimization precision by subdividing a time area or increasing the number of distribution points in the area; and obtaining a feasible solution when all the self intervals meet the precision requirement.
(2) Criteria for subdividing the grid or adding distribution points
The subdivision time (h) or the increased number (p) of distribution points is determined by the curvature of the sub-time domain interval, if the curvature is larger than a set value, the interval is more oscillatory, the time region needs to be further divided, otherwise, the number of distribution points is increased.
Defining the relative curvature of the jth sub-temporal region as:
wherein
Respectively representing the maximum curvature and the average curvature of the matching points in the interval, and the curvature of the matching points is calculated in the following way:
(3) Point placement number determination
If the interval is smooth, the number of the configuration points is increased, and the configuration points are determined by the following formula:
in the formula, ceil represents rounding upwards, and A is an adjustable parameter larger than 0; the above equation shows that the larger the sub-interval residual error is, the more the number of dispensing points to be added is, and the less the number of dispensing points is increased in the contrary manner.
(4) Time zone number and location determination
The number of divided time regions is determined by the following formula:
b is an adjustable parameter; the above equation shows that the larger the residual error is, the more the trajectory oscillates, and more sub-time regions need to be divided. The positions of the time division regions are determined according to the average curvature, the larger the curvature is, the more the track vibrates, and the time period is narrower to ensure the optimization precision. Curvature density coefficient constant c defining jth time region j Satisfies the following conditions:
the ith subdivision point in the jth time zone should be arranged at the collocation point
Such that:
on the basis of the adjustment, the problem of trajectory planning after self-adaptive pseudo-spectrum method dispersion is solved, as shown in fig. 7, the change situation of the height and the speed in the reentry process is shown, and the height and speed constraint of the terminal can be met; as shown in fig. 8, the variation of the ballistic inclination angle and the ballistic declination angle during reentry is shown; as shown in fig. 9, the change of the attack angle and the roll angle during reentry is shown; as shown in fig. 10, process constraint satisfaction during reentry is demonstrated, and the graphical curves show that the process constraints can be satisfied; as shown in fig. 11 and 12, it is shown that the reentry procedure can avoid one or more no-fly regions and meet the reentry trajectory requirement.
Aiming at the high uncertainty, high nonlinearity, high complexity and high constraint environment of the reentry process of the hypersonic flight vehicle, the reentry minimum heat flow requirement and the flight control forbidding area avoiding requirement, the invention converts the original complex track planning problem into the nonlinear programming problem which can be solved on the basis of the discrete state variable and the control variable of the self-adaptive pseudo-spectral method, thereby realizing the multi-constraint requirement and the flight control forbidding area avoiding requirement of the reentry track planning and finally obtaining the smaller pneumatic heating effect.
Finally, it should be noted that: the above examples are only intended to illustrate the technical solution of the present invention, and not to limit it; although the present invention has been described in detail with reference to the foregoing embodiments, it should be understood by those of ordinary skill in the art that: the technical solutions described in the foregoing embodiments may still be modified, or some technical features may be equivalently replaced; and such modifications or substitutions do not depart from the spirit and scope of the corresponding technical solutions of the embodiments of the present invention.
Claims (8)
1. A hypersonic aircraft flight-forbidden zone avoidance guidance method based on a pseudo-spectrum method is characterized by comprising the following steps:
s1, establishing an unpowered reentry guidance kinetic equation of the hypersonic aircraft;
s2, establishing an environment atmospheric density and gravity acceleration model;
s3, analyzing the hypersonic aircraft pneumatic model based on the S1 and the S2;
s4, establishing process constraint models of the hypersonic aircraft, such as heat flow, dynamic pressure, overload, no-fly zones and the like;
s5, discretizing reentry process state variables and control variables based on a Radau pseudo-spectrum method principle according to the model and the constraint;
s6, for the constructed optimization problem, integrating and solving the initial state variable and the control variable estimation value by adopting a four-order Ronge-Kutta method;
and S7, based on the discrete state and the control variable, establishing a nonlinear programming problem by taking the minimum heat flow and the track oscillation as optimization targets, and performing self-adaptive correction on the discrete state according to the precision requirement in the solving process.
2. The hypersonic aircraft flight-forbidden region avoidance guidance method based on the pseudo-spectrum method as claimed in claim 1, wherein the equation established in step S1 is as follows:
wherein r is the geocentric distance, λ is the longitude, φ is the latitude, S e The angle represents a range angle, V is a flight speed, sigma is an inclination angle, theta is a trajectory inclination angle psi and a heading deflection angle, is an included angle between the local horizontal projection of a velocity vector and the true north direction, and is changed into positive clockwise rotation; l, D represent aerodynamic lift and drag, m 0 G represents the earth gravity acceleration of the aircraft mass and the height of the aircraft:
wherein R is 0 Representing the radius of the earth, S ref M, ρ represents the aircraft reference area, mass and atmospheric density, C L ,C D Representing the lift coefficient and the drag coefficient.
3. The hypersonic aircraft flight-forbidden region avoidance guidance method based on the pseudo-spectrum method as claimed in claim 2, wherein the equation established in step S2 is as follows:
ρ＝ρ 0 e -h/7200 (10)
wherein ρ 0 ＝1.225kg/m 3 Is standard atmospheric density, g 0 ＝9.81m/s 2 Is standard earth gravitational acceleration, h represents aircraft flight altitude, R 0 And r represents the earth radius and the aircraft-to-geocentric distance.
4. The hypersonic aircraft flight-forbidden region avoidance guidance method based on the pseudo-spectrum method as claimed in claim 3, wherein the method for analyzing the hypersonic aircraft pneumatic model based on S1 and S2 in the step S3 is as follows:
assuming that the moment of the reentry gliding missile is in a balanced state and has no lateral force influence in the whole process, the aerodynamic force and the aerodynamic moment of the aircraft in the reentry process are mainly determined by the following parameters:
supposing that the moment of the reentry gliding missile is in a balanced state and has no lateral force influence in the whole process, the aerodynamic force and the aerodynamic moment suffered by the aircraft in the reentry process are mainly determined by the following parameters:
f and M respectively represent aerodynamic force and aerodynamic moment borne by the aircraft, and Re represents the Reynolds number to measure the magnitude of the friction force; ma is the mach number, alpha, beta,
representing the current airflow angle and the airflow angular velocity of the current gliding missile;
characterizing the rotation angular speed of the monomer of the gliding missile;
representing the current rudder deflection angle of the gliding missile; reference area S of wing ref Obtaining a function expression of pneumatic parameters after normalization processing on the current dynamic pressure q and the reference characteristic length l; shape, scale and power respectively represent the appearance characteristics, size and engine characteristics of the aircraft;
lift coefficient of current gliding missile:
current gliding missile resistance coefficient:
the above two physical quantities are determined by the following parameters:
C L |C D ＝f{M,Re,α,β,shape,poweron/off} (15)
after reasonable simplification, the final form is:
C L |C D ＝f{M,α,β} (16)
under the condition of hypersonic flight, the following aerodynamic coefficient fitting results are finally shown by combining aerodynamic data:
C L ＝C L00 +C L10 α+C L01 M+C L20 α 2 +C L11 αM+C L02 M 2 +C L21 α 2 M+C L12 αM 2 +C L03 M 3 +C L22 α 2 M 2 +C L13 αM 3 +C L04 M 4
C D ＝C D00 +C D10 α+C D01 M+C D20 α 2 +C D11 αM+C D02 M 2 +C D21 α 2 M+C D12 αM 2 +C D03 M 3 +C D22 α 2 M 2 +C D13 αM 3 +C D04 M 4 (17)
in the formula C Lij Representing the lift polynomial fitting coefficient, C Dij Representing the resistance polynomial fit coefficients.
5. The hypersonic aircraft flight-forbidden region avoidance guidance method based on the pseudo-spectrum method as claimed in claim 4, wherein the method for establishing the process constraint models of the hypersonic aircraft in the step S4, such as heat flow, dynamic pressure, overload, flight-forbidden region and the like, comprises the following steps:
the reentry process of the hypersonic aerocraft is mainly limited by heat flow, dynamic pressure and overload
In the formula (I), the compound is shown in the specification,
is the heat flow density at the stagnation point;
maximum allowable heat flux density; rho, rho 0 Respectively representing the corresponding atmospheric density of actual altitude and the standard atmospheric density at sea level, C1 and R d Respectively, coefficients relating to the aircraft; v is used for carrying out the treatment of the waste gas,
respectively representing the current flight speed and the normalized speed; q. q.s max N is the total overload constraint for the maximum dynamic pressure allowed; n is max Maximum total overload allowed;
the no-fly zone is regarded as an infinite-height cylinder, and the track of the aircraft can only be avoided from the left side and the right side without considering the situation of avoiding from the upper side or the lower side; let λ be m ,φ m ,λ z ,φ z ,R z ,R e Respectively representing the current longitude and latitude of the aircraft, the longitude and latitude of the center of the no-fly zone, the radius of the no-fly zone and the radius of the earth, and then satisfying the path constraint shown in the formula:
6. the hypersonic aircraft flight-forbidden region avoidance guidance method based on the pseudo-spectral method as claimed in claim 5, wherein the method in step S5 is as follows:
dispersing a control variable and a state variable on a series of Legendre-Gauss-Rafau points through a Radau pseudo-spectrum method, approximating the control variable and the state variable through an interpolation polynomial, approximating a kinetic differential equation through derivation of the state variable, converting a continuous optimal control problem into a series of NLP problems constrained by algebra, and solving the problems by adopting a corresponding numerical method; the method mainly comprises the following steps:
(1) Time domain transformation;
(2) Approximating the state variable and the control variable polynomial based on time domain transformation;
(3) Based on time domain transformation, carrying out constraint transformation on a differential equation;
(4) Based on time domain transformation, carrying out constraint transformation on the terminal;
(5) Based on the time domain transformation, a performance index function approximation is obtained.
7. The hypersonic aircraft flight-forbidden zone avoidance guidance method based on the pseudo-spectrum method is characterized in that the method in the step S6 is as follows:
the Ronge-Kutta integral equation is as follows:
in the formula y n ,y n+1 Respectively representing the current time and the next time of the dependent variable of the functional relationship f, x n ,h d Respectively representing the value of the independent variable at the current moment and the updating step length;
and obtaining state variables and control variable estimated values capable of roughly completing reentry tasks by presetting a roll angle control quantity, a preset attack angle profile and course angle deviation corridor guidance so as to provide initial value guess for subsequent nonlinear programming problem solution.
8. The hypersonic aircraft flight-forbidden zone avoidance guidance method based on the pseudo-spectrum method as claimed in claim 7, wherein the method of step S7 is as follows:
obtaining an NLP problem after dispersing state variables, control variables and constraint conditions, setting an objective function as minimum heat flow and trajectory oscillation, starting to solve an optimization problem, and considering self-adaptive correction of a pseudo-spectrum orthogonal node in the solving process, wherein the method mainly comprises the following steps:
(1) Calculating approximation precision;
(2) Designing a judgment standard for subdividing grids or increasing distribution points according to the approximation precision calculation;
(3) Determining the number of the dispensing points according to the judgment result of the judgment standard;
(4) And according to the judgment result of the judgment standard, determining the number and the position of the divided time regions.
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202211192338.2A CN115686059A (en) | 2022-09-28 | 2022-09-28 | Hypersonic aircraft flight-forbidden region avoidance guidance method based on pseudo-spectrum method |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202211192338.2A CN115686059A (en) | 2022-09-28 | 2022-09-28 | Hypersonic aircraft flight-forbidden region avoidance guidance method based on pseudo-spectrum method |
Publications (1)
| Publication Number | Publication Date |
|---|---|
| CN115686059A true CN115686059A (en) | 2023-02-03 |
Family
ID=85065093
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN202211192338.2A Pending CN115686059A (en) | 2022-09-28 | 2022-09-28 | Hypersonic aircraft flight-forbidden region avoidance guidance method based on pseudo-spectrum method |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN115686059A (en) |
Cited By (1)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN117850248A (en) * | 2024-03-07 | 2024-04-09 | 中国人民解放军国防科技大学 | Spacecraft hidden maneuver trajectory planning method based on Gaussian pseudo-spectrum method |
-
2022
- 2022-09-28 CN CN202211192338.2A patent/CN115686059A/en active Pending
Cited By (2)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN117850248A (en) * | 2024-03-07 | 2024-04-09 | 中国人民解放军国防科技大学 | Spacecraft hidden maneuver trajectory planning method based on Gaussian pseudo-spectrum method |
| CN117850248B (en) * | 2024-03-07 | 2024-05-28 | 中国人民解放军国防科技大学 | Spacecraft hidden maneuver trajectory planning method based on Gaussian pseudo-spectrum method |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN110412874B (en) | Multi-missile cooperative guidance law design method for maneuvering target and time delay communication | |
| CN114003050B (en) | Active defense guidance method of three-body countermeasure strategy based on differential game | |
| Fonod et al. | Estimation enhancement by cooperatively imposing relative intercept angles | |
| CN108534614A (en) | A kind of real-time Predictor-corrector guidance method of three-dimensional omnidirectional | |
| CN114020019A (en) | Guidance method and device for aircraft | |
| CN110822994A (en) | Linear pseudo-spectrum spreading control guidance method with falling angle constraint | |
| CN115686059A (en) | Hypersonic aircraft flight-forbidden region avoidance guidance method based on pseudo-spectrum method | |
| Fonod et al. | Blinding guidance against missiles sharing bearings-only measurements | |
| CN111176315B (en) | Variant cruise missile control method based on L1 adaptive control | |
| Xu et al. | A trajectory design method for coupling aircraft radar cross-section characteristics | |
| Gaudet et al. | Terminal adaptive guidance for autonomous hypersonic strike weapons via reinforcement learning | |
| Yuru et al. | Impact time control guidance against maneuvering targets based on a nonlinear virtual relative model | |
| Farooq et al. | Trajectory optimization for air-to-surface missiles with imaging radars | |
| CN116795130A (en) | Obstacle avoidance guidance method for maneuvering target | |
| CN114610057B (en) | Design method for maneuver burst prevention strategy of high Mach aircraft | |
| Zhang et al. | Closed-form time-to-go estimation for proportional navigation guidance considering drag | |
| Fu et al. | Air combat assignment problem based on bayesian optimization algorithm | |
| Chen et al. | On-line trajectory generation of midcourse cooperative guidance for multiple interceptors | |
| CN112346474A (en) | Design method of differential game guidance law with limited time convergence | |
| Xu et al. | Optimal Design of Cooperative Penetration Trajectories for Multiaircraft | |
| Gour et al. | A novel IEG strategy for realistically modeled seeker-less interceptors | |
| Taur | Nonlinear guidance and navigation of a tactical missile with high heading error | |
| Nie et al. | Adaptive tracking algorithm based on 3D variable turn model | |
| Zhang et al. | Pseudospectral method for autonomous attack trajectory planning of a fixed-wing UCAV | |
| Chen et al. | A composite optimal predictive guidance law for hypersonic target interception |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| PB01 | Publication | ||
| PB01 | Publication | ||
| SE01 | Entry into force of request for substantive examination | ||
| SE01 | Entry into force of request for substantive examination |