CN112861249A - Gliding trajectory degradation solution along with speed change based on attack angle and resistance acceleration - Google Patents

Gliding trajectory degradation solution along with speed change based on attack angle and resistance acceleration Download PDF

Info

Publication number
CN112861249A
CN112861249A CN202011532176.3A CN202011532176A CN112861249A CN 112861249 A CN112861249 A CN 112861249A CN 202011532176 A CN202011532176 A CN 202011532176A CN 112861249 A CN112861249 A CN 112861249A
Authority
CN
China
Prior art keywords
gliding
equation
formula
solution
aircraft
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
CN202011532176.3A
Other languages
Chinese (zh)
Other versions
CN112861249B (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.)
23 Units Of Chinese People's Liberation Army 96901 Force
Original Assignee
23 Units Of Chinese People's Liberation Army 96901 Force
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 23 Units Of Chinese People's Liberation Army 96901 Force filed Critical 23 Units Of Chinese People's Liberation Army 96901 Force
Priority to CN202011532176.3A priority Critical patent/CN112861249B/en
Publication of CN112861249A publication Critical patent/CN112861249A/en
Application granted granted Critical
Publication of CN112861249B publication Critical patent/CN112861249B/en
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Classifications

    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F30/00Computer-aided design [CAD]
    • G06F30/10Geometric CAD
    • G06F30/15Vehicle, aircraft or watercraft design
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10Complex mathematical operations
    • G06F17/16Matrix or vector computation, e.g. matrix-matrix or matrix-vector multiplication, matrix factorization
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F30/00Computer-aided design [CAD]
    • G06F30/20Design optimisation, verification or simulation
    • G06F30/28Design optimisation, verification or simulation using fluid dynamics, e.g. using Navier-Stokes equations or computational fluid dynamics [CFD]
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F2113/00Details relating to the application field
    • G06F2113/08Fluids
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F2119/00Details relating to the type or aim of the analysis or the optimisation
    • G06F2119/14Force analysis or force optimisation, e.g. static or dynamic forces
    • YGENERAL 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
    • Y02TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
    • Y02TCLIMATE CHANGE MITIGATION TECHNOLOGIES RELATED TO TRANSPORTATION
    • Y02T90/00Enabling technologies or technologies with a potential or indirect contribution to GHG emissions mitigation

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • Geometry (AREA)
  • Mathematical Analysis (AREA)
  • Mathematical Optimization (AREA)
  • Mathematical Physics (AREA)
  • Pure & Applied Mathematics (AREA)
  • Computational Mathematics (AREA)
  • General Engineering & Computer Science (AREA)
  • Computer Hardware Design (AREA)
  • Evolutionary Computation (AREA)
  • Computing Systems (AREA)
  • Algebra (AREA)
  • Data Mining & Analysis (AREA)
  • Fluid Mechanics (AREA)
  • Automation & Control Theory (AREA)
  • Aviation & Aerospace Engineering (AREA)
  • Databases & Information Systems (AREA)
  • Software Systems (AREA)
  • Control Of Position, Course, Altitude, Or Attitude Of Moving Bodies (AREA)
  • Navigation (AREA)

Abstract

The invention belongs to the technical field of guidance, and relates to a gliding trajectory degradation solution along with speed change based on an attack angle and a resistance acceleration. The gliding trajectory degradation solution along with the change of the speed based on the attack angle and the resistance acceleration comprises the following steps: establishing a dynamic model of the gliding aircraft; establishing an approximate equation of a dynamic model of the gliding aircraft containing three-dimensional motion; and solving the analytic solution of the first approximation equation. Compared with the traditional method, the method has the advantages that the method takes the attack angle and the resistance acceleration as control variables, establishes the three-dimensional gliding aircraft dynamics approximate equation, obviously improves the precision, can quickly solve the large-range and three-dimensional maneuvering trajectory, and can meet the requirements of on-line trajectory planning and quick guidance solving of the missile-borne computer in engineering realization.

Description

Gliding trajectory degradation solution along with speed change based on attack angle and resistance acceleration
Technical Field
The invention belongs to the technical field of guidance, and particularly relates to a gliding trajectory degradation along with speed change based on an attack angle and a resistance acceleration.
Background
The guidance planning of the hypersonic trajectory is complex and the calculation amount is large, the traditional numerical calculation method for solving the differential equation is long in calculation time, and the requirements for online trajectory planning and guidance quick calculation are difficult to meet. Because ballistic planning needs to execute a large number of iterative operations, the analytic solutions have great value in the aspect of online ballistic planning, but most of glide ballistic analytic solutions in the prior art are difficult to meet the engineering application requirements of online ballistic planning and guidance quick solving due to the defects in the aspects of calculation speed and precision when processing three-dimensional maneuvering ballistics with large-scale transverse motion, so that a price-reducing solution which can quickly solve the large-scale three-dimensional maneuvering ballistics and can meet the engineering application requirements in precision is necessary to be researched.
Disclosure of Invention
The invention aims to provide a gliding trajectory degradation solution of an attack angle and a resistance acceleration along with the change of speed, the attack angle and the resistance acceleration are used as control variables, the requirements of high-speed gliding trajectory rapid guidance solution and online trajectory planning are met, a large-range three-dimensional maneuvering trajectory can be rapidly solved compared with a traditional analysis solution, and the calculation speed and the calculation precision can meet the engineering application requirements of online guidance rapid solution.
The technical scheme adopted by the invention is as follows:
1. the gliding trajectory degradation solution along with the change of the speed based on the attack angle and the resistance acceleration is characterized in that:
the method comprises the following steps:
step 1: establishing a dynamics model of a gliding aircraft with a velocity v as an independent variable
Defining: initial longitude, initial latitude and initial azimuth angle of gliding aircraft in geocentric coordinate system are respectively lambda0、φ0And alpha0
Defining: the offset geocentric coordinate system rotates the geocentric coordinate system around the z, y and x axes by a rotating angle lambda respectively0、-φ0And alpha0Obtaining a coordinate system;
defining: the longitude, the latitude and the heading angle of the gliding aircraft under the offset geocentric coordinate system are respectively lambda, phi and psi;
introducing dimensionless heights
Figure BSA0000228278770000011
Wherein
Figure BSA0000228278770000012
re=Re+Hm,HmIs the average flying height R of the gliding aircraft in the gliding flight sectioneIs the radius of the earth, reThe distance between the center of mass of the aircraft and the geocentric; v. the0And vfThe speeds of the aircraft at the glide starting point and the glide ending point respectively;
establishing a dynamic model of the gliding aircraft with the speed v as an independent variable:
Figure BSA0000228278770000021
wherein A isDFor resisting acceleration, LDyIs the longitudinal lift-to-drag ratio, LDzIs the transverse lift-drag ratio, omegaeIs the angular velocity of rotation of the earth, mu is the gravitational constant of the earth;
Figure BSA0000228278770000022
Cσ≈2vωe(sinφ-sinψtanθcosφ)
Cθ≈-2vωecosφcosψ
step 2: approximate equation for establishing dynamics model of gliding aircraft based on Newton iteration method
Let xvWhen v, the variable of formula (1) is replaced and rewritten to y' (x)v)=f(xv,y),
Wherein the content of the first and second substances,
Figure BSA0000228278770000023
remember yiIs y' (x)v)=f(xvY) and setting the i +1 th iteration value as yi+1=yi+δyiTo, for
Figure BSA0000228278770000031
A first order Taylor approximation is performed on theta, phi, psi, lambda to expand the function y' (x)v)=f(xvY) is constructed as a newton iteration in function space, i.e.:
Figure BSA0000228278770000032
wherein
Figure BSA0000228278770000033
xvini
Figure BSA0000228278770000034
θini、φini、ψini、λiniInitial conditions for compatibility;
knowing the value y of the ith iterationiThen, the solution y of the (i + 1) th step can be obtained by using the formula (2)i+1
Defining an equation formed by using the formula (2) at the Nth time as an approximate equation of the formula (1) at the Nth time, and solving the approximate equation into an approximate solution of the Nth time; definition of
Figure BSA0000228278770000035
θ0,φ0,ψ0,λ0Is a zero-order approximate solution of formula (1);
Figure BSA0000228278770000036
θ1,φ1,ψ1,λ1is a first order approximation solution of formula (1);
when the offset geocentric system equation is adopted, the gliding aircraft approximately flies to the target along the equator of the offset geocentric coordinate system, and the gliding flight section has phi approximately equal to 0, theta approximately equal to 0, psi approximately equal to 0,
Figure BSA0000228278770000037
The characteristic of (1) is that the zero-order approximate solution is:
Figure BSA0000228278770000038
substituting the equation (3) into the equation (2) to obtain a first order approximation equation of the glide flight segment as follows:
Figure BSA0000228278770000041
wherein
Figure BSA0000228278770000042
Figure BSA0000228278770000043
As can be seen from the equation (4), the Jacobian matrix has the characteristic of block sparsity, and theta can be solved first1And
Figure BSA0000228278770000044
and phi1And psi1Then solve for λ1
Due to phi1And psi1Is the solution of the lateral motion variable,
Figure BSA0000228278770000045
and theta1Is a solution of a longitudinal motion variable, phi1And psi1To pair
Figure BSA0000228278770000046
And theta1Insensitive, equation (4) is further simplified as:
Figure BSA0000228278770000047
step 3: establishing the relation between the three-dimensional section and the motion state
Using angle of attack and resistance acceleration as control variables, order
Figure BSA0000228278770000048
Wherein, c0,c1,c2Are all controlled variable resistance acceleration ADThe design parameters of (1); α is the angle of attack, and f (v) is the angle of attack control function designed with velocity as the independent variable.
Defining: the flight height of the gliding aircraft is h, and the first derivative of the aerodynamic drag acceleration and the aerodynamic drag speed with respect to the energy is obtained by:
Figure BSA0000228278770000051
finishing to obtain:
Figure BSA0000228278770000052
wherein: cDIs the drag coefficient of the gliding aircraft, hsThe coefficients of the exponential formula being a function of atmospheric density;
neglecting the effect of the drag coefficient derivative, the approximate altitude of the gliding aircraft can be found:
Figure BSA0000228278770000053
and (3) solving a second derivative related to energy for the resistance acceleration to obtain a longitudinal lift-drag ratio:
LDy=a(AD″-b)
wherein the content of the first and second substances,
Figure BSA0000228278770000054
and (3) utilizing a pneumatic coefficient table or a pneumatic interpolation function to obtain a transverse lift-drag ratio:
Figure BSA0000228278770000055
wherein:
Figure BSA0000228278770000056
the sign of which is determined by judging the azimuth deviation;
step 4: analytic solution to solve first order approximation equation
Substituting formula (6) for formula (3) according to the definition of the offset geocentric coordinate system, and calculating from v0By v, we get:
Figure BSA0000228278770000057
the first two formulas of the combined vertical type (5) can obtain
Figure BSA0000228278770000061
Order to
Figure BSA0000228278770000062
Left-multiplying equation (7) by M (v, v)0) After arrangement, the differential multiplication rule is applied reversely, and then integration and arrangement are carried out to obtain:
Figure BSA0000228278770000063
using M (v, v)0) Solving its inverse matrix [ M (v, v)0)]-1Formula (8) is substituted for
Figure BSA0000228278770000064
Finish to obtain phi1And psi1The analytical formula (2):
Figure BSA0000228278770000065
lagrange interpolation polynomial in interval [ v ] by taking root of Legendre of degree n as sampling nodef,v0]Approximate λ (v, v) on segment0)、cos(λ(xv,v0))m3(xv)、sin(λ(xv,v0))m3(xv) Obtaining:
Figure BSA0000228278770000066
Figure BSA0000228278770000067
Figure BSA0000228278770000068
substituting formula (9) to obtain phi1And psi1Analytic solution of (2):
Figure BSA0000228278770000071
wherein: c. Cp0(k)、cp1(k)、cp2(k)Is the coefficient of the interpolation polynomial;
by calculated phi1And psi1Substituting formula (5) to obtain λ1Analytic solution of (2):
Figure BSA0000228278770000072
the invention has the beneficial effects that:
the invention establishes a gliding trajectory degradation solution of the attack angle and the resistance acceleration along with the change of the speed, realizes the rapid calculation of the high-altitude gliding three-dimensional trajectory, has good calculation speed and precision in the aspect of calculating the large-range three-dimensional maneuvering trajectory compared with the traditional analytic solution, and can meet the application requirement of the online guidance planning engineering of the large-range three-dimensional gliding aircraft.
Drawings
None.
Detailed Description
The technical scheme of the invention is further specifically described by combining specific embodiments.
1. The gliding trajectory degradation solution along with the change of the speed based on the attack angle and the resistance acceleration is characterized in that:
the method comprises the following steps:
step 1: establishing a dynamics model of a gliding aircraft with a velocity v as an independent variable
Defining: initial longitude, initial latitude and initial azimuth angle of gliding aircraft in geocentric coordinate system are respectively lambda0、φ0And alpha0
Defining: the offset geocentric coordinate system rotates the geocentric coordinate system around the z, y and x axes by a rotating angle lambda respectively0、-φ0And alpha0Obtaining a coordinate system;
defining: the longitude, the latitude and the heading angle of the gliding aircraft under the offset geocentric coordinate system are respectively lambda, phi and psi;
taking into account that the flying height during gliding flight is a small quantity relative to the radius of the earth, dimensionless heights are introduced
Figure BSA0000228278770000073
Wherein
Figure BSA0000228278770000074
re=Re+Hm,HmFor the average flight height, R, of the gliding aircraft in the gliding flight sectioneIs the radius of the earth, reThe distance between the center of mass of the aircraft and the geocentric; v. the0And vfThe speeds of the aircraft at the glide starting point and the glide ending point respectively;
the gliding flight section is that after the gliding aircraft enters again, the gliding aircraft has a small ballistic inclination angle and a long distanceA flight segment away from the landing area toward the land. Establishing a dynamics model of the gliding aircraft with the velocity v as an independent variable (introducing dimensionless altitude into a classical model of the hypersonic dynamics)
Figure BSA0000228278770000081
To obtain):
Figure BSA0000228278770000082
wherein A isDFor resisting acceleration, LDyIs the longitudinal lift-to-drag ratio, LDzIs the transverse lift-drag ratio, omegaeIs the angular velocity of rotation of the earth, mu is the gravitational constant of the earth;
Figure BSA0000228278770000083
Cσ≈2vωe(sinφ-sinψtanθcosφ)
Cθ≈-2vωecosφcosψ
step 2: approximate equation for establishing dynamics model of gliding aircraft based on Newton iteration method
Let xvWhen v, the variable of formula (1) is replaced and rewritten to y' (x)v)=f(xv,y),
Wherein the content of the first and second substances,
Figure BSA0000228278770000091
and (3) a Newton iteration method is applied to approximate the multi-dimensional nonlinear trajectory equation into a plurality of low-dimensional variable coefficient linear differential equations so as to realize high-precision rapid calculation of the flight trajectory. Remember yiIs y' (x)v)=f(xvY) and setting the i +1 th iteration value as yi+1=yi+δyi. To solve equation (2), a zero-order approximation solution of equation (1) needs to be provided before starting the iterative computation. The closer the zero-order solution is to the true solution, the number of iterations required to achieve the same accuracy requirementThe less. For the problem of aircraft guidance, a larger ballistic estimation error is usually allowed, if a zero-time solution is closer to a true solution, the requirement can be met by one-time solution, and the method is suitable for the problem of aircraft guidance
Figure BSA0000228278770000092
A first order Taylor approximation is performed on theta, phi, psi, lambda to expand the function y' (x)v)=f(xvY) is constructed as a newton iteration in function space, i.e.:
Figure BSA0000228278770000093
wherein
Figure BSA0000228278770000094
xvini
Figure BSA0000228278770000095
θini、φini、ψini、λiniInitial conditions for compatibility;
knowing the value y of the ith iterationiThen, the solution y of the (i + 1) th step can be obtained by using the formula (2)i+1
Defining an equation formed by using the formula (2) at the Nth time as an approximate equation of the formula (1) at the Nth time, and solving the approximate equation into an approximate solution of the Nth time; definition of
Figure BSA0000228278770000096
θ0,φ0,ψ0,λ0Is a zero-order approximate solution of formula (1);
Figure BSA0000228278770000097
θ1,φ1,ψ1,λ1is a first order approximation solution of formula (1);
c used for calculating the Jacobian matrix in order to obtain the Jacobian matrix which is accurate and has the blocking sparse characteristic of decoupling the equationθ、CσIn the expression, only ω is ignorede 2And J2, trying to retain ωeThe first order item of (2). When the offset geocentric system equation is adopted, the gliding aircraft approximately flies to the target along the equator of the offset geocentric coordinate system, and the gliding flight section has phi ≈ 0, theta ≈ 0, psi ≈ 0,
Figure BSA0000228278770000101
The characteristic of (1) is that the zero-order approximate solution is:
Figure BSA0000228278770000102
substituting the equation (3) into the equation (2) to obtain a first order approximation equation of the glide flight segment as follows:
Figure BSA0000228278770000103
wherein
Figure BSA0000228278770000104
Figure BSA0000228278770000105
Equation (4) retains the basic characteristics of equation (1) and contains a first term of the earth rotation rate, so that the accuracy of the 1-time approximate solution is obviously improved compared with the accuracy of the 0-time approximate solution, and the 1-time approximate solution is used for guidance calculation to obtain better guidance performance. If desired, equation (2) can be used multiple times to obtain a more accurate solution.
As can be seen from the equation (4), the Jacobian matrix has the characteristic of block sparsity, and theta can be solved first1And
Figure BSA0000228278770000106
and phi1And psi1Then solve for λ1
Due to theta1And psi1Is in the transverse directionThe solution of the motion variable is then determined,
Figure BSA0000228278770000107
and theta1Is a solution of a longitudinal motion variable, phi1And psi1To pair
Figure BSA0000228278770000108
And theta1Insensitive, equation (4) is further simplified as:
Figure BSA0000228278770000111
step 3: establishing the relation between the three-dimensional section and the motion state
On the premise of meeting the requirements of trajectory planning, the angle of attack and the resistance acceleration are taken as control variables to order
Figure BSA0000228278770000112
Wherein, c0,c1,c2Are all controlled variable resistance acceleration ADThe design parameters of (1); α is the angle of attack, and f (v) is the angle of attack control function designed with velocity as the independent variable.
Defining: the flight height of the gliding aircraft is h, and the first derivative of the aerodynamic drag acceleration and the aerodynamic drag speed with respect to the energy is obtained by:
Figure BSA0000228278770000113
finishing to obtain:
Figure BSA0000228278770000114
wherein: cDIs the drag coefficient of the gliding aircraft, hsThe coefficients of the exponential formula being a function of atmospheric density;
neglecting the effect of the drag coefficient derivative, the approximate altitude of the gliding aircraft can be found:
Figure BSA0000228278770000115
and (3) solving a second derivative related to energy for the resistance acceleration to obtain a longitudinal lift-drag ratio:
LDy=a(AD″-b)
wherein the content of the first and second substances,
Figure BSA0000228278770000121
and (3) utilizing a pneumatic coefficient table or a pneumatic interpolation function to obtain a transverse lift-drag ratio:
Figure BSA0000228278770000122
wherein:
Figure BSA0000228278770000123
the sign of which is determined by judging the azimuth deviation;
step 4: analytic solution to solve first order approximation equation
Since the zeroth order generalized longitude equation is independent, first derivation is made here about λ0Is defined by the offset geocentric coordinate system, formula (6) is substituted for formula (3), and the formula (6) is derived from v0By v, we get:
Figure BSA0000228278770000124
to solve for phi1,ψ1The first two formulas of the combined vertical type (5) can be obtained
Figure BSA0000228278770000125
Order to
Figure BSA0000228278770000126
Left-multiplying equation (7) by M (v, v)0) After arrangement, the differential multiplication rule is applied reversely, and then integration and arrangement are carried out to obtain:
Figure BSA0000228278770000127
using M (v, v)0) Solving its inverse matrix [ M (v, v)0)]-1Formula (8) is substituted for
Figure BSA0000228278770000128
Finish to obtain phi1And psi1The analytical formula (2):
Figure BSA0000228278770000131
lagrange interpolation polynomial in interval [ v ] by taking root of Legendre of degree n as sampling nodef,v0]Approximate λ (v, v) on segment0)、cos(λ(xv,v0))m3(xv)、sin(λ(xv,v0))m3(xv) Obtaining:
Figure BSA0000228278770000132
Figure BSA0000228278770000133
Figure BSA0000228278770000134
substituting formula (9) to obtain phi1And psi1Analytic solution of (2):
Figure BSA0000228278770000135
wherein: c. Cp0(k)、cp1(k)、cp2(k)Is the coefficient of the interpolation polynomial;
by calculated phi1And psi1Substituting formula (5) to obtain λ1Analytic solution of (2):
Figure BSA0000228278770000136
the above detailed description of the present invention is only used for illustrating the present invention and is not limited to the technical solutions described in the embodiments of the present invention, and it should be understood by those skilled in the art that the present invention can be modified or substituted equally to achieve the same technical effects; as long as the use requirements are met, the method is within the protection scope of the invention.

Claims (1)

1. The gliding trajectory degradation solution along with the change of the speed based on the attack angle and the resistance acceleration is characterized in that:
the method comprises the following steps:
step 1: establishing a dynamics model of a gliding aircraft with a velocity v as an independent variable
Defining: initial longitude, initial latitude and initial azimuth angle of gliding aircraft in geocentric coordinate system are respectively lambda0、φ0And alpha0
Defining: the offset geocentric coordinate system rotates the geocentric coordinate system around the z, y and x axes by a rotating angle lambda respectively0、-φ0And alpha0Obtaining a coordinate system;
defining: the longitude, the latitude and the heading angle of the gliding aircraft under the offset geocentric coordinate system are respectively lambda, phi and psi;
introducing dimensionless heights
Figure FSA0000228278760000011
Wherein
Figure FSA0000228278760000012
re=Re+Hm,HmFor the average flight height, R, of the gliding aircraft in the gliding flight sectioneIs the radius of the earth, reThe distance between the center of mass of the aircraft and the geocentric; v. the0And vfThe speeds of the aircraft at the glide starting point and the glide ending point respectively;
establishing a dynamic model of the gliding aircraft with the speed v as an independent variable:
Figure FSA0000228278760000013
wherein A isDFor resisting acceleration, LDyIs the longitudinal lift-to-drag ratio, LDzIs the transverse lift-drag ratio, omegaeIs the angular velocity of rotation of the earth, mu is the gravitational constant of the earth;
Figure FSA0000228278760000014
Cσ≈2vωe(sinφ-sinψtanθcosφ)
Cθ≈-2vωecosφcosψ
step 2: approximate equation for establishing dynamics model of gliding aircraft based on Newton iteration method
Let xvWhen v, the variable of formula (1) is replaced and rewritten to y' (x)v)=f(xv,y),
Wherein the content of the first and second substances,
Figure FSA0000228278760000021
remember yiIs y' (x)v)=f(xvY) and setting the i +1 th iteration value as yi+1=yi+δyiTo, for
Figure FSA0000228278760000022
A first order Taylor approximation is performed on theta, phi, psi, lambda to expand the function y' (x)v)=f(xvY) is constructed as a newton iteration in function space, i.e.:
Figure FSA0000228278760000023
yi+1(vini)=yini
wherein
Figure FSA0000228278760000024
xvini
Figure FSA0000228278760000025
θini、θini、ψini、λiniInitial conditions for compatibility;
knowing the value y of the ith iterationiThen, the solution y of the (i + 1) th step can be obtained by using the formula (2)i+1
Defining an equation formed by using the formula (2) at the Nth time as an approximate equation of the formula (1) at the Nth time, and solving the approximate equation into an approximate solution of the Nth time; definition of
Figure FSA0000228278760000026
θ0,φ0,ψ0,λ0Is a zero-order approximate solution of formula (1);
Figure FSA0000228278760000027
θ1,φ1,ψ1,λ1is a first order approximation solution of formula (1);
when the offset geocentric system equation is adopted, the gliding aircraft approximately flies to the target along the equator of the offset geocentric coordinate system, and the gliding flight section has phi approximately equal to 0, theta approximately equal to 0, psi approximately equal to 0,
Figure FSA0000228278760000028
The characteristic of (1) is that the zero-order approximate solution is:
Figure FSA0000228278760000031
substituting the equation (3) into the equation (2) to obtain a first order approximation equation of the glide flight segment as follows:
Figure FSA0000228278760000032
wherein
Figure FSA0000228278760000033
Figure FSA0000228278760000034
As can be seen from the equation (4), the Jacobian matrix has the characteristic of block sparsity, and theta can be solved first1And
Figure FSA0000228278760000035
and phi1And psi1Then solve for λ1
Due to phi1And psi1Is the solution of the lateral motion variable,
Figure FSA0000228278760000036
and theta1Is a solution of a longitudinal motion variable, phi1And psi1To pair
Figure FSA0000228278760000037
And theta1Insensitive, equation (4) is further simplified as:
Figure FSA0000228278760000038
step 3: establishing the relation between the three-dimensional section and the motion state
Using angle of attack and resistance acceleration as control variables, order
Figure FSA0000228278760000041
α=f(v)
Wherein, c0,c1,c2Are all controlled variable resistance acceleration ADThe design parameters of (1); α is the angle of attack, and f (v) is the angle of attack control function designed with velocity as the independent variable.
Defining: the flight height of the gliding aircraft is h, and the first derivative of the aerodynamic drag acceleration and the aerodynamic drag speed with respect to the energy is obtained by:
Figure FSA0000228278760000042
finishing to obtain:
Figure FSA0000228278760000043
wherein: cDIs the drag coefficient of the gliding aircraft, hsThe coefficients of the exponential formula being a function of atmospheric density;
neglecting the effect of the drag coefficient derivative, the approximate altitude of the gliding aircraft can be found:
Figure FSA0000228278760000044
and (3) solving a second derivative related to energy for the resistance acceleration to obtain a longitudinal lift-drag ratio:
LDy=a(AD″-b)
wherein the content of the first and second substances,
Figure FSA0000228278760000045
and (3) utilizing a pneumatic coefficient table or a pneumatic interpolation function to obtain a transverse lift-drag ratio:
Figure FSA0000228278760000046
wherein:
Figure FSA0000228278760000051
the sign of which is determined by judging the azimuth deviation;
step 4: analytic solution to solve first order approximation equation
Substituting formula (6) for formula (3) according to the definition of the offset geocentric coordinate system, and calculating from v0By v, we get:
Figure FSA0000228278760000052
the first two formulas of the combined vertical type (5) can obtain
Figure FSA0000228278760000053
Order to
Figure FSA0000228278760000054
Left-multiplying equation (7) by M (v, v)0) After arrangement, the differential multiplication rule is applied reversely, and then integration and arrangement are carried out to obtain:
Figure FSA0000228278760000055
using M (v, v)0) Solving its inverse matrix [ M (v, v)0)]-1Formula (8) is substituted for
Figure FSA0000228278760000056
Finish to obtain phi1And psi1The analytical formula (2):
Figure FSA0000228278760000057
lagrange interpolation polynomial in interval [ v ] by taking root of Legendre of degree n as sampling nodef,v0]Approximate λ (v, v) on segment0)、cos(λ(xv,v0))m3(xv)、sin(λ(xv,v0))m3(xv) Obtaining:
Figure FSA0000228278760000058
Figure FSA0000228278760000059
Figure FSA0000228278760000061
substituting formula (9) to obtain phi1And psi1Analytic solution of (2):
Figure FSA0000228278760000062
wherein: c. Cp0(k)、cp1(k)、cp2(k)Is the coefficient of the interpolation polynomial;
by calculated phi1And psi1Substituting formula (5) to obtain λ1Analytic solution of (2):
Figure FSA0000228278760000063
CN202011532176.3A 2020-12-21 2020-12-21 Method for calculating degradation solution of gliding trajectory along with speed change based on attack angle and resistance Active CN112861249B (en)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202011532176.3A CN112861249B (en) 2020-12-21 2020-12-21 Method for calculating degradation solution of gliding trajectory along with speed change based on attack angle and resistance

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202011532176.3A CN112861249B (en) 2020-12-21 2020-12-21 Method for calculating degradation solution of gliding trajectory along with speed change based on attack angle and resistance

Publications (2)

Publication Number Publication Date
CN112861249A true CN112861249A (en) 2021-05-28
CN112861249B CN112861249B (en) 2023-03-31

Family

ID=75996263

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202011532176.3A Active CN112861249B (en) 2020-12-21 2020-12-21 Method for calculating degradation solution of gliding trajectory along with speed change based on attack angle and resistance

Country Status (1)

Country Link
CN (1) CN112861249B (en)

Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN107941087A (en) * 2017-10-18 2018-04-20 北京航空航天大学 A kind of superb steady gliding reentry guidance method of high lift-drag ratio based on resistance profiles
CN111306989A (en) * 2020-03-12 2020-06-19 北京航空航天大学 Hypersonic velocity reentry guidance method based on steady glide trajectory analytic solution
CN111488646A (en) * 2020-03-02 2020-08-04 北京航空航天大学 Analytic solving method for hypersonic steady gliding trajectory under rotating earth

Patent Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN107941087A (en) * 2017-10-18 2018-04-20 北京航空航天大学 A kind of superb steady gliding reentry guidance method of high lift-drag ratio based on resistance profiles
CN111488646A (en) * 2020-03-02 2020-08-04 北京航空航天大学 Analytic solving method for hypersonic steady gliding trajectory under rotating earth
CN111306989A (en) * 2020-03-12 2020-06-19 北京航空航天大学 Hypersonic velocity reentry guidance method based on steady glide trajectory analytic solution

Non-Patent Citations (3)

* Cited by examiner, † Cited by third party
Title
侯明善: "防区外投放制导炸弹滑翔段垂直面最优制导", 《兵工学报》 *
孙建波等: "滑翔飞行器线性伪谱模型预测控制三维轨迹规划", 《红外与激光工程》 *
王洁瑶等: "助推-滑翔弹道高精度滑翔射程解析估算方法", 《宇航学报》 *

Also Published As

Publication number Publication date
CN112861249B (en) 2023-03-31

Similar Documents

Publication Publication Date Title
CN111783307B (en) Hypersonic aircraft state estimation method
CN109116860B (en) Nonlinear robust control method for three-rotor unmanned aerial vehicle
CN105184109B (en) Disturb trajectory motors in boost phase penetration state deviation analytic method under graviational interaction
CN109614633A (en) A kind of composite rotor craft non-linear modeling method and Calculate Ways
CN107330152B (en) Efficient pneumatic balancing method suitable for rotor craft
Zhu et al. Impact time and angle control guidance independent of time-to-go prediction
CN109885074B (en) Finite time convergence attitude control method for quad-rotor unmanned aerial vehicle
CN111290278B (en) Hypersonic aircraft robust attitude control method based on prediction sliding mode
CN108681331A (en) A kind of Attitude tracking control method of Near Space Flying Vehicles
CN104536448B (en) Backstepping based control method for unmanned-plane attitude system
CN112198885A (en) Unmanned aerial vehicle control method capable of meeting autonomous landing requirement of maneuvering platform
CN111008488B (en) Propeller unmanned aerial vehicle launching process reaction torque modeling method
CN112861250B (en) Method for calculating degradation solution of gliding trajectory along with energy change based on attack angle and resistance
CN116187199A (en) Unsteady aerodynamic modeling method integrating intelligent technology
CN112507467B (en) Method for calculating descending order solution of glide trajectory along with speed change based on resistance and lift-drag ratio
CN116627156B (en) Four-rotor unmanned aerial vehicle attitude disturbance rejection control method
CN112861249B (en) Method for calculating degradation solution of gliding trajectory along with speed change based on attack angle and resistance
Meng et al. Characteristic model based control of the X-34 reusable launch vehicle in its climbing phase
CN112507465B (en) Gliding trajectory energy change degradation solution calculation method based on resistance and lift-drag ratio
CN110989397A (en) Aircraft accident search simulation method and system
CN112861251B (en) Glide trajectory degradation solution calculation method with speed change considered by earth rotation rate
CN112507466A (en) Gliding trajectory energy change-dependent degradation solution considering earth rotation rate
CN110562492B (en) Method for quickly generating Mars atmospheric entrance track of detector
Abdallah et al. Modelling and simulation of an anti-tank guided missile
Kim et al. Flight dynamics analyses of a propeller-driven airplane (I): aerodynamic and inertial modeling of the propeller

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