CN111306989A - Hypersonic velocity reentry guidance method based on steady glide trajectory analytic solution - Google Patents

Hypersonic velocity reentry guidance method based on steady glide trajectory analytic solution Download PDF

Info

Publication number
CN111306989A
CN111306989A CN202010170138.1A CN202010170138A CN111306989A CN 111306989 A CN111306989 A CN 111306989A CN 202010170138 A CN202010170138 A CN 202010170138A CN 111306989 A CN111306989 A CN 111306989A
Authority
CN
China
Prior art keywords
angle
terminal
taem
attack
representing
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.)
Granted
Application number
CN202010170138.1A
Other languages
Chinese (zh)
Other versions
CN111306989B (en
Inventor
陈万春
张晚晴
余文斌
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Beihang University
Original Assignee
Beihang University
Priority date (The priority date 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 date listed.)
Filing date
Publication date
Application filed by Beihang University filed Critical Beihang University
Priority to CN202010170138.1A priority Critical patent/CN111306989B/en
Publication of CN111306989A publication Critical patent/CN111306989A/en
Application granted granted Critical
Publication of CN111306989B publication Critical patent/CN111306989B/en
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • F—MECHANICAL ENGINEERING; LIGHTING; HEATING; WEAPONS; BLASTING
    • F41—WEAPONS
    • F41G—WEAPON SIGHTS; AIMING
    • F41G3/00—Aiming or laying means
    • F41G3/22—Aiming or laying means for vehicle-borne armament, e.g. on aircraft

Abstract

The invention relates to a hypersonic reentry guidance method based on a steady gliding trajectory analytic solution, which comprises the following steps: 1. establishing an aircraft dynamics model; 2. establishing a constraint model; 3. an initial descent stage guidance method; 4. a steady glide phase guidance method; 5. a height adjustment section guidance method; has the advantages that: the longitudinal and transverse courses are predicted and corrected based on the three-dimensional trajectory analytic solution, so that the calculation amount is greatly reduced, the calculation speed is improved, and the method is convenient to apply to an airborne system in real time. The transverse movement is determined by twice tilting reversal points, the position of the first reversal point is corrected to control the terminal traverse error, and the position of the second reversal point is corrected to meet the terminal speed constraint. Compared with the traditional method for controlling the transverse motion by using the course error threshold, the method has the advantages of reducing the number of times of tilting and reversing, improving the energy utilization efficiency and reducing the requirements on a control system. The height adjusting section utilizes the idea of proportional guidance, so that the terminal roll angle is converged to 0, and the target hitting effect is enhanced.

Description

Hypersonic velocity reentry guidance method based on steady glide trajectory analytic solution
Technical Field
The invention provides a hypersonic velocity reentry guidance method based on steady glide trajectory analytic solution, belonging to the technical field of space technology and weapons
Background
The design of the reentry guidance law is a core technology of whether the hypersonic reentry gliding aircraft can finish the accurate hitting task. The ultimate goal of reentry guidance is to achieve successful guidance of the hypersonic reentry vehicle to a specified target while satisfying path constraints and terminal constraints.
The reentry aircraft is subjected to strict path constraint in the reentry process, however, the traditional algorithm is always difficult to process the path constraint, and difficulties and challenges are brought to the online generation of the trajectory and the online calculation of the guidance law. Meanwhile, the extremely fast glide flight speed puts forward higher requirements on the calculation time, so that the time constraint generated by the online guidance law instruction is very strong, and for the task of target relocation, the algorithm is required to have stronger adaptivity and rapidity so as to adapt to the requirements of the large-maneuvering flight task.
Disclosure of Invention
Aiming at the problems of reentry guidance of the hypersonic aircraft, a high-precision hypersonic reentry guidance method based on a three-dimensional trajectory analytic solution is provided. Firstly, an auxiliary geocentric rotation coordinate system is constructed, a new kinetic model of the longitudinal course, the transverse course and the course angle with respect to the energy is established on the basis of the auxiliary geocentric rotation coordinate system, the model is favorable for decoupling and linearizing the highly nonlinear and strongly coupled kinetic model to obtain a three-dimensional ballistic analytic solution, and a guidance method is designed on the basis of the analytic solution. The guidance method utilizes a longitudinal analytic solution to plan a longitudinal reference profile in real time, utilizes a transverse analytic solution to correct a first inversion point to eliminate a transverse error, utilizes trajectory simulation for several times to correct a second inversion point to ensure terminal speed constraint and corrects the minimum value of an attack angle to ensure terminal height constraint.
The invention discloses a hypersonic reentry guidance method based on a steady glide trajectory analytic solution, which is characterized by comprising the following steps: the method comprises the following steps:
the method comprises the following steps: establishing an aircraft dynamics model;
the hypersonic flight vehicle generally refers to an aircraft with flight Mach number larger than 5 and capable of realizing hypersonic flight in an atmospheric layer and a trans-atmospheric layer. Under the background of the rotating earth, a standard atmosphere model and an inverse square gravitational field model are adopted to establish a six-degree-of-freedom kinetic equation of the hypersonic aerocraft as follows:
Figure BDA0002408897380000021
Figure BDA0002408897380000022
Figure BDA0002408897380000023
Figure BDA0002408897380000024
Figure BDA0002408897380000025
Figure BDA0002408897380000026
Figure BDA0002408897380000027
where d/dt represents the derivative over time, t is time, h is altitude, R is range, λ is longitude, φ is latitude, V is the velocity of the aircraft relative to the earth, γ is the ballistic inclination, ψ is the aircraft heading angle, based on local north, σ is the roll angle, R is the altitude, and R is the altitude, whereeIs the radius of the earth, and takes the value of 6378.137km, omegaeIs the angular velocity of rotation of the earth, g represents the acceleration of gravity; l ═ ρ V2SCLThe lift acceleration is given by/2 m, and D is equal to rho V2SCD(2 m) is the acceleration of resistance, CL、CDRespectively lift coefficient and drag coefficient.
Step two: establishing a constraint model;
the reentry aircraft needs to meet the process constraints of heat flux density, dynamic pressure, overload and the like in the flight process, and the process constraints are expressed as follows:
Figure BDA0002408897380000028
Figure BDA0002408897380000031
Figure BDA0002408897380000032
parameter p0Denotes nominal atmospheric density, e denotes natural constant, h denotes altitude, hSIndicating nominal height, n overload, and m mass.
Wherein, KQIs the heat flux density coefficient;
Figure BDA0002408897380000033
representing a maximum allowable heat flow density value; n represents an allowable overload value; n ismaxRepresents a maximum allowable overload value; q represents an allowable dynamic pressure value; q. q.smaxIndicating the maximum allowable dynamic pressure value.
Furthermore, due to the capability constraints of the flight control system, the rate of change of the control quantities (angle of attack α and roll angle σ) also satisfies the constraints:
Figure BDA0002408897380000034
in addition to process constraints, the aircraft should also satisfy terminal constraints, with: vTAEM＝2000m/s,HTAEM＝25km,Rtm＝100km,σTAEM≈0°,γTAEM≈0°andΔψTAEM0 deg. where the subscript TAEM denotes the terminal amount. VTAEM、HTAEM、σTAEM、γTAEM、ΔψTAEMRespectively representing the terminal speed, the terminal height, the terminal roll angle, the terminal trajectory inclination angle and the terminal course angle error. RtmIndicating the shot-eye distance.
Step three: an initial descent stage guidance method;
according to the characteristics of reentry trajectory, the reentry trajectory is divided into three parts: an initial descent section, a steady glide section and a height adjustment section. The initial descent segment is designed to allow a more gradual transition of the trajectory to a smooth glide segment. The whole flight trajectory is mostly in a stable gliding section, which not only determines the range, speed and height of the flight terminal, but also determines whether the aircraft can meet various process constraints. The height adjustment segment is the last phase of the trajectory, as to whether the trajectory can strictly meet the terminal constraints. The flying speed and flying height of the aircraft are low in this section, and severe requirements are provided for a guidance method for designing a height adjusting section.
In the initial descent, the variation of the controlled variable has little effect on the flight, so the fixed values of angle of attack and roll angle are designed here to avoid altitude variationsmaxFly 20deg and zero roll angle. Rate of change of ballistic inclination
Figure BDA0002408897380000035
By time, it is meant that lift can meet the smooth glide requirements so that the angle of attack can be smoothly transitioned from the maximum angle of attack to the reference angle of attack profile αbsl。
At this stage, the command angle of attack and roll are designed as follows:
Figure BDA0002408897380000041
σcmd＝0 (12)
wherein the content of the first and second substances,
Δγ＝γ-γSG(13)
Figure BDA0002408897380000042
wherein, γSGRepresenting the angle of inclination, k, of the smooth gliding trajectoryγThe feedback factor is expressed and here takes the value 5. when Δ γ is 0, the aircraft enters the smooth glide phase. parameter αcmdIs an instruction angle of attack; sigmacmdTilting to commandAn angle; delta gamma is the deviation of the current trajectory inclination angle and the steady glide trajectory inclination angle; delta gamma1As the rate of change of ballistic inclination
Figure BDA0002408897380000043
Deviation value of trajectory inclination angle at 0 from the trajectory inclination angle of smooth glide.
Step four: a steady glide phase guidance method;
ballistic planning in the steady glide phase is the basis of the proposed guidance algorithm, which requires approaching the aircraft to the target while guaranteeing process constraints. The algorithm determines the longitudinal range and lateral range profiles based on a new three-dimensional analytical solution. Here, a new ballistic analytical solution is first briefly introduced.
In order to meet the characteristic that the longitudinal stroke of the reentry trajectory is far greater than the transverse stroke, two coordinate systems are firstly established, and the meanings of the longitudinal stroke and the transverse stroke can be intuitively represented.
Defining an auxiliary earth-centered rotational reference system (AGR coordinate system) rotating with the earth: the origin is at the geocentric E,
Figure BDA0002408897380000044
the shaft points to the initial position of the aircraft,
Figure BDA0002408897380000045
the axis is in the plane of a great circle passing through the aircraft and the target point and perpendicular to the aircraft
Figure BDA0002408897380000046
The shaft is provided with a plurality of axial holes,
Figure BDA0002408897380000047
the axes may be determined according to the right hand rule.
Meanwhile, a local generalized north east down coordinate system (GNED coordinate system) is defined: the vertical projection point of the origin o from the mass center M of the aircraft to the ground;
Figure BDA0002408897380000048
the axis is vertically downwards directed to the center of the earth,
Figure BDA0002408897380000049
axis-oriented AGR coordinate system
Figure BDA00024088973800000410
The direction of the "north" direction is,
Figure BDA00024088973800000411
the axes are determined by the right hand rule.
Define the longitudinal extent as
Figure BDA00024088973800000412
Define the course as
Figure BDA00024088973800000413
Wherein
Figure BDA00024088973800000414
And
Figure BDA00024088973800000415
respectively representing a generalized longitude and a generalized latitude in an AGR coordinate system.
Defining the energy as
Figure BDA0002408897380000051
Defining the longitudinal lift-to-drag ratio Lcos sigma/D ═ L1D, transverse lift-drag ratio Lsin sigma/D ═ L2and/D. The generalized longitude and latitude coordinate system in the AGR coordinate axis is beneficial to linearizing and decoupling a complex dynamic model, and the following three-dimensional ballistic analysis solution is obtained.
Figure BDA0002408897380000052
Figure BDA0002408897380000053
Figure BDA0002408897380000054
Wherein the content of the first and second substances,
Figure BDA0002408897380000055
R*＝Re+H；
Figure BDA0002408897380000056
hij,i＝1,2；j＝0,1,α1and pijI is 1, 2; j is 0,1,2,3,4 are constant coefficients. E, E0Respectively representing the current energy and the initial energy value; x is the number ofERepresenting an energy integral variable; l is1Representing lift values in a longitudinal plane; x is the number ofC0Representing an initial traverse value;
Figure BDA0002408897380000057
representing an initial course angle error under a GNED coordinate system;
Figure BDA0002408897380000058
indicating the heading angle error in the GNED coordinate system.
Since the analytic solution has high precision, the analytic solution can be applied to guidance method derivation to reduce calculation amount and improve universality.
1. Design of reference angle of attack
Because the analytic solution is designed with energy as the argument, to make it more convenient to utilize the analytic solution, all reference profiles in the guidance method also use energy as the argument, the reference angle of attack profile αbslThe design is as follows:
Figure BDA0002408897380000061
wherein E isα＝-5.55×107J/kg is located near the intersection point of the smooth gliding section and the final altitude adjusting section, in order to exert the maximum capability of the aircraft, the attack angle is designed to be α1The quadratic function for which the angle of attack is designed to be energy is designed so that the angle of attack can be from α1Gently transition to α2. In this systemIn the method, α is preset for the first time2The ballistic simulation is later applied to correct it to strictly meet the terminal height constraint, 6 deg. When the design of the attack angle section is finished, the corresponding reference lift-drag ratio section (L/D)bslAs determined accordingly.
2. Baseline lift-to-drag ratio profile design
In order to meet the requirement of range, the longitudinal lift-drag ratio profile is designed as follows
Figure BDA0002408897380000062
Wherein, let the terminal longitudinal analysis solve xDfEqual to the remaining range, then obtain L1/D1. Design L1/D2＝L/D(ETAEM) In which ETAEMIs terminal energy, so when xE＝ETAEMWhen (L)1/D)bsl|E＝E_TAEM＝(LcosσTAEM/D)bsl＝L/D(ETAEM) Can make the terminal tilt angle sigmaTAEMIs 0. By designing the form of equation (20), the roll angle can be controlled to be approximately constant during flight.
Therefore, to solve for L1/D1First, the remaining range x needs to be calculatedDfIn the GER coordinate system, η is defined as a vector pointing from the center of the earth to the aircraft
Figure BDA0002408897380000063
Vector directed to aircraft from earth center
Figure BDA0002408897380000064
The included angle therebetween. The remaining range can be expressed as:
xDf＝Reη-STAEM(21)
wherein S isTAEMIndicating the distance to the target at the end of the reentry segment.
Solving for L1/D1The process of (a) can be described as follows:
when E is more than or equal to EαWhen there is
xD(Eα,E)+xD(ETAEM,Eα)＝xDf(22)
Wherein the content of the first and second substances,
Figure BDA0002408897380000071
wherein the content of the first and second substances,
Figure BDA0002408897380000072
Figure BDA0002408897380000073
and is provided with a plurality of groups of the following components,
Figure BDA0002408897380000074
wherein the content of the first and second substances,
Figure BDA0002408897380000075
c2＝(a1-c1d0)/d1(28)
c3＝a0-c2d0(29)
Figure BDA0002408897380000076
Figure BDA0002408897380000077
Figure BDA0002408897380000078
the formula (22) is solved to obtain
Figure BDA0002408897380000079
Wherein the content of the first and second substances,
Figure BDA0002408897380000081
Figure BDA0002408897380000082
when E < EαIn the meantime, the aircraft is mainly in the altitude adjustment phase and does not update L any more1/D1The value of (c).
3. Reference lateral lift-drag ratio and reference roll angle design
The modulus of the transverse lift-drag ratio can be well planned from the previous (L/D)bslAnd longitudinal lift-to-drag ratio (L)1/D)bslIs determined as
Figure BDA0002408897380000083
In order to satisfy the terminal traverse constraint, the transverse lift-drag ratio profile is designed as follows:
Figure BDA0002408897380000084
wherein E isBR1And EBR2Respectively, the first and second roll reversal points. sgn is a sign variable representing the initial roll reversal direction. After obtaining the transverse lift-drag ratio section, referring to the roll angle section sigmabslIt was also determined that this is expressed as follows:
Figure BDA0002408897380000085
wherein, (L/D)realShowing the actual lift-to-drag ratio profile.
4. Roll reversal point correction algorithm
In this algorithm, the roll reversal point is designed based on a course analytic solution to eliminate dead-end course errors. Here two tilting reversal points are arranged. The first reversal point E is first designed and correctedBR1. First order second reversal point EBR2＝EαThen the terminal traverse xC(ETAEME) is only EBR1As a function of (c). As can be seen from equation (17), the terminal traverse xCfCan be represented as EBR1The piecewise function of (c) is as follows:
Figure BDA0002408897380000086
f, G is a custom function, and the expression is as follows
Figure BDA0002408897380000091
Figure BDA0002408897380000092
Wherein, the parameter xE2,xE1Expressing upper and lower energy bounds; k is a constant of 1 to 5; p1Representing polynomial fitting coefficients.
Obtaining E from equation (38)BR1Derivative of (2) can be obtained
Figure BDA0002408897380000093
Wherein E ═ EBR1Upper mark for right limit of+Is represented by E ═ EBR1Upper mark for left limit of department-And (4) showing. Because L is2/DbslWhen E is equal to EBR1Is not continuous in place, so there is
Figure BDA0002408897380000094
i is 0,1,2,3, 4. In order to make the terminal traverse error 0, x needs to be satisfiedCf(EBR1) 0. The essence of this algorithm is therefore to solve the non-linear equation xCf(EBR1) Root 0, the newton method is as follows:
Figure BDA0002408897380000095
parameter(s)
Figure BDA0002408897380000096
And
Figure BDA0002408897380000097
and not the same as the above-mentioned case,
Figure BDA0002408897380000098
represents EBR1To the power of k in the first order,
Figure BDA0002408897380000099
e representing the kth iterationBR1The value is obtained. Solution stop set to
Figure BDA00024088973800000910
Thus xCf(EBR1) With EBR1Monotonously changing, the solution of this non-linear equation can be obtained by a few simple iterations.
5. Design of command attack angle and inclination angle
To ensure the tracking accuracy of the reference profile, ballistic damping control techniques are employed to suppress long-period, weakly damped oscillations in the presence of reentry ballisticscmdAnd the command roll angle sigmacmdThe following were used:
αcmd＝αbsl+cosσbslKγ(γSG-γ) (43)
Figure BDA00024088973800000911
wherein, KγIs a feedback coefficient, and takes the value of 5; gamma raySGIs an approximate solution to the angle of inclination of the smooth glide trajectory, which can be expressed as follows
Figure BDA0002408897380000101
Wherein S is the aircraft reference area; cLIs the coefficient of lift; ρ is the air density.
The steady glide phase needs to satisfy the process constraint, so the process constraint of equations (8) - (10) is converted into the roll angle constraint, which is expressed as follows
Figure BDA0002408897380000102
Wherein L ismaxFor maximum lift, the expression is as follows:
Figure BDA0002408897380000103
Figure BDA0002408897380000104
Figure BDA0002408897380000105
wherein HminThe lowest boundary of the corresponding glide height can be obtained through process constraint; constant coefficient kσ-50. Thus, to satisfy process constraints, the roll angle needs to satisfy the following conditions:
|σcmd|≤σmax(50)
step five: a height adjustment section guidance method;
1. and correcting the command attack angle and the second inversion point.
(1) Command angle of attack α2Correction algorithm
To meet the terminal height constraint, a ballistic simulation is used to modify the commanded angle of attack α2. Terminal height H obtained from trajectory simulationfAnd terminal constraint height HTAEMThe magnitude of the deviation therebetween, to further correct the resulting final commanded angle of attack α2. In the following variables, the subscript f denotes the terminal value calculated from the ballistic simulation, and the subscript TAEM denotes the terminal constraint value.
Due to gammaf,σfAnd
Figure BDA0002408897380000106
may be approximated as 0 and thus may be respectively falseLet cos gammaf＝1,cosσ f1 and
Figure BDA0002408897380000107
therefore, according to equation (6):
Figure BDA0002408897380000111
wherein, CL(est)(αbsl(Ef) Q) represents the lift coefficient, q, corresponding to the terminal angle of attack in a ballistic simulationfIndicating the terminal dynamic pressure.
According to terminal constraints, have
Figure BDA0002408897380000112
Wherein the content of the first and second substances,
Figure BDA0002408897380000113
indicating the corrected command angle of attack, phif、ψf、Vf、HfRespectively representing terminal latitude, course angle, speed and height obtained by trajectory simulation calculation; q. q.sTAEM、φTAEM、ψTAEM、VTAEM、HTAEMRespectively representing the dynamic pressure value, the latitude value, the course angle value, the speed value and the height value of the terminal constraint.
At phif＝φTAEMAnd psif＝ψTAEMUnder the assumption of (1), linearizing
Figure BDA0002408897380000114
By making equation (51) equal to equation (52), it can be obtained:
Figure BDA0002408897380000115
wherein E isfRepresenting a terminal energy value obtained by ballistic simulation calculation;
Figure BDA0002408897380000116
representing the derivative of the lift coefficient with respect to angle of attack, α solved for in equation (53)2The control device can be taken into a formula (19) to be used as a reference attack angle for guidance control.
(2) Second roll reversal point EBR2Correction algorithm
And in the height adjusting stage, guidance is performed by adopting a proportional guidance law. Through the analysis of the guidance algorithm, if E is reducedBR2The guidance law command will generate a larger roll angle to eliminate a larger course error generated during the last roll reversal, so that the longitudinal lift-drag ratio is reduced, and the terminal speed V is reducedf. Thus, VfCan be regarded as relating to EBR2Is a monotonic function of (a). To guarantee the terminal constraints, this problem can be considered as solving a non-linear equation Vf(EBR2)＝VTAEMThe problem of (2). The solution is carried out by adopting a secant method, and V is simulated by a plurality of trajectoriesfMaking a prediction and adjusting E according to the deviation of the predicted value from the terminal constraintBR2The value of (c).
Figure BDA0002408897380000117
EBR2The condition for ending the iterative correction is that the calculated deviation satisfies
Figure BDA0002408897380000121
2. Reference roll angle
After the second roll reversal, the aircraft enters a height adjustment phase. At this stage, the reference angle of attack is still determined by equation (19), and the reference roll angle is controlled by the proportional guidance law. Angular rate of view
Figure BDA0002408897380000122
Can be expressed as
Figure BDA0002408897380000123
Then, the lateral maneuvering acceleration command a obtained by the proportional guidanceL2As follows
Figure BDA0002408897380000124
Wherein k isPNIs a transverse acceleration maneuvering coefficient and takes the value as
Figure BDA0002408897380000125
Figure BDA0002408897380000126
The longitudinal stroke at the time of the second tilt reversal is shown. Assuming that the aircraft is approximately gliding smoothly, the resultant acceleration in the longitudinal plane is 0. Thus, the longitudinal maneuvering acceleration command aL1Is composed of
Figure BDA0002408897380000127
So that the reference roll angle can be obtained as
Figure BDA0002408897380000128
3. Command angle of attack and roll angle
To ensure terminal velocity constraints, the energy-independent residual glide distance reference profile is tracked by adjusting the angle of attack
Figure BDA0002408897380000129
The commanded angle of attack and roll angle may be expressed as follows:
Figure BDA00024088973800001210
Figure BDA00024088973800001211
wherein k isαThe value of (A) can be adjusted in time according to the task, and can be selected as 2 pi/(1.8 × 10)7)；
Figure BDA00024088973800001212
Representing the remaining range, the value can be obtained by the last previous ballistic simulation.
To meet the process constraints, the roll angle still needs to satisfy the following requirements, as before:
|σcmd|≤σmax(61)
through the five steps, the hypersonic reentry guidance method based on the steady gliding trajectory analytic solution can be obtained.
The invention relates to a hypersonic reentry guidance method based on steady glide trajectory analytic solution, which has the advantages that:
(1) the longitudinal and transverse courses are predicted and corrected based on the three-dimensional trajectory analytic solution, so that the calculation amount is greatly reduced, the calculation speed is improved, and the method is convenient to apply to an airborne system in real time.
(2) The transverse movement is determined by twice tilting reversal points, the position of the first reversal point is corrected to control the terminal traverse error, and the position of the second reversal point is corrected to meet the terminal speed constraint. Compared with the traditional method for controlling the transverse motion by using the course error threshold, the method reduces the tilting and reversing times, improves the energy utilization efficiency and reduces the requirement on a control system.
(3) The height adjusting section utilizes the idea of proportional guidance, so that the terminal roll angle is converged to 0, and the target hitting effect is enhanced.
Drawings
Fig. 1 is a schematic view of the equatorial rotation coordinate system (GER) of the earth's center and the eastern coordinate system (NED) of the north and east of the local area.
Fig. 2 is a schematic diagram of an auxiliary earth rotation reference system (AGR) and an auxiliary north east down coordinate system (GNED).
FIG. 3 is a latitude and longitude plot for five target states.
Fig. 4 is a height-energy curve for five target states.
Fig. 5 is a velocity-energy curve for five target states.
Fig. 6 is an angle of attack versus energy curve for five target states.
Fig. 7 is a roll angle versus energy curve for five target states.
Figure 8 is a ballistic dip-energy curve for five target states.
FIG. 9 is a heading error-energy curve for five target states.
In the above figures, the symbols and symbols are as follows:
in the context of figure 1 of the drawings,
Figure BDA0002408897380000131
representing the equatorial rotation coordinate system (GER), oxyzRepresents the local North-East down coordinate system (NED), North represents the North direction, East represents the East direction. h is altitude, R is range, λ is longitude, φ is latitude, V is velocity of the vehicle relative to the earth, γ is ballistic inclination, ψ is vehicle heading angle, σ is roll angle, based on local north. In the context of figure 2, it is shown,
Figure BDA0002408897380000132
represents an assisted earth-rotation reference frame (AGR),
Figure BDA0002408897380000133
representing the north, east and down-the-earth coordinate system (NED), Axis represents the earth polar Axis. M is the current position point of the aircraft, and o is the projection point of the current position point of the aircraft on the earth surface. T is the target point, oTIs the projected point of the target point on the surface of the earth. OmegaexIndicating the rotational angular velocity of the earth
Figure BDA0002408897380000141
Component of the axis, ωeyIndicating the rotational angular velocity of the earth
Figure BDA0002408897380000142
Component of the axis, ωzxIndicating the rotational angular velocity of the earth
Figure BDA0002408897380000143
The component of the axis. In FIG. 3, Longitude represents Longitude, Latitude represents Latitude, and T1-T5 represent examples of five different objects. In FIG. 4, Altitude represents height, Energy represents energy, H represents height curve obtained by trajectory simulation, HminRepresenting a minimum height boundary that satisfies the process constraints. In fig. 5, Velocity represents altitude and Energy represents Energy. In FIG. 6, Angle of attach represents Angle of Attack and Energy represents Energy. In FIG. 7, Bank Angle represents roll Angle and Energy represents Energy. In FIG. 8, Flight-pathAngle represents roll angle and Energy represents Energy. In fig. 9, Heading Angle represents a Heading Angle, and Energy represents Energy.
Detailed Description
The invention will be further explained in detail with reference to the drawings and the embodiments.
Aiming at the problems of reentry guidance of the hypersonic aircraft, a high-precision hypersonic reentry guidance method based on a three-dimensional trajectory analytic solution is provided. Firstly, an auxiliary geocentric rotation coordinate system is constructed, a new kinetic model of the longitudinal course, the transverse course and the course angle with respect to the energy is established on the basis of the auxiliary geocentric rotation coordinate system, the model is favorable for decoupling and linearizing the highly nonlinear and strongly coupled kinetic model to obtain a three-dimensional ballistic analytic solution, and a guidance method is designed on the basis of the analytic solution. The guidance method utilizes a longitudinal analytic solution to plan a longitudinal reference profile in real time, utilizes a transverse analytic solution to correct a first inversion point to eliminate a transverse error, utilizes trajectory simulation for several times to correct a second inversion point to ensure terminal speed constraint and corrects the minimum value of an attack angle to ensure terminal height constraint.
The invention discloses a hypersonic reentry guidance method based on a steady glide trajectory analytic solution, which is characterized by comprising the following steps: the method comprises the following steps:
the method comprises the following steps: establishing an aircraft dynamics model;
as shown in fig. 1, the geocentric equatorial rotation coordinate system (GER coordinate system) is established: origin at geocentric E, zeThe axis is perpendicular to the earth's equatorial surface, pointing to the north pole; x is the number ofeAxis and yeThe axes are in the equatorial plane and perpendicular to each other. The coordinate system rotates with the earth, which rotates around zeThe angular velocity of rotation of the shaft is the rotational angular velocity omega of the earthe。
North east local coordinate system (NED coordinate system): defining a vertical projection point of an origin o to the ground at the mass center M of the aircraft; the x axis points to the local north, the y axis points to the local east, and the z axis points to the geocentric vertically and downwardly.
Under the background of the rotating earth, a standard atmosphere model and an inverse square gravitational field model are adopted to establish a six-degree-of-freedom kinetic equation of the hypersonic aerocraft as follows:
Figure BDA0002408897380000151
Figure BDA0002408897380000152
Figure BDA0002408897380000153
Figure BDA0002408897380000154
Figure BDA0002408897380000155
Figure BDA0002408897380000156
Figure BDA0002408897380000157
where d/dt represents the derivative over time, t is time, h is altitude, R is range, λ is longitude, φ is latitude, V is the velocity of the aircraft relative to the earth, γ is the ballistic inclination, ψ is the aircraft heading angle, based on local north, σ is the roll angle, R is the altitude, and R is the altitude, whereeIs the radius of the earth, and takes the value of 6378.137km, omegaeIs the angular velocity of rotation of the earth, L ═ ρ V2SCLThe lift acceleration is given by/2 m, and D is equal to rho V2SCDThe/2 m is resistance acceleration; s is the aircraft reference area; cLIs the coefficient of lift; ρ is the air density.
Step two: establishing a constraint model;
the reentry aircraft needs to meet the process constraints of heat flux density, dynamic pressure, overload and the like in the flight process, and the process constraints are expressed as follows:
Figure BDA0002408897380000161
Figure BDA0002408897380000162
Figure BDA0002408897380000163
wherein, KQIs the heat flux density coefficient;
Figure BDA0002408897380000164
representing a maximum allowable heat flow density value; n ismaxRepresents a maximum allowable overload value; q. q.smaxIn addition, due to capacity constraints of the flight control system, the rate of change of the control quantities (angle of attack α and roll angle σ) also satisfies the constraint:
Figure BDA0002408897380000165
in addition to process constraints, the aircraft should also satisfy terminal constraints, with: vTAEM＝2000m/s,HTAEM＝25km,Rtm＝100km,σTAEM≈0°,γTAEM≈0°andΔψTAEM≈0°.
Step three: an initial descent stage guidance method;
according to the characteristics of reentry trajectory, the reentry trajectory is divided into three parts: an initial descent section, a steady glide section and a height adjustment section. The initial descent segment is designed to allow a more gradual transition of the trajectory to a smooth glide segment. The whole flight trajectory is mostly in a stable gliding section, which not only determines the range, speed and height of the flight terminal, but also determines whether the aircraft can meet various process constraints. The height adjustment segment is the last phase of the trajectory, as to whether the trajectory can strictly meet the terminal constraints. The flying speed and flying height of the aircraft are low in this section, and severe requirements are provided for a guidance method for designing a height adjusting section.
In the initial descent, the variation of the controlled variable has little effect on the flight, so the fixed values of angle of attack and roll angle are designed here to avoid altitude variationsmaxFly 20deg and zero roll angle. Rate of change of ballistic inclination
Figure BDA0002408897380000166
At this time, it means that the lift force can satisfy the requirement of smooth gliding, so that the attack angle can be smoothly transited from the maximum attack angle to the reference attack angle αbsl。
At this stage, the command angle of attack and roll are designed as follows:
Figure BDA0002408897380000167
σcmd＝0 (12)
wherein the content of the first and second substances,
Δγ＝γ-γSG(13)
Figure BDA0002408897380000171
wherein, γSGRepresenting the angle of inclination, k, of the smooth gliding trajectoryγRepresenting the feedback coefficient and here taking the value 5. When Δ γ is 0, the aircraft enters a smooth glide phase.
Step four: a steady glide phase guidance method;
ballistic planning in the steady glide phase is the basis of the proposed guidance algorithm, which requires approaching the aircraft to the target while guaranteeing process constraints. The algorithm determines the longitudinal range and lateral range profiles based on a new three-dimensional analytical solution. Here, a new ballistic analytical solution is first briefly introduced.
As shown in fig. 2, in order to meet the characteristic that the reentry trajectory has a greater longitudinal range than transverse range, two coordinate systems are first established, and the meanings of the longitudinal range and the transverse range can be visually expressed.
Defining an auxiliary earth-centered rotational reference system (AGR coordinate system) rotating with the earth: the origin is at the geocentric E,
Figure BDA0002408897380000172
the shaft points to the initial position of the aircraft,
Figure BDA0002408897380000173
the axis is in the plane of a great circle passing through the aircraft and the target point and perpendicular to the aircraft
Figure BDA0002408897380000174
The shaft is provided with a plurality of axial holes,
Figure BDA0002408897380000175
the axes may be determined according to the right hand rule.
Meanwhile, a local generalized north east down coordinate system (GNED coordinate system) is defined: the vertical projection point of the origin o from the mass center M of the aircraft to the ground;
Figure BDA0002408897380000176
the axis is vertically downwards directed to the center of the earth,
Figure BDA0002408897380000177
axis-oriented AGR coordinate system
Figure BDA0002408897380000178
The direction of the "north" direction is,
Figure BDA0002408897380000179
the axes are determined by the right hand rule.
Define the longitudinal extent as
Figure BDA00024088973800001710
Define the course as
Figure BDA00024088973800001711
Wherein
Figure BDA00024088973800001712
And
Figure BDA00024088973800001713
respectively representing a generalized longitude and a generalized latitude in an AGR coordinate system.
Defining the energy as
Figure BDA00024088973800001714
The generalized longitude and latitude coordinate system in the AGR coordinate axis is beneficial to linearizing and decoupling a complex dynamic model, and the following three-dimensional ballistic analysis solution is obtained.
Figure BDA0002408897380000181
Figure BDA0002408897380000182
Figure BDA0002408897380000183
Wherein the content of the first and second substances,
Figure BDA0002408897380000184
R*＝Re+H；
Figure BDA0002408897380000185
hij,i＝1,2；j＝0,1,α1and pijI is 1, 2; j is 0,1,2,3,4 are constant coefficients. Since the analytic solution has high precision, the analytic solution can be applied to guidance method derivation to reduce calculation amount and improve universality.
6. Design of reference angle of attack
Because the analytic solution is designed by taking energy as an independent variable, all reference profiles in the guidance method also adopt energy as the independent variable so as to more conveniently utilize the analytic solution. The reference angle of attack profile is designed as:
Figure BDA0002408897380000186
wherein E isα＝-5.55×107J/kg is located near the intersection point of the smooth gliding section and the final altitude adjusting section, in order to exert the maximum capability of the aircraft, the attack angle is designed to be α1The quadratic function for which the angle of attack is designed to be energy is designed so that the angle of attack can be from α1Gently transition to α2In this guidance method α is preset for the first time2The ballistic simulation is later applied to correct it to strictly meet the terminal height constraint, 6 deg. When the design of the attack angle section is finished, the corresponding lift-drag ratio section L/DbslAs determined accordingly.
7. Baseline lift-to-drag ratio profile design
In order to meet the requirement of range, the longitudinal lift-drag ratio profile is designed as follows
Figure BDA0002408897380000191
Wherein, the length range resolution is equal to the remaining range, so as to obtain L1/D1. Design L1/D2＝L/D(ETAEM) In order to make the terminal tilt angle 0. By designing the form of equation (20), the roll angle can be controlled to be approximately constant during flight.
Therefore, to solve for L1/D1First, the remaining range x needs to be calculatedDfIn the GER coordinate system, η is defined as a vector pointing from the center of the earth to the aircraft
Figure BDA0002408897380000192
Vector directed to aircraft from earth center
Figure BDA0002408897380000193
The included angle therebetween. ThenThe remaining range can be expressed as:
xDf＝Reη-STAEM(21)
wherein S isTAEMIndicating the distance to the target at the end of the reentry segment.
Solving for L1/D1The process of (a) can be described as follows:
when E is more than or equal to EαWhen there is
xD(Eα,E)+xD(ETAEM,Eα)＝xDf(22)
Wherein the content of the first and second substances,
Figure BDA0002408897380000194
wherein the content of the first and second substances,
Figure BDA0002408897380000195
Figure BDA0002408897380000196
and is provided with a plurality of groups of the following components,
Figure BDA0002408897380000197
wherein the content of the first and second substances,
Figure BDA0002408897380000198
c2＝(a1-c1d0)/d1(28)
c3＝a0-c2d0(29)
Figure BDA0002408897380000201
Figure BDA0002408897380000202
Figure BDA0002408897380000203
the formula (22) is solved to obtain
Figure BDA0002408897380000204
Wherein the content of the first and second substances,
Figure BDA0002408897380000205
Figure BDA0002408897380000206
when E < EαIn the meantime, the aircraft is mainly in the altitude adjustment phase and does not update L any more1/D1The value of (c).
8. Reference lateral lift-drag ratio and reference roll angle design
The module value of the transverse lift-drag ratio can be well planned by the previously planned L/DbslAnd a longitudinal lift-to-drag ratio L1/DbsIs determined as
Figure BDA0002408897380000207
In order to satisfy the terminal traverse constraint, the transverse lift-drag ratio profile is designed as follows:
Figure BDA0002408897380000208
wherein E isBR1And EBR2Respectively, the first and second roll reversal points. sgn is a sign variable representing the initial roll reversal direction. After obtaining the transverse lift-drag ratio profile, the roll angle profile is also determined as follows:
Figure BDA0002408897380000211
wherein, L/DrealShowing the actual lift-to-drag ratio profile.
9. Roll reversal point correction algorithm
In this algorithm, the roll reversal point is designed based on a course analytic solution to eliminate dead-end course errors. Here two tilting reversal points are arranged. The first reversal point E is first designed and correctedBR1. First order second reversal point EBR2＝EαThen the terminal traverse xC(ETAEME) is only EBR1As a function of (c). As can be seen from equation (17), the terminal traverse can be represented as EBR1The piecewise function of (c) is as follows:
Figure BDA0002408897380000212
wherein the content of the first and second substances,
Figure BDA0002408897380000213
Figure BDA0002408897380000214
obtaining E from equation (38)BR1Derivative of (2) can be obtained
Figure BDA0002408897380000215
Wherein E ═ EBR1Upper mark for right limit of+Is represented by E ═ EBR1Upper mark for left limit of department-And (4) showing. Because L is2/DbslWhen E is equal to EBR1Is not continuous in place, so there is
Figure BDA0002408897380000216
i is 0,1,2,3, 4. In order to make the terminal traverse error 0, x needs to be satisfiedCf(EBR1) 0. The essence of this algorithm is therefore to solve the non-linear equation xCf(EBR1) Root 0, the newton method is as follows:
Figure BDA0002408897380000217
solution stop set to
Figure BDA0002408897380000218
Thus xCf(EBR1) With EBR1Monotonously changing, the solution of this non-linear equation can be obtained by a few simple iterations.
10. Design of command attack angle and inclination angle
To ensure the tracking accuracy of the reference profile, ballistic damping control techniques are employed herein to suppress long-period, weakly damped oscillations present in the reentry trajectory. The design command attack angle and roll angle are as follows:
αcmd＝αbsl+cosσbslKγ(γSG-γ) (43)
Figure BDA0002408897380000221
wherein, KγIs a feedback coefficient, and takes the value of 5; gamma raySGIs an approximate solution to the angle of inclination of the smooth glide trajectory, which can be expressed as follows
Figure BDA0002408897380000222
The steady glide phase needs to satisfy the process constraint, so the process constraint of equations (8) - (10) is converted into the roll angle constraint, which is expressed as follows
Figure BDA0002408897380000223
Wherein the content of the first and second substances,
Figure BDA0002408897380000224
Figure BDA0002408897380000225
Figure BDA0002408897380000226
wherein HminThe lowest boundary of the corresponding glide height can be obtained through process constraint; constant coefficient kσ-50. Thus, to satisfy process constraints, the roll angle needs to satisfy the following conditions:
|σcmd|≤σmax(50)
step five: a height adjustment section guidance method;
4. and correcting the command attack angle and the second inversion point.
(3) Command angle of attack α2Correction algorithm
To meet the terminal height constraint, a ballistic simulation is used to modify the commanded angle of attack α2. Terminal height H obtained from trajectory simulationfAnd terminal constraint height HTAEMThe magnitude of the deviation therebetween, to further correct the resulting final commanded angle of attack α2. In the following variables, the subscript f denotes the terminal value calculated from the ballistic simulation, and the subscript TAEM denotes the terminal constraint value.
Due to gammaf,σfAnd
Figure BDA0002408897380000231
can be approximated to 0, and thus can assume cos γ, respectivelyf＝1,cosσ f1 and
Figure BDA0002408897380000232
therefore, according to equation (6):
Figure BDA0002408897380000233
wherein, CL(est)(αbsl(Ef) Q) represents the lift coefficient, q, corresponding to the terminal angle of attack in a ballistic simulationfIndicating the terminal dynamic pressure.
According to terminal constraints, have
Figure BDA0002408897380000234
Wherein the content of the first and second substances,
Figure BDA0002408897380000235
representing the corrected commanded angle of attack.
At phif＝φTAEMAnd psif＝ψTAEMUnder the assumption of (1), linearizing
Figure BDA0002408897380000236
By making equation (51) equal to equation (52), it can be obtained:
Figure BDA0002408897380000237
α solved by equation (53)2The control device can be taken into a formula (19) to be used as a reference attack angle for guidance control.
(4) Second roll reversal point EBR2Correction algorithm
And in the height adjusting stage, guidance is performed by adopting a proportional guidance law. Through the analysis of the guidance algorithm, if E is reducedBR2The guidance law command will generate a larger roll angle to eliminate a larger course error generated during the last roll reversal, so that the longitudinal lift-drag ratio is reduced, and the terminal speed V is reducedf. Thus, VfCan be regarded as relating to EBR2Is a monotonic function of (a). To guarantee the terminal constraints, this problem can be considered as solving a non-linear equation Vf(EBR2)＝VTAEMThe problem of (2). The solution is carried out by adopting a secant method, and V is simulated by a plurality of trajectoriesfMaking a prediction and adjusting E according to the deviation of the predicted value from the terminal constraintBR2The value of (c).
Figure BDA0002408897380000241
EBR2The condition for ending the iterative correction is that the calculated deviation satisfies
Figure BDA0002408897380000242
5. Reference roll angle
After the second roll reversal, the aircraft enters a height adjustment phase. At this stage, the reference angle of attack is still determined by equation (19), and the reference roll angle is controlled by the proportional guidance law. The line-of-sight angular rate may be expressed as
Figure BDA0002408897380000243
Then, the lateral maneuver acceleration command obtained by the proportional steering is as follows
Figure BDA0002408897380000244
Wherein the content of the first and second substances,
Figure BDA0002408897380000245
Figure BDA0002408897380000246
the longitudinal stroke at the time of the second tilt reversal is shown. Assuming that the aircraft is approximately gliding smoothly, the resultant acceleration in the longitudinal plane is 0. Thus, there are
Figure BDA0002408897380000247
So that the reference roll angle can be obtained as
Figure BDA0002408897380000248
6. Command angle of attack and roll angle
To ensure terminal velocity constraints, the energy-independent residual glide distance reference profile is tracked by adjusting the angle of attack
Figure BDA0002408897380000249
Command angle of attack and roll angleTo represent the following:
Figure BDA00024088973800002410
Figure BDA00024088973800002411
wherein k isαThe value of (A) can be adjusted in time according to the task, and can be selected as 2 pi/(1.8 × 10)7)；
Figure BDA00024088973800002412
The value of (c) can be obtained by the last previous ballistic simulation.
To meet the process constraints, the roll angle still needs to satisfy the following requirements, as before:
|σcmd|≤σmax(61)
through the five steps, the hypersonic reentry guidance method based on the steady gliding trajectory analytic solution can be obtained.
To demonstrate the universality of the guidance method, 5 examples of different target positions are considered here. Five different target points are indicated by T1 to T5, respectively, which are set at (-77 °, 0 °) (T1), (-34 °, -32.5 °) (T2), (5 °, -50 °) (T3), (65 °, -33 °) (T4), and (118 °, 0 °) (T5), respectively. The initial conditions are all h0＝80km,V0＝7200m/s,λ0＝0°,θ045 ° and γ0＝0°。ψ0The pointing direction is determined by the target position. The terminal constraint conditions are all VTAEM＝2000m/s,HTAEM25km and R tm100 km. The simulation results are shown in fig. 3 to 9.
In fig. 3, a red circle represents a range of 100km from the target. As is apparent from the figure, all the trajectories fall on the circles, and the guidance method is proved to meet the terminal distance constraint. Table 1 shows the state quantities of range, speed, altitude, roll angle, trajectory inclination and course error of other terminals.
TABLE 1 simulation results
Figure BDA0002408897380000251
It can be seen from fig. 4 that all but two westward-fly examples satisfy the terminal process constraints. This is because when the aircraft flies westward, the coriolis force is applied downward, so that the glide height is reduced and the heat flux density is increased. Furthermore, from fig. 4-9 and table 1, we can readily see that all examples satisfy Vf≈VTAEM＝2000m/s,Hf≈HTAEM＝25km,Rtm＝100km,σf≈0°,γf0DEG and [ Delta ] phifTerminal constraint of 0 deg..
In conclusion, the method is deduced through the steps, namely the hypersonic velocity reentry guidance method based on the steady glide trajectory analytic solution, case simulation results show that the method can meet various constraints required by reentry guidance, and the method is designed based on the trajectory analytic solution, has the characteristics of high precision and small calculated amount, and can be applied to an airborne system of an aircraft in real time.

Claims (2)

1. A hypersonic velocity reentry guidance method based on steady glide trajectory analytic solution is characterized in that: the method comprises the following steps:
the method comprises the following steps: establishing an aircraft dynamics model;
under the background of the rotating earth, a standard atmosphere model and an inverse square gravitational field model are adopted to establish a six-degree-of-freedom kinetic equation of the hypersonic aerocraft as follows:
Figure FDA0002408897370000011
Figure FDA0002408897370000012
Figure FDA0002408897370000013
Figure FDA0002408897370000014
Figure FDA0002408897370000015
Figure FDA0002408897370000016
Figure FDA0002408897370000017
where d/dt represents the derivative over time, t is time, h is altitude, R is range, λ is longitude, φ is latitude, V is the velocity of the aircraft relative to the earth, γ is the ballistic inclination, ψ is the aircraft heading angle, based on local north, σ is the roll angle, R is the altitude, and R is the altitude, whereeIs the radius of the earth, and takes the value of 6378.137km, omegaeIs the angular velocity of rotation of the earth, g represents the acceleration of gravity; l ═ ρ V2SCLThe lift acceleration is given by/2 m, and D is equal to rho V2SCD(2 m) is the acceleration of resistance, CL、CDRespectively a lift coefficient and a drag coefficient;
step two: establishing a constraint model;
the reentry vehicle needs to meet the process constraints of heat flux density, dynamic pressure and overload during flight, as follows:
Figure FDA0002408897370000021
Figure FDA0002408897370000022
Figure FDA0002408897370000023
parameter p0Denotes nominal atmospheric density, e denotes natural constant, h denotes altitude, hSRepresenting nominal height, n representing overload, m representing mass;
wherein, KQIs the heat flux density coefficient;
Figure FDA0002408897370000024
representing a maximum allowable heat flow density value; n represents an allowable overload value; n ismaxRepresents a maximum allowable overload value; q represents an allowable dynamic pressure value; q. q.smaxRepresents the maximum allowable dynamic pressure value;
furthermore, due to the capability constraints of the flight control system, the rate of change of the control quantity also satisfies the constraint:
Figure FDA0002408897370000025
Figure FDA0002408897370000026
in addition to process constraints, the aircraft should also satisfy terminal constraints, with: vTAEM＝2000m/s,HTAEM＝25km,Rtm＝100km,σTAEM≈0°,γTAEM≈0°andΔψTAEM0 ° where the subscript TAEM denotes the terminal amount; vTAEM、HTAEM、σTAEM、γTAEM、ΔψTAEMRespectively representing the terminal speed, the terminal height, the terminal roll angle, the terminal trajectory inclination angle and the terminal course angle error; rtmRepresenting the shot-eye distance;
step three: an initial descent stage guidance method;
during the initial descent, the maximum angle of attack α is selectedmaxFlying at 20deg and zero roll angle; rate of change of ballistic inclination
Figure FDA0002408897370000027
By time, it is meant that lift meets the smooth glide requirements, so that the angle of attack transitions smoothly from the maximum angle of attack to the reference angle of attack profile αbsl；
At this stage, the command angle of attack and roll are designed as follows:
Figure FDA0002408897370000028
σcmd＝0 (12)
wherein the content of the first and second substances,
Δγ＝γ-γSG(13)
Figure FDA0002408897370000031
wherein, γSGRepresenting the angle of inclination, k, of the smooth gliding trajectoryγRepresenting the feedback coefficient and taking 5, when delta gamma is 0, the aircraft enters the steady glide phase, and the parameter αcmdIs an instruction angle of attack; sigmacmdIs a commanded roll angle; delta gamma is the deviation of the current trajectory inclination angle and the steady glide trajectory inclination angle; delta gamma1As the rate of change of ballistic inclination
Figure FDA0002408897370000032
A deviation value of the trajectory inclination angle when it is 0 from the steady glide trajectory inclination angle;
step four: a steady glide phase guidance method;
ballistic planning in the steady glide phase is the basis of the proposed guidance algorithm, and the aircraft needs to approach the target while process constraints are ensured;
firstly, establishing the following two coordinate systems, and intuitively obtaining the meanings of representing the longitudinal course and the transverse course;
defining an auxiliary earth center rotation reference system (AGR coordinate system) rotating along with the earth: the origin is at the geocentric E,
Figure FDA0002408897370000033
the shaft points to the initial position of the aircraft,
Figure FDA0002408897370000034
with the axis in over-flightIn the plane of the great circle of the target point and perpendicular to
Figure FDA0002408897370000035
The shaft is provided with a plurality of axial holes,
Figure FDA0002408897370000036
the axis is determined according to the right-hand rule;
meanwhile, a local generalized north-east coordinate system, i.e., GNED coordinate system, is defined: the vertical projection point of the origin o from the mass center M of the aircraft to the ground;
Figure FDA0002408897370000037
the axis is vertically downwards directed to the center of the earth,
Figure FDA0002408897370000038
axis-oriented AGR coordinate system
Figure FDA0002408897370000039
The direction of the "north" direction is,
Figure FDA00024088973700000310
the axis is determined by the right hand rule;
define the longitudinal extent as
Figure FDA00024088973700000311
Define the course as
Figure FDA00024088973700000312
Wherein
Figure FDA00024088973700000313
And
Figure FDA00024088973700000314
respectively representing a generalized longitude and a generalized latitude in an AGR coordinate system;
defining the energy as
Figure FDA00024088973700000315
Defining the longitudinal lift-to-drag ratio Lcos sigma/D ═ L1D, transverse lift-drag ratio Lsin sigma/D ═ L2D; the generalized longitude and latitude coordinate system in the AGR coordinate axis is beneficial to linearizing and decoupling a complex dynamic model, and the following three-dimensional ballistic analysis solution is obtained;
Figure FDA0002408897370000041
Figure FDA0002408897370000042
Figure FDA0002408897370000043
wherein the content of the first and second substances,
Figure FDA0002408897370000044
R*＝Re+H；
Figure FDA0002408897370000045
hij,i＝1,2；j＝0,1,α1and pijI is 1, 2; j is 0,1,2,3 and 4 are constant coefficients; e, E0Respectively representing the current energy and the initial energy value; x is the number ofERepresenting an energy integral variable; l is1Representing lift values in a longitudinal plane; x is the number ofC0Representing an initial traverse value;
Figure FDA0002408897370000046
representing an initial course angle error under a GNED coordinate system;
Figure FDA0002408897370000047
representing the course angle error under the GNED coordinate system;
4.1 datum Angle of attack design
All reference profiles in the guidance method also use energy as an independent variable, a reference attack angle profile αbslThe design is as follows:
Figure FDA0002408897370000048
wherein E isα＝-5.55×107J/kg is located near the boundary point between the smooth gliding section and the final height adjusting section, and in order to exert the maximum capability of the aircraft, the attack angle is designed to be α1Designing the quadratic function of the angle of attack as energy to make the angle of attack from α deg1Gently transition to α2In this guidance method α is preset for the first time2If the terminal height is 6deg, trajectory simulation is applied to correct the terminal height to strictly meet the terminal height constraint; when the design of the attack angle section is finished, the corresponding reference lift-drag ratio section (L/D)bslIs also determined accordingly;
4.2 benchmark lift-drag ratio profile design
In order to meet the requirement of range, the longitudinal lift-drag ratio profile is designed as follows
Figure FDA0002408897370000051
Wherein, let the terminal longitudinal analysis solve xDfEqual to the remaining range, i.e. obtaining L1/D1(ii) a Design L1/D2＝L/D(ETAEM) In which ETAEMIs terminal energy, so when xE＝ETAEMWhen (L)1/D)bsl|E＝E_TAEM＝(LcosσTAEM/D)bsl＝L/D(ETAEM) So that the tilt angle sigma of the terminalTAEMIs 0; by means of the form of the design formula (20), the control roll angle is approximately constant during the flight;
therefore, to solve for L1/D1First, the remaining range x needs to be calculatedDfIn the GER coordinate system, η is defined as the vector pointing from the center of the earth to the aircraft
Figure FDA0002408897370000058
Vector directed to aircraft from earth center
Figure FDA0002408897370000053
The included angle between them; the remaining range is then expressed as:
xDf＝Reη-STAEM(21)
wherein S isTAEMIndicating the distance to the target at the end of the reentry segment;
solving for L1/D1The process of (a) is described as follows:
when E is more than or equal to EαWhen there is
xD(Eα,E)+xD(ETAEM,Eα)＝xDf(22)
Wherein the content of the first and second substances,
Figure FDA0002408897370000054
wherein the content of the first and second substances,
Figure FDA0002408897370000055
Figure FDA0002408897370000056
and is provided with a plurality of groups of the following components,
Figure FDA0002408897370000057
wherein the content of the first and second substances,
Figure FDA0002408897370000061
c2＝(a1-c1d0)/d1(28)
c3＝a0-c2d0(29)
Figure FDA0002408897370000062
Figure FDA0002408897370000063
Figure FDA0002408897370000064
the formula (22) is solved to obtain
Figure FDA0002408897370000065
Wherein the content of the first and second substances,
Figure FDA0002408897370000066
Figure FDA0002408897370000067
when E < EαIn the meantime, the aircraft is mainly in the altitude adjustment phase and does not update L any more1/D1A value of (d);
4.3 reference transverse lift-drag ratio and reference roll angle design
The modulus of the transverse lift-drag ratio is well planned from the previous (L/D)bsAnd longitudinal lift-to-drag ratio (L)1/D)bsIs determined as
Figure FDA0002408897370000068
In order to satisfy the terminal traverse constraint, the transverse lift-drag ratio profile is designed as follows:
Figure FDA0002408897370000069
wherein E isBR1And EBR2Respectively representing a first and a second roll reversal point; sgn is a sign variable representing the initial roll reversal direction; after obtaining the transverse lift-drag ratio section, referring to the roll angle section sigmabslIt was also determined that this is expressed as follows:
Figure FDA0002408897370000071
wherein, (L/D)realRepresenting the actual lift-to-drag ratio profile;
4.4 roll reversal Point correction Algorithm
Designing a roll reversal point based on a traverse analytic solution to eliminate a terminal traverse error; arranging two tilting reversal points; the first reversal point E is first designed and correctedBR1(ii) a First order second reversal point EBR2＝EαThen the terminal traverse xC(ETAEME) is only EBR1A function of (a); as seen from the formula (17), the terminal traverse xCfIs represented by EBR1The piecewise function of (c) is as follows:
Figure FDA0002408897370000072
f, G is a custom function, and the expression is as follows
Figure FDA0002408897370000073
Figure FDA0002408897370000074
Wherein, the parameter xE2,xE1Expressing upper and lower energy bounds; k is a constant of 1 to 5; p1Representing a polynomial fitting coefficient;
obtaining E from equation (38)BR1Derivative of (2) can be obtained
Figure FDA0002408897370000075
Wherein E ═ EBR1Upper mark for right limit of+Is represented by E ═ EBR1Upper mark for left limit of department-Represents; because L is2/DbslWhen E is equal to EBR1Is not continuous in place, so there is
Figure FDA0002408897370000076
In order to make the terminal traverse error 0, x needs to be satisfiedCf(EBR1) 0; the essence of this algorithm is therefore to solve the non-linear equation xCf(EBR1) Root 0, the newton method is as follows:
Figure FDA0002408897370000081
parameter(s)
Figure FDA0002408897370000082
And
Figure FDA0002408897370000083
and not the same as the above-mentioned case,
Figure FDA0002408897370000084
represents EBR1To the power of k in the first order,
Figure FDA0002408897370000085
e representing the kth iterationBR1A value; solution stop set to
Figure FDA0002408897370000086
(1 m); thus xCf(EBR1) With EBR1Monotonously changing, and the solution of the nonlinear equation is obtained by simple iterations;
4.5 Command attack Angle and roll Angle design
In order to ensure the tracking precision of the reference profile, a trajectory damping control technology is adopted to inhibit the oscillation of long-period and weak damping existing in a reentry trajectory; design ofCommand angle of attack αcmdAnd the command roll angle sigmacmdThe following were used:
αcmd＝αbsl+cosσbslKγ(γSG-γ) (43)
Figure FDA0002408897370000087
wherein, KγIs a feedback coefficient, and takes the value of 5; gamma raySGIs an approximate solution to the angle of inclination of the smooth glide trajectory, as shown below
Figure FDA0002408897370000088
Wherein S is the aircraft reference area; cLIs the coefficient of lift; ρ is the air density;
the smooth glide phase needs to satisfy the process constraint, so the process constraint of equation (8) -10 is converted into the roll angle constraint, which is expressed as follows
Figure FDA0002408897370000089
Wherein L ismaxFor maximum lift, the expression is as follows:
Figure FDA00024088973700000810
Figure FDA00024088973700000811
Figure FDA00024088973700000812
wherein HminThe lowest boundary corresponding to the glide height is obtained through process constraint; constant coefficient kσ-50; thus, to satisfy process constraints, the roll angle needs to satisfy the following conditions:
|σcmd|≤σmax(50)
step five: a height adjustment section guidance method;
5.1 correcting the command attack angle and the second reversal point;
(1) command angle of attack α2Correction algorithm
To meet the terminal height constraint, a ballistic simulation is used to modify the commanded angle of attack α2(ii) a Terminal height H obtained from trajectory simulationfAnd terminal constraint height HTAEMThe magnitude of the deviation therebetween, to further correct the resulting final commanded angle of attack α2(ii) a In the following variables, subscript f represents a terminal value obtained by ballistic simulation calculation, and subscript TAEM represents a terminal constraint value;
due to gammaf,σfAnd
Figure FDA0002408897370000091
is approximately 0, so cos γ is assumed separatelyf＝1,cosσf1 and
Figure FDA0002408897370000092
therefore, it is obtained from equation (6):
Figure FDA0002408897370000093
wherein, CL(est)(αbsl(Ef) Q) represents the lift coefficient, q, corresponding to the terminal angle of attack in a ballistic simulationfRepresenting terminal dynamic pressure; according to terminal constraints, have
Figure FDA0002408897370000094
Wherein the content of the first and second substances,
Figure FDA0002408897370000095
indicating the corrected command angle of attack, phif、ψf、Vf、HfRespectively representing trajectory simulation calculationThe terminal latitude, course angle, speed and height; q. q.sTAEM、φTAEM、ψTAEM、VTAEM、HTAEMRespectively representing a dynamic pressure value, a latitude value, a course angle value, a speed value and a height value of terminal constraint; at phif＝φTAEMAnd psif＝ψTAEMUnder the assumption of (1), linearizing
Figure FDA0002408897370000096
By making equation (51) equal to equation (52), we obtain:
Figure FDA0002408897370000097
wherein E isfRepresenting a terminal energy value obtained by ballistic simulation calculation;
Figure FDA0002408897370000101
representing the derivative of the lift coefficient with respect to angle of attack, α solved for in equation (53)2Taking the belt-in type (19) as a reference attack angle to conduct guidance control;
(2) second roll reversal point EBR2Correction algorithm
In the height adjustment stage, guidance is designed by adopting a proportional guidance law; through the analysis of the guidance algorithm, if E is reducedBR2The guidance law command will generate a larger roll angle to eliminate a larger course error generated during the last roll reversal, so that the longitudinal lift-drag ratio is reduced, and the terminal speed V is reducedf(ii) a Thus, VfIs regarded as relating to EBR2A monotonic function of (a); to guarantee the terminal constraints, this problem is treated as solving the nonlinear equation Vf(EBR2)＝VTAEMThe problem of the solution of (a); solving by adopting a secant method, and simulating V by a plurality of trajectoriesfMaking a prediction and adjusting E according to the deviation of the predicted value from the terminal constraintBR2A value of (d);
Figure FDA0002408897370000102
EBR2the condition for ending the iterative correction is that the calculated deviation satisfies
Figure FDA0002408897370000103
5.2 reference roll angle
After the second roll reversal, the aircraft enters a height adjustment phase; at this stage, the reference angle of attack is still determined by equation (19), and the reference roll angle is controlled by the proportional guidance law; angular rate of view
Figure FDA0002408897370000104
Is shown as
Figure FDA0002408897370000105
Then, the lateral maneuvering acceleration command a obtained by the proportional guidanceL2As follows
Figure FDA0002408897370000106
Wherein k isPNIs a transverse acceleration maneuvering coefficient and takes the value as
Figure FDA0002408897370000107
Figure FDA0002408897370000108
Represents the longitudinal stroke at the time of the second tilt reversal; assuming that the aircraft glides approximately smoothly, the resultant acceleration in the longitudinal plane is 0; thus, the longitudinal maneuvering acceleration command aL1Is composed of
Figure FDA0002408897370000109
So that the reference roll angle is
Figure FDA00024088973700001010
5.3 Command Angle of attack and roll
To ensure terminal velocity constraints, the energy-independent residual glide distance reference profile is tracked by adjusting the angle of attack
Figure FDA0002408897370000111
The commanded angle of attack and roll angle are expressed as follows:
Figure FDA0002408897370000112
Figure FDA0002408897370000113
wherein k isαThe value of (2) is selected as 2 pi/(1.8 × 10) according to task timely adjustment7)；
Figure FDA0002408897370000114
Representing the remaining range, obtained by the last trajectory simulation before;
to meet the process constraints, the roll angle still needs to satisfy the following requirements, as before:
|σcmd|≤σmax(61)
through the five steps, the hypersonic reentry guidance method based on the steady glide trajectory analytic solution is obtained.
2. The hypersonic re-entry guidance method based on the smooth gliding trajectory analytic solution as claimed in claim 1, characterized in that: the hypersonic flight vehicle is an aircraft with a flight Mach number larger than 5 and capable of realizing hypersonic flight in an atmosphere layer and a trans-atmosphere layer.
CN202010170138.1A 2020-03-12 2020-03-12 Hypersonic velocity reentry guidance method based on steady glide trajectory analytic solution Active CN111306989B (en)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202010170138.1A CN111306989B (en) 2020-03-12 2020-03-12 Hypersonic velocity reentry guidance method based on steady glide trajectory analytic solution

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202010170138.1A CN111306989B (en) 2020-03-12 2020-03-12 Hypersonic velocity reentry guidance method based on steady glide trajectory analytic solution

Publications (2)

Publication Number Publication Date
CN111306989A true CN111306989A (en) 2020-06-19
CN111306989B CN111306989B (en) 2021-12-28

Family

ID=71149646

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202010170138.1A Active CN111306989B (en) 2020-03-12 2020-03-12 Hypersonic velocity reentry guidance method based on steady glide trajectory analytic solution

Country Status (1)

Country Link
CN (1) CN111306989B (en)

Cited By (12)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN111674573A (en) * 2020-08-11 2020-09-18 北京控制与电子技术研究所 Non-parallel gravitational field deep space impact control method and system based on proportional guidance
CN112009698A (en) * 2020-08-28 2020-12-01 北京空间技术研制试验中心 Return flight energy management method for hypersonic cruise aircraft
CN112162567A (en) * 2020-09-09 2021-01-01 北京航空航天大学 Avoidance guidance method suitable for online no-fly zone of aircraft
CN112507466A (en) * 2020-12-21 2021-03-16 中国人民解放军96901部队23分队 Gliding trajectory energy change-dependent degradation solution considering earth rotation rate
CN112507465A (en) * 2020-12-21 2021-03-16 中国人民解放军96901部队23分队 Gliding trajectory degradation solution with energy change based on resistance and longitudinal lift-drag ratio
CN112507467A (en) * 2020-12-21 2021-03-16 中国人民解放军96901部队23分队 Gliding trajectory degradation solution with speed change based on resistance and longitudinal lift-drag ratio
CN112861249A (en) * 2020-12-21 2021-05-28 中国人民解放军96901部队23分队 Gliding trajectory degradation solution along with speed change based on attack angle and resistance acceleration
CN112861250A (en) * 2020-12-21 2021-05-28 中国人民解放军96901部队23分队 Gliding trajectory degradation solution along with energy change based on attack angle and resistance acceleration
CN113375634A (en) * 2021-04-30 2021-09-10 北京临近空间飞行器系统工程研究所 Altitude measurement method based on atmospheric model and aircraft normal overload combination
CN113741509A (en) * 2021-07-30 2021-12-03 北京航空航天大学 Hypersonic glide aircraft downward pressing section energy management method
CN114265433A (en) * 2021-12-24 2022-04-01 北京航空航天大学 Transverse guidance method and system for meeting large-transverse-stroke maneuver
CN114675545A (en) * 2022-05-26 2022-06-28 中国人民解放军火箭军工程大学 Hypersonic aircraft reentry cooperative guidance method based on reinforcement learning

Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US8489258B2 (en) * 2009-03-27 2013-07-16 The Charles Stark Draper Laboratory, Inc. Propulsive guidance for atmospheric skip entry trajectories
CN104035335A (en) * 2014-05-27 2014-09-10 北京航空航天大学 High accuracy longitudinal and cross range analytical prediction method based smooth gliding reentry guidance method
CN104656450A (en) * 2015-01-20 2015-05-27 北京航空航天大学 Method for designing hypersonic aircraft to re-enter into trajectory in smooth gliding mode
CN105550402A (en) * 2015-12-07 2016-05-04 北京航空航天大学 Attack angle or inclination angle frequency conversion based design method for hypersonic steady maneuver gliding trajectory
CN108036676A (en) * 2017-12-04 2018-05-15 北京航空航天大学 A kind of autonomous reentry guidance method of full directive based on three-dimensional resolution Value of Reentry Vehicle

Patent Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US8489258B2 (en) * 2009-03-27 2013-07-16 The Charles Stark Draper Laboratory, Inc. Propulsive guidance for atmospheric skip entry trajectories
CN104035335A (en) * 2014-05-27 2014-09-10 北京航空航天大学 High accuracy longitudinal and cross range analytical prediction method based smooth gliding reentry guidance method
CN104656450A (en) * 2015-01-20 2015-05-27 北京航空航天大学 Method for designing hypersonic aircraft to re-enter into trajectory in smooth gliding mode
CN105550402A (en) * 2015-12-07 2016-05-04 北京航空航天大学 Attack angle or inclination angle frequency conversion based design method for hypersonic steady maneuver gliding trajectory
CN108036676A (en) * 2017-12-04 2018-05-15 北京航空航天大学 A kind of autonomous reentry guidance method of full directive based on three-dimensional resolution Value of Reentry Vehicle

Cited By (21)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN111674573A (en) * 2020-08-11 2020-09-18 北京控制与电子技术研究所 Non-parallel gravitational field deep space impact control method and system based on proportional guidance
CN112009698B (en) * 2020-08-28 2021-11-23 北京空间技术研制试验中心 Return flight energy management method for hypersonic cruise aircraft
CN112009698A (en) * 2020-08-28 2020-12-01 北京空间技术研制试验中心 Return flight energy management method for hypersonic cruise aircraft
CN112162567A (en) * 2020-09-09 2021-01-01 北京航空航天大学 Avoidance guidance method suitable for online no-fly zone of aircraft
CN112861249B (en) * 2020-12-21 2023-03-31 中国人民解放军96901部队23分队 Method for calculating degradation solution of gliding trajectory along with speed change based on attack angle and resistance
CN112507466B (en) * 2020-12-21 2023-05-09 中国人民解放军96901部队23分队 Method for calculating glide trajectory energy-dependent degradation solution considering earth rotation rate
CN112861249A (en) * 2020-12-21 2021-05-28 中国人民解放军96901部队23分队 Gliding trajectory degradation solution along with speed change based on attack angle and resistance acceleration
CN112861250A (en) * 2020-12-21 2021-05-28 中国人民解放军96901部队23分队 Gliding trajectory degradation solution along with energy change based on attack angle and resistance acceleration
CN112507467B (en) * 2020-12-21 2023-05-12 中国人民解放军96901部队23分队 Method for calculating descending order solution of glide trajectory along with speed change based on resistance and lift-drag ratio
CN112507465A (en) * 2020-12-21 2021-03-16 中国人民解放军96901部队23分队 Gliding trajectory degradation solution with energy change based on resistance and longitudinal lift-drag ratio
CN112507467A (en) * 2020-12-21 2021-03-16 中国人民解放军96901部队23分队 Gliding trajectory degradation solution with speed change based on resistance and longitudinal lift-drag ratio
CN112507465B (en) * 2020-12-21 2023-04-14 中国人民解放军96901部队23分队 Gliding trajectory energy change degradation solution calculation method based on resistance and lift-drag ratio
CN112507466A (en) * 2020-12-21 2021-03-16 中国人民解放军96901部队23分队 Gliding trajectory energy change-dependent degradation solution considering earth rotation rate
CN112861250B (en) * 2020-12-21 2023-03-28 中国人民解放军96901部队23分队 Method for calculating degradation solution of gliding trajectory along with energy change based on attack angle and resistance
CN113375634A (en) * 2021-04-30 2021-09-10 北京临近空间飞行器系统工程研究所 Altitude measurement method based on atmospheric model and aircraft normal overload combination
CN113741509A (en) * 2021-07-30 2021-12-03 北京航空航天大学 Hypersonic glide aircraft downward pressing section energy management method
CN113741509B (en) * 2021-07-30 2024-02-27 北京航空航天大学 Hypersonic gliding aircraft hold-down section energy management method
CN114265433A (en) * 2021-12-24 2022-04-01 北京航空航天大学 Transverse guidance method and system for meeting large-transverse-stroke maneuver
CN114265433B (en) * 2021-12-24 2024-02-02 北京航空航天大学 Transverse guidance method and system for meeting large-journey maneuver
CN114675545B (en) * 2022-05-26 2022-08-23 中国人民解放军火箭军工程大学 Hypersonic aircraft reentry cooperative guidance method based on reinforcement learning
CN114675545A (en) * 2022-05-26 2022-06-28 中国人民解放军火箭军工程大学 Hypersonic aircraft reentry cooperative guidance method based on reinforcement learning

Also Published As

Publication number Publication date
CN111306989B (en) 2021-12-28

Similar Documents

Publication Publication Date Title
CN111306989B (en) Hypersonic velocity reentry guidance method based on steady glide trajectory analytic solution
CN109508030B (en) Collaborative analysis reentry guidance method considering multi-forbidden flight area constraint
CN107941087B (en) A kind of superb steady gliding reentry guidance method of high lift-drag ratio based on resistance profiles
CN109740198B (en) Analytic prediction-based three-dimensional reentry guidance method for gliding aircraft
Lu et al. Adaptive terminal guidance for hypervelocity impact in specified direction
Xie et al. Highly constrained entry trajectory generation
CN106227972A (en) A kind of optimization method of the steady glide trajectories of hypersonic aircraft
CN106371312B (en) Lift formula based on fuzzy controller reenters prediction-correction method of guidance
Li et al. Improved artificial potential field based lateral entry guidance for waypoints passage and no-fly zones avoidance
CN109062241B (en) Autonomous full-shot reentry guidance method based on linear pseudo-spectrum model predictive control
CN104656450B (en) Method for designing hypersonic aircraft to re-enter into trajectory in smooth gliding mode
CN109270960A (en) Online Optimal Feedback reentry guidance method based on Radau puppet spectrometry
Yu et al. Analytical entry guidance for coordinated flight with multiple no-fly-zone constraints
CN107861517A (en) The online trajectory planning method of guidance of great-jump-forward reentry vehicle based on linear pseudo- spectrum
CN110015446B (en) Semi-analytic Mars entry guidance method
CN106444822A (en) Space vector field guidance based stratospheric airship&#39;s trajectory tracking control method
Leavitt et al. Performance of evolved acceleration guidance logic for entry (EAGLE)
CN112558475A (en) Reentry guidance method considering flight time based on convex optimization
CN113835442B (en) Hypersonic gliding aircraft linear pseudo-spectrum reentry guidance method and hypersonic gliding aircraft linear pseudo-spectrum reentry guidance system
CN113741509B (en) Hypersonic gliding aircraft hold-down section energy management method
CN107102547A (en) A kind of RLV landing phase Guidance Law acquisition methods based on sliding mode control theory
CN106200664A (en) A kind of adapt to attitude control method the most out of control
Yomchinda et al. Autonomous control and path planning for autorotation of unmanned helicopters
CN105354380A (en) Perturbation factor effect-compensated method for rapidly correcting glide trajectory
CN113093789B (en) Planning method for avoiding trajectory of aircraft no-fly zone based on path point optimization

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
GR01 Patent grant
GR01 Patent grant