CN116502570B - Analysis method for longitudinal and transverse coupling motion stability of ultra-high-speed aircraft - Google Patents

Abstract

The application belongs to the field of aerodynamics and flight mechanics, and discloses a method for analyzing longitudinal and transverse coupling motion stability of an ultra-high-speed aircraft, which comprises the following steps: aiming at the motion characteristics of the aircraft to be researched, a pitching channel motion equation is established; obtaining a rolling moment static derivative and a rolling moment dynamic derivative through CFD calculation; obtaining a functional expression of a rolling moment static derivative relative to a pitch angle through polynomial fitting; establishing a rolling channel dynamics equation considering inertial coupling and pneumatic/motion coupling effects; solving the stability limit of the equation by a Hill infinite determinant; solving a critical frequency of pitching motion which causes the rolling motion to diverge through a stability limit; the above steps are repeated by selecting different motion amplitudes and average pitch angles. The application solves the problem of low efficiency of the existing ultra-high speed aircraft in multi-channel coupling instability analysis, and is suitable for motion stability analysis considering inertial coupling and pneumatic/motion coupling effects.

Description

Analysis method for longitudinal and transverse coupling motion stability of ultra-high-speed aircraft

Technical Field

The application relates to the field of aerodynamics and flight mechanics, in particular to a method for analyzing the stability of longitudinal and transverse coupling movement of an ultra-high-speed aircraft.

Background

The ultra-high speed aircraft is affected by the engine body configuration and the high-altitude high-speed flight environment, the lateral stability is poor, instability caused by multi-channel coupling motion is easy to occur in the maneuvering process, for example, the multi-channel coupling instability of the American HTV-2 aircraft in 2010 causes test flight failure. When the ultra-high speed aircraft performs pitching rapid maneuver, stronger inertial coupling and pneumatic/motion coupling effects exist between the pitching channel and the rolling channel, so that pitching motion interferes with the rolling channel, thereby affecting the stability of the rolling channel and possibly leading to instability of the aircraft. But currently there is less research on analysis of the stability of the coupled motion of ultra-high speed aircraft in the longitudinal/lateral (i.e. pitch/roll) direction.

Conventional stability analysis methods generally divide the motion of an aircraft into a baseline motion and a perturbed motion based on small perturbation assumptions. The disturbance motion of the aircraft is described by separate models of longitudinal (pitch) and lateral heading (roll and yaw), respectively, a set of differential equations describing the motion of the aircraft is linearized, and the linear stability theory is used as a theoretical analysis method. However, in the process that the pitch angle of the ultra-high-speed aircraft is widely changed, the pitch channel and the roll channel have strong inertial coupling and pneumatic/motion coupling effects, at the moment, the motion equation is nonlinear, and the traditional theory method is not applicable any more, so that a new stability analysis method is required to be sought.

In addition to theoretical analysis, aircraft stability predictions can be developed, typically based on CFD calculations and wind tunnel test methods. However, aiming at the multi-channel coupling motion working condition of the ultra-high speed aircraft, a great amount of calculation resources are required to be consumed for directly developing CFD calculation, the calculation period is long, and the problems of complex technology and high cost exist in performing the wind tunnel multi-degree-of-freedom dynamic test.

In order to meet engineering application requirements of ultra-high speed aircraft technology, development of an analysis method for longitudinal/transverse coupling motion stability of an ultra-high speed aircraft is needed currently so as to realize rapid and efficient stability prediction.

Disclosure of Invention

The application aims to provide a method for analyzing the stability of longitudinal and transverse coupling movement of an ultra-high speed aircraft, which solves the problem of low efficiency of the existing ultra-high speed aircraft in multichannel coupling instability analysis.

In order to achieve the above object, the present application provides a technical solution: a method for analyzing the stability of longitudinal and transverse coupling movement of a super-high-speed aircraft comprises the following steps:

s1, aiming at the motion characteristics of an aircraft to be researched, given an incoming flow working condition and a pitch angle motion range of the aircraft, selecting a plurality of sample points { theta ] in the pitch angle motion range _i Then selecting the sinusoidal oscillation motion amplitude A and the average pitch angle theta of the pitch channel ₀ Thereby establishing a pitch channel motion equation θ (t);

s2, obtaining a sample point { theta ] through CFD calculation _i Static derivative of roll moment at }θ ₀ Rolling moment derivative at +.>

S3, obtaining a functional expression of the static derivative of the rolling moment relative to the pitch angle through polynomial fittingSubstituting θ (t) into +.>Obtaining the functional expression->

S4, A and theta are calculated ₀ 、θ(t)、Substituting roll channel dynamics equation f taking inertial coupling and pneumatic/kinematic coupling effects into account ₁ And will equation f ₁ Equation f reduced to Mathieu equation form ₂ ；

S5, solving an equation f through a Hill infinite determinant method ₂ Stability margin f (lambda, epsilon) |of (2) _λ＝1 ＝0、f(λ,ε)| _λ＝0.25 ＝0；

S6, through f (lambda, epsilon) | _λ＝1 ＝0、f(λ,ε)| _λ＝0.25 Solution of critical frequency f=0 _c | _λ＝0 、 When the forced pitching frequency f meets f>f _c | _λ＝0 、/> In any one of the above, the rolling motion diverges under the combined action of inertial coupling and pneumatic/kinematic coupling;

s7, selecting different motion amplitudes A and average pitch angles theta ₀ And repeating the steps S1-S6 to obtain the longitudinal/transverse coupling motion stability of the aircraft in the given motion parameter range.

Further, the pitch channel motion equation θ (t) in step S1 is as follows:

θ(t)＝θ ₀ +A sin(ωt) (1)

wherein θ is the pitch angle of the aircraft, θ ₀ The average pitch angle, A is the amplitude, ω is the pitch oscillation angular frequency, and t is the time.

Further, the rolling moment static derivative in step S2Roll moment derivative +.>The calculation method of (2) is as follows:

wherein q is the incoming flow pressure, S _r For aircraft reference area, l _r For the reference length of the aircraft,is the pitch angle theta _i Rolling static derivative at ∈ ->Is the average pitch angle theta ₀ Roll derivative at I _X For winding the aircraft body axis o-x _b Is a rotational inertia of the bearing.

Further, the polynomial fitting method in step S3 uses the following fitting function:

wherein b ₁ 、b ₂ Is the fitting coefficient.

Further, the roll channel dynamics equation f in step S4 ₁ The following are provided:

wherein the method comprises the steps ofA pneumatic/kinematic coupling moment term;

is an inertial coupling moment term;

is a rolling damping moment term;

I _y 、I _z respectively, about the aircraft body axis o-y _b 、o-z _b Is a rotational inertia of the bearing.

Further, in step S4, the roll channel dynamics equation f ₂ The following are provided:

wherein the method comprises the steps of

Further, in step S5, equation f is solved by Hill infinite determinant ₂ Stability margin f (lambda, epsilon) |of (2) _λ＝1 ＝0、f(λ,ε)| _λ＝0.25 The calculation method for=0 is as follows:

for stability limits around λ=1, the solution of the equation with a period of 2π is expressed as follows:

s in _n For each coefficient of the complex fourier series, formula (8) is substituted into formula (6) and is considered Then there are:

formula (9) can be simplified as

At e only ^inτ All coefficients of (2) are zero, i.e.:

to make the constitution of { S } _n The infinite system of equations of } has a non-zero solution, the value of the corresponding coefficient matrix determinant has to be zero, where λ+.n for any n is satisfied ² ；

The form of the coefficient matrix determinant is:

wherein:

similarly, for a stability limit around λ=0.25, the solution of the equation with a period of 4pi is expressed as follows:

the substitution of formula (14) into formula (6) is as follows:

formula (15) can be simplified as

Only atAll coefficients of (2) are zero, equation (9) can be true, namely:

the form of the coefficient matrix determinant is:

wherein:

equation (12) and equation (18) are expressed as:

f(λ,ε,μ,δ)＝0 (20)

taking δ and μ as constants, equation (20) is a hidden function for λ, ε, and the curve f (λ, ε) of ε with λ around λ=1 and λ=0.25 is obtained by the "fimplitit" command of Matlab _λ＝1 ＝0、f(λ,ε)| _λ＝0.25 =0, i.e. stability margin.

Further, the step S6 is performed by f (λ, ε) | _λ＝1 ＝0、f(λ,ε)| _λ＝0.25 Solution of critical frequency f=0 _c | _λ＝0 、The calculation method of (2) is as follows:

for a given moment of inertia parameter a _x Pitch motion parameter θ ₀ A to establish a function epsilon=epsilon (λ) of epsilon with respect to λ, i.e.:

drawing epsilon=epsilon (lambda) in a stability map as a working condition curve, wherein the working condition curve sequentially passes through stability limits near lambda=0, lambda=0.25 and lambda=1 to form 5 intersection points, and the coordinates of the intersection points are set as (lambda) _i ,ε _i ) From formula (7):

thereby obtaining the critical frequency f corresponding to the stability critical point _i Critical frequency f _i Sequentially defined as f in the order of the corresponding lambda value from smaller to larger _c | _λ＝0 、

Compared with the prior art, the application has the beneficial effects that:

aiming at the problem of analysis of longitudinal/transverse coupling motion stability of the ultra-high speed aircraft, the scheme provides a method for constructing a kinetic equation comprising inertial coupling and pneumatic/motion coupling effects and carrying out stability analysis based on Mathieu equation stability theory, so that calculation time and cost can be reduced while prediction accuracy is ensured, qualitative, quantitative, rapid and efficient prediction of the longitudinal/transverse coupling motion stability of the ultra-high speed aircraft is realized, and guidance is provided for relevant theoretical research and engineering application of the ultra-high speed aircraft.

Drawings

FIG. 1 is a flow chart of the present embodiment;

FIG. 2 is a schematic illustration of an ultra-high speed aircraft in this embodiment;

FIG. 3 is a graph showing the comparison of the fitted function curve and the CFD result of the ultra-high speed aircraft in the present embodiment;

FIG. 4 is a graph showing the stability of the ultra-high speed aircraft in the present embodiment;

FIG. 5 is a graph comparing the stability of the working condition f-A of the ultra-high speed aircraft with the CFD result.

Detailed Description

The application is described in further detail below with reference to the attached drawings and embodiments:

examples

As shown in FIG. 1, in the analysis method of the longitudinal and transverse coupling motion stability of the ultra-high-speed aircraft, a kinetic equation comprising inertial coupling and pneumatic/motion coupling effects is firstly required to be constructed, and then stability is carried out based on Mathieu equation stability theory. The embodiment specifically comprises the following steps:

s1, the embodiment takes the ultra-high-speed aircraft shown in fig. 2 as an object to study the longitudinal/transverse coupling motion stability of the ultra-high-speed aircraft. Given that the incoming flow working condition is Mach number 6, the flying height is 28 km, the pitch angle movement range is-5-15 degrees, and a sample point { theta ] is selected in the pitch angle movement range _i = { -5 °,0 °,5 °,10 °,15 ° }. Given pitch channel sinusoidal oscillation motion amplitude a=10°, average pitch angle θ ₀ =5°, thereby establishing a pitch channel motion equation θ (t). The pitch channel motion equation θ (t) is as follows:

θ(t)＝θ ₀ +A sinωt (1)

wherein θ is the pitch angle of the aircraft, θ ₀ The average pitch angle, A is the amplitude, ω is the pitch oscillation angular frequency, and t is the time.

S2, obtaining a sample point { theta ] through CFD calculation _i Static derivative of roll moment at }θ ₀ Rolling moment derivative at +.>Wherein the static derivative of the roll moment +.>Roll moment derivative +.>The calculation method of (2) is as follows:

wherein q is the incoming flow pressure, S _r For aircraft reference area, l _r For the reference length of the aircraft,is the pitch angle theta _i Rolling static derivative at ∈ ->Is the average pitch angle theta ₀ Roll derivative at I _X For winding the aircraft body axis o-x _b The moment of inertia of the body axis o-x _b The definition is shown in fig. 2.

S3, obtaining a functional expression of the static derivative of the rolling moment relative to the pitch angle through polynomial fittingFitting by least squares to obtain a fitted function curve and CFD result pair such as shown in FIG. 3, and obtaining a fitting coefficient b ₁ ＝-2.9951,b ₂ = -67.93698, determining coefficient R ² = 0.998817, the use of a second order polynomial gives a better fit. Substituting θ (t) into +.>Obtaining the functional expression->The fitting function used in the polynomial fitting method is as follows:

wherein b ₁ 、b ₂ Is the fitting coefficient.

S4, A is,θ ₀ 、θ(t)、Substituting roll channel dynamics equation f taking inertial coupling and pneumatic/kinematic coupling effects into account ₁ And will equation f ₁ Equation f reduced to Mathieu equation form ₂ The method comprises the steps of carrying out a first treatment on the surface of the Wherein roll channel dynamics equation f ₁ The following are provided:

wherein the method comprises the steps ofA pneumatic/kinematic coupling moment term;

gamma is an inertial coupling moment term;

is a rolling damping moment term;

I _y 、I _z respectively, about the aircraft body axis o-y _b 、o-z _b The moment of inertia of the body axis o-y _b 、o-z _b The definition is shown in fig. 2.

Roll channel dynamics equation f ₂ The following are provided:

wherein the method comprises the steps of

S5, solving an equation f through a Hill infinite determinant method ₂ Stability margin f (lambda, epsilon) |of (2) _λ＝1 ＝0、f(λ,ε)| _λ＝0.25 =0, and the stability margin was plotted in the λ - ε plane, and the stability profile was obtained as shown in fig. 4.

Wherein the stability margin f (lambda, epsilon) | _λ＝1 ＝0、f(λ,ε)| _λ＝0.25 The calculation method for=0 is as follows:

as shown in fig. 4, as known from the Floquet theory, the stability limit of the Mathieu equation can divide the λ -epsilon plane into a plurality of stable regions and unstable regions, and there is one stability limit starting from λ=0 on the λ axis, and two stability limits starting from the vicinity of λ=0.25 and λ=1, respectively. Wherein the stability equation starting at λ=0 has a periodic solution at the stability limit, wherein the stability limit starting at λ=0.25 has a solution of period 4pi and the stability limit starting at λ=1 has a solution of period 2pi. The periodic solution may be represented by a complex fourier series. The stability margin of the Mathieu equation can be obtained by the hilt infinite determinant method. The specific method comprises the following steps:

for stability limits around λ=1, the solution of the equation with a period of 2π is expressed as follows:

s in _n For each coefficient of the complex fourier series, formula (8) is substituted into formula (6) and is considered Then there are:

formula (9) can be simplified as

At e only ^inτ All coefficients of (2) are zero, i.e.:

to make the expression (11) be constructed with respect to { S ] _n The infinite system of equations of } has a non-zero solution, the value of the corresponding coefficient matrix determinant has to be zero, where λ+.n for any n is satisfied ² ；

The form of the coefficient matrix determinant is:

wherein:

similarly, for a stability limit around λ=0.25, the solution of the equation with a period of 4pi is expressed as follows:

the substitution of formula (14) into formula (6) is as follows:

formula (15) can be simplified as

Only atAll coefficients of (2) are zero, equation (9) can be true, namely:

the form of the coefficient matrix determinant is:

wherein:

the Hill determinant of the form (12) and the form (18) has infinite orders, namely N in the form should tend to infinity, but in the actual solving process, the infinite determinant is truncated by taking a finite N order, and N can be taken as 5.

Equation (12) and equation (18) are expressed as:

f(λ,ε,μ,δ)＝0 (20)

taking δ and μ as constants, equation (20) is a hidden function for λ, ε, and the curve f (λ, ε) of ε with λ around λ=1 and λ=0.25 is obtained by the "fimplitit" command of Matlab _λ＝1 ＝0、f(λ,ε)| _λ＝0.25 =0, i.e. stability margin.

S6, through f (lambda, epsilon) | _λ＝1 ＝0、f(λ,ε)| _λ＝0.25 Solution of critical frequency f=0 _c | _λ＝0 、 When the forced pitching frequency f meets f>f _c | _λ＝0 、/> In any case, the rolling motion diverges under the combined effect of inertial coupling and pneumatic/kinematic coupling.

Wherein f (lambda, epsilon) is passed _λ＝1 ＝0、f(λ,ε)| _λ＝0.25 Solution of critical frequency f=0 _c | _λ＝0 、 The calculation method of (2) is as follows:

for a given moment of inertia parameter a _x Pitch motion parameter θ ₀ Specific conditions of a, establish a function epsilon=epsilon (λ) of epsilon with respect to λ, namely:

as shown in fig. 4, epsilon=epsilon (λ) is plotted in the stability map as a working condition curve, which sequentially passes through stability limits around λ=0, λ=0.25, and λ=1 to form 5 intersection points, and coordinates of the intersection points are set to (λ _i ,ε _i ) From formula (7):

thereby obtaining the critical frequency f corresponding to the stability critical point _i Critical frequency f _i Sequentially defined as f in the order of the corresponding lambda value from smaller to larger _c | _λ＝0 、The coordinates of the intersection point of the operating mode curve and the stability limit and the corresponding frequencies are shown in table 1 below.

TABLE 1 intersection point coordinates of working condition curves and stability margins and corresponding frequencies

S7, selecting A=5°, θ ₀ =5°, repeating steps S1 to S6 to obtain the critical frequency f under the working condition _c | _λ＝0 、Based on the critical frequency results for the a=10° and a=5° operating conditions, an f-a stability plan can be plotted and plotted against CFD results such as shown in fig. 5. The shaded area in the figure is the unstable area predicted by theoretical analysis, each point in the figure represents a CFD working condition result, "·" represents that the CFD result is rolling motion convergence under the working condition, "×" represents that the CFD result is rolling motion divergence under the working condition, ". "means that the CFD results in limit cycle motion, i.e., critical steady state, under this condition. As can be seen in fig. 5, the numerical simulation results of the working conditions outside the shadow area are all convergent, and the numerical simulation results of the working conditions inside the shadow area are basically divergent, which indicates that the numerical simulation results are basically consistent with the theoretical prediction, which indicates that the stability prediction method disclosed by the application is effective in the working conditions studied herein, and the longitudinal/transverse coupling motion stability of the ultra-high speed aircraft in the given motion parameter range can be obtained by the method.

The foregoing is merely exemplary embodiments of the present application, and detailed technical solutions or features that are well known in the art have not been described in detail herein. It should be noted that, for those skilled in the art, several variations and modifications can be made without departing from the technical solution of the present application, and these should also be regarded as the protection scope of the present application, which does not affect the effect of the implementation of the present application and the practical applicability of the patent. The protection scope of the present application is subject to the content of the claims, and the description of the specific embodiments and the like in the specification can be used for explaining the content of the claims.

Claims (5)

1. A method for analyzing the stability of longitudinal and transverse coupling movement of a super-high-speed aircraft is characterized by comprising the following steps: the method comprises the following steps:

s1, aiming at the motion characteristics of an aircraft to be researched, given an incoming flow working condition and a pitch angle motion range of the aircraft, selecting a plurality of sample points { theta ] in the pitch angle motion range _i Then selecting the sinusoidal oscillation motion amplitude A and the average pitch angle theta of the pitch channel ₀ Thereby establishing a pitch channel motion equation θ (t);

s2, obtaining a sample point { theta ] through CFD calculation _i Static derivative of roll moment at }θ ₀ Rolling moment derivative at +.>

S3, obtaining a functional expression of the static derivative of the rolling moment relative to the pitch angle through polynomial fittingSubstituting θ (t) intoObtaining the functional expression->

S4, A and theta are calculated ₀ 、θ(t)、Substituting roll channel dynamics equation f taking inertial coupling and pneumatic/kinematic coupling effects into account ₁ And will equation f ₁ Equation f reduced to Mathieu equation form ₂ ；

S5, solving an equation f through a Hill infinite determinant method ₂ Is set at a constant limit f (lambda,ε)| _λ＝1 ＝0、f(λ，ε)| _λ＝0.25 ＝0；

s6, through f (lambda, epsilon) | _λ＝1 ＝0、f(λ，ε)| _λ＝0.25 Solution of critical frequency f=0 _c | _λ＝0 、 When the forced pitching frequency f satisfies +.> In any one of the above, the rolling motion diverges under the combined action of inertial coupling and pneumatic/kinematic coupling;

s7, selecting different motion amplitudes A and average pitch angles theta ₀ Repeating the steps S1-S6 to obtain the longitudinal/transverse coupling motion stability of the aircraft within the given motion parameter range;

the polynomial fitting method in step S3 uses the following fitting function:

wherein b ₁ 、b ₂ Fitting coefficients;

roll channel dynamics equation f in step S4 ₁ The following are provided:

wherein the method comprises the steps ofA pneumatic/kinematic coupling moment term;

is an inertial coupling moment term;

is a rolling damping moment term;

I _y 、I _z respectively, about the aircraft body axis o-y _b 、o-z _b Is a rotational inertia of (a);

roll channel dynamics equation f in step S4 ₂ The following are provided:

wherein the method comprises the steps of

2. The method for analyzing the longitudinal and transverse coupling motion stability of the ultra-high-speed aircraft according to claim 1, wherein the method comprises the following steps of: the pitch channel motion equation θ (t) in step S1 is as follows:

θ(t)＝θ ₀ +Asin(ωt) (1)

wherein θ is the pitch angle of the aircraft, θ ₀ The average pitch angle, A is the amplitude, ω is the pitch oscillation angular frequency, and t is the time.

3. The method for analyzing the longitudinal and transverse coupling motion stability of the ultra-high-speed aircraft according to claim 1, wherein the method comprises the following steps of: static derivative of roll moment in step S2Roll moment derivative +.>The calculation method of (2) is as follows:

wherein q is the incoming flow pressure, S _r For aircraft reference area, l _r For the reference length of the aircraft,is the pitch angle theta _i Rolling static derivative at ∈ ->Is the average pitch angle theta ₀ Roll derivative at I _X For winding the aircraft body axis o-x _b Is a rotational inertia of the bearing.

4. The method for analyzing the longitudinal and transverse coupling motion stability of the ultra-high-speed aircraft according to claim 1, wherein the method comprises the following steps of: solving equation f by Hill infinite determinant in step S5 ₂ Stability margin f (lambda, epsilon) |of (2) _λ＝1 ＝0、f(λ，ε)| _λ＝0.25 The calculation method for=0 is as follows:

for stability limits around λ=1, the solution of the equation with a period of 2π is expressed as follows:

s in _n For each coefficient of the complex fourier series, formula (8) is substituted into formula (6) and is considered Then there are:

formula (9) can be simplified as

At e only ^inτ All coefficients of (2) are zero, i.e.:

to make the constitution of { S } _n The infinite system of equations of } has a non-zero solution, the value of the corresponding coefficient matrix determinant has to be zero, where λ+.n for any n is satisfied ² ；

The form of the coefficient matrix determinant is:

wherein:

similarly, for a stability limit around λ=0.25, the solution of the equation with a period of 4pi is expressed as follows:

the substitution of formula (14) into formula (6) is as follows:

formula (15) can be simplified as

Only atAll coefficients of (2) are zero, equation (9) can be true, namely:

the form of the coefficient matrix determinant is:

wherein:

equation (12) and equation (18) are expressed as:

f(λ，ε，μ，δ)＝0 (20)

taking δ and μ as constants, equation (20) is a hidden function for λ, ε, and the curve f (λ, ε) of ε with λ around λ=1 and λ=0.25 is obtained by the "fimplitit" command of Matlab _λ＝1 ＝0、f(λ，ε)| _λ＝0.25 =0, i.e. stability margin.

5. The method for analyzing the longitudinal and transverse coupling motion stability of the ultra-high-speed aircraft according to claim 1, wherein the method comprises the following steps of: in step S6, f (lambda, epsilon) |is passed _λ＝1 ＝0、f(λ，ε)| _λ＝0.25 Solution of critical frequency =0The calculation method of (2) is as follows:

for a given moment of inertia parameter a _x Pitch motion parameter θ ₀ A to establish a function epsilon=epsilon (λ) of epsilon with respect to λ, i.e.:

drawing epsilon=epsilon (lambda) in a stability map as a working condition curve, wherein the working condition curve sequentially passes through stability limits near lambda=0, lambda=0.25 and lambda=1 to form 5 intersection points, and the coordinates of the intersection points are set as (lambda) _i ，ε _i ) From formula (7):

thereby obtaining the critical frequency f corresponding to the stability critical point _i Critical frequency f _i Sequentially defined as f in the order of the corresponding lambda value from smaller to larger _c | _λ＝0 、