CN117494611A - Supersonic flutter analysis method based on viscous local flow piston theory - Google Patents
Supersonic flutter analysis method based on viscous local flow piston theory Download PDFInfo
- Publication number
- CN117494611A CN117494611A CN202311542986.0A CN202311542986A CN117494611A CN 117494611 A CN117494611 A CN 117494611A CN 202311542986 A CN202311542986 A CN 202311542986A CN 117494611 A CN117494611 A CN 117494611A
- Authority
- CN
- China
- Prior art keywords
- piston
- flow field
- grid
- matrix
- wall
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Pending
Links
- 238000004458 analytical method Methods 0.000 title claims abstract description 38
- 239000011159 matrix material Substances 0.000 claims abstract description 49
- 238000000034 method Methods 0.000 claims abstract description 17
- 230000000694 effects Effects 0.000 claims abstract description 8
- 239000013598 vector Substances 0.000 claims description 15
- 238000010586 diagram Methods 0.000 claims description 14
- 238000013016 damping Methods 0.000 claims description 9
- 238000006243 chemical reaction Methods 0.000 claims description 6
- 238000006073 displacement reaction Methods 0.000 claims description 5
- 238000000605 extraction Methods 0.000 claims description 2
- 238000012546 transfer Methods 0.000 claims description 2
- 230000009466 transformation Effects 0.000 claims description 2
- 230000007704 transition Effects 0.000 claims description 2
- 238000004364 calculation method Methods 0.000 abstract description 13
- 238000012937 correction Methods 0.000 description 5
- 230000008878 coupling Effects 0.000 description 5
- 238000010168 coupling process Methods 0.000 description 5
- 238000005859 coupling reaction Methods 0.000 description 5
- 230000008859 change Effects 0.000 description 4
- 230000004048 modification Effects 0.000 description 4
- 238000012986 modification Methods 0.000 description 4
- 230000008901 benefit Effects 0.000 description 3
- 230000007547 defect Effects 0.000 description 3
- 238000013461 design Methods 0.000 description 3
- 238000004422 calculation algorithm Methods 0.000 description 2
- 238000012407 engineering method Methods 0.000 description 2
- 238000005192 partition Methods 0.000 description 2
- 238000004088 simulation Methods 0.000 description 2
- 230000009286 beneficial effect Effects 0.000 description 1
- 230000007423 decrease Effects 0.000 description 1
- 238000011161 development Methods 0.000 description 1
- 229910003460 diamond Inorganic materials 0.000 description 1
- 239000010432 diamond Substances 0.000 description 1
- 238000005516 engineering process Methods 0.000 description 1
- 230000006872 improvement Effects 0.000 description 1
- 230000003993 interaction Effects 0.000 description 1
- 239000000463 material Substances 0.000 description 1
- 230000002265 prevention Effects 0.000 description 1
- 238000011160 research Methods 0.000 description 1
- 239000007787 solid Substances 0.000 description 1
- 230000035899 viability Effects 0.000 description 1
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F30/00—Computer-aided design [CAD]
- G06F30/20—Design optimisation, verification or simulation
- G06F30/28—Design optimisation, verification or simulation using fluid dynamics, e.g. using Navier-Stokes equations or computational fluid dynamics [CFD]
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
- G06F17/10—Complex mathematical operations
- G06F17/11—Complex mathematical operations for solving equations, e.g. nonlinear equations, general mathematical optimization problems
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
- G06F17/10—Complex mathematical operations
- G06F17/16—Matrix or vector computation, e.g. matrix-matrix or matrix-vector multiplication, matrix factorization
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F30/00—Computer-aided design [CAD]
- G06F30/10—Geometric CAD
- G06F30/15—Vehicle, aircraft or watercraft design
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F30/00—Computer-aided design [CAD]
- G06F30/20—Design optimisation, verification or simulation
- G06F30/23—Design optimisation, verification or simulation using finite element methods [FEM] or finite difference methods [FDM]
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F2111/00—Details relating to CAD techniques
- G06F2111/10—Numerical modelling
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F2119/00—Details relating to the type or aim of the analysis or the optimisation
- G06F2119/14—Force analysis or force optimisation, e.g. static or dynamic forces
Landscapes
- Engineering & Computer Science (AREA)
- Physics & Mathematics (AREA)
- General Physics & Mathematics (AREA)
- Theoretical Computer Science (AREA)
- Mathematical Physics (AREA)
- Mathematical Analysis (AREA)
- Mathematical Optimization (AREA)
- Pure & Applied Mathematics (AREA)
- Geometry (AREA)
- Computational Mathematics (AREA)
- General Engineering & Computer Science (AREA)
- Data Mining & Analysis (AREA)
- Computer Hardware Design (AREA)
- Evolutionary Computation (AREA)
- Algebra (AREA)
- Computing Systems (AREA)
- Databases & Information Systems (AREA)
- Software Systems (AREA)
- Fluid Mechanics (AREA)
- Automation & Control Theory (AREA)
- Aviation & Aerospace Engineering (AREA)
- Operations Research (AREA)
- Aerodynamic Tests, Hydrodynamic Tests, Wind Tunnels, And Water Tanks (AREA)
Abstract
The invention belongs to the technical field of flutter analysis of supersonic velocity domain aircrafts, and provides a supersonic velocity flutter analysis method based on a viscous local flow piston theory. According to the method, an N-S equation is solved through CFD values, a large attack angle and a complex appearance interference flow field are obtained, and then the effective aerodynamic appearance of the aircraft is obtained according to a vorticity criterion; discretizing the surface of the aircraft to obtain a high-fidelity three-dimensional pneumatic grid suitable for a piston theory; interpolating flow field parameters and structural mode shape at the effective appearance of the CFD onto a piston grid by using a thin plate spline interpolation method, calculating a generalized aerodynamic force influence coefficient matrix and assembling a structural dynamics equation in a state space form; and carrying out flutter analysis by analyzing characteristic roots of the state space equation. The analysis method can analyze any complex aerodynamic shape, consider the viscosity effect of real gas, and consider the calculation precision and efficiency.
Description
Technical Field
The invention belongs to the technical field of flutter analysis of supersonic speed domain aircrafts, and particularly relates to a supersonic speed flutter analysis method based on a viscous local flow piston theory.
Background
Due to the remarkable advantage of high flying speed, the supersonic aircraft can effectively realize the characteristics of remote rapid deployment, high-efficiency burst prevention, high viability and the like, and has extremely high strategic value in the field of military aviation. The new generation supersonic and hypersonic aircrafts under development are widely applied to light materials and are designed by adopting large thin-wall structures, so that the inherent vibration frequency of the aircrafts is low, and the problem of fluid-solid coupling is more remarkable. Therefore, the flutter analysis has important significance for the design of supersonic and hypersonic aircrafts.
One of the main key technologies of supersonic flutter analysis is to determine unsteady aerodynamic forces, and the previous algorithms can be divided into two categories: the CFD/CSD coupling-based high-precision numerical simulation method and the approximate theory-based engineering method comprise a piston theory, a unified lifting surface theory and the like. The numerical simulation method has the defects of large calculated amount and low calculation efficiency, and is not suitable for preliminary design of the aircraft; and the engineering method is difficult to analyze the complex appearance and the large attack angle of the wing body assembly and the like. The Chinese patent publication No. CN104133933B discloses a hypersonic aeroelastic analysis method based on the piston theory, but the method needs to simplify the wing body assembly into a two-dimensional aerodynamic grid, and cannot consider the thickness effect of wings and a fuselage and the ternary effect of a low aspect ratio wing surface.
Based on the problems of the two methods, the local flow piston theory combining the high-precision CFD and the piston theory is developed in the last 90 th century, and the calculation efficiency and the precision can be considered. In 6 th month 1995, yang Bingyuan and Song Weili were first introduced to the local flow piston theory and flutter analysis in the "vibration and impact" volume 14, 2 nd publication "calculation of high angle of attack airfoil supersonic flutter with the local flow piston theory". In 9 2005, zhang Weiwei, she Zhengyin and the like published paper of "mechanics journal" volume 37, 5 ", based on the research of the aeroelastic calculation method of the local flow piston theory, the local flow piston theory is deduced again through the momentum theorem, and the theoretical basis of the popularization of the local flow piston theory from two-dimensional airfoil surfaces to three-dimensional complex shapes is laid. The Chinese patent publication No. CN106508028B describes a supersonic and hypersonic flutter analysis method which is based on the local flow piston theory and can be applied to complex shapes, but does not consider the influence of a viscous boundary layer on the effective aerodynamic shape, and only derives a generalized aerodynamic force calculation method of a symmetrical airfoil three-dimensional airfoil and a rotator body.
Aiming at the defects of the prior art, the invention realizes aerodynamic force calculation of any complex aerodynamic shape through a high-fidelity three-dimensional aerodynamic grid, introduces effective shape correction to consider the viscosity effect of real gas, and adopts a Kriging method to interpolate a near-wall flow field so as to improve the efficiency of the effective shape correction. On the basis, a supersonic flutter analysis method based on a viscous local flow piston theory is provided.
Disclosure of Invention
Aiming at the defects in the prior art, the invention provides a supersonic flutter analysis method based on a viscous local flow piston theory, which solves an N-S equation through CFD values to obtain a large attack angle and complex profile interference flow field, and further obtains the effective aerodynamic profile of an aircraft according to a vorticity criterion; discretizing the surface of the aircraft to obtain a high-fidelity three-dimensional pneumatic grid suitable for a piston theory; interpolating flow field parameters and structural mode shape at the effective appearance of the CFD onto a piston grid by using a thin plate spline interpolation method, calculating a generalized aerodynamic force influence coefficient matrix and assembling a structural dynamics equation in a state space form; and carrying out flutter analysis by analyzing characteristic roots of the state space equation. The analysis method can analyze any complex aerodynamic shape, consider the viscosity effect of real gas, and consider the calculation precision and efficiency.
The technical scheme of the invention is as follows:
a supersonic flutter analysis method based on a viscous local flow piston theory comprises the following steps:
step S1, a structural finite element model is established aiming at a target aircraft, and translational modes, generalized masses M and generalized rigidity K of all nodes on the surface of the aircraft are analyzed and extracted;
s2, aiming at a target aircraft, establishing a CFD model of an external flow field, solving an N-S equation, and deriving near-wall flow field parameters;
s3, interpolating the near-wall flow field to obtain flow field parameters of grid points of each wall, and calculating the effective appearance of the aircraft under the viscous effect according to the vorticity criterion;
s4, establishing an interpolation function according to a thin plate spline method to realize parameter transfer among different grids;
s5, calculating a generalized aerodynamic force influence coefficient matrix according to a local flow piston theory;
and S6, assembling a state space equation, and analyzing the flutter characteristics of the aircraft through the eigenvalues of the state transition matrix.
Further, the translational modal expression in the step S1 is as follows:
where i denotes the number of the aircraft surface node, j denotes the modal order,representing the translational component of the mode shape of the j-order mode of the surface node, u xij 、u yij And u zij Respectively representing the jth order modal array form of the ith aircraft surface node of the structure finite element model in the directions of the x axis, the y axis and the z axis of the structure finite element global coordinate systemA translational component of the direction;
the conversion is carried out under a piston grid coordinate system:
wherein L is ps Representing a transformation array of the structural finite element model coordinate system to the piston grid coordinate system.
Further, the near-wall flow field parameter expression in the step S2 is as follows:
P o ＝[p o ρ o v xo v yo v zo ξ o a o ]
wherein [ v ] xo ,v yo ,v zo ]For the gas flow velocity at flow field node o, p o 、ρ o 、ξ o 、a o The pressure, gas density, flow field vorticity and local sound velocity at the flow field node o, respectively.
Further, the step S3 specifically includes:
step S3-1, setting the coordinate of the wall grid point m of the CFD model as x m The normal vector of the wall surface is n m Extracting all flow field nodes above the wall grid point m in the near-wall flow field to be used as a data point set of flow field interpolation;
the extraction mode of the data point set is as follows:
wherein U is m Represents a data point set, d represents a flow field node in the data point set, x d Representing the coordinates of a flow field node d, I 3 Representing a 3 x 3 identity matrix, epsilon representing a preset threshold;
step S3-2, obtaining flow field parameters of the wall grid point m at any distance S along the normal direction by using a common kriking method, wherein the flow field parameters are expressed as follows:
wherein x is s Representing the coordinates at s, lambda d Representing weights, P d Representing the flow field parameter vector at the CFD near-wall flow field node d, the components of which define the near-wall flow field parameter P in step S2 o The same;
weight vector λ= [ λ ] d ]Obtained according to the following formula:
wherein e= [111]Is a column vector of all 1's, matrixThe covariance matrix between the known nodes is adopted, and the vector r is the covariance between the known nodes and the nodes to be solved;
s3-3, obtaining the effective appearance and the distance S between the wall surfaces according to the vorticity criterion c Expressed as:
wherein Ma ∞ 、c ∞ 、ρ ∞ 、T ∞ 、μ ∞ Mach number, sound velocity, density, temperature and dynamic viscosity, T, of the incoming flow far ahead respectively w T' is Anderson reference temperature, x m Representing the distance from the wall grid point m to the leading edge of the wing, L m Representing the reference length, i.e. the chord length of the wing passing through the wall grid point m, a 1 、a 2 Are effective shape fitting coefficients, and a is taken for the arc wing 1 ＝12.05,a 2 = -0.8, taking a for diamond wing 1 ＝9.05,a 2 ＝0.19，ξ c Is the vorticity criterion at the wall grid point m, ζ (x s ) The vorticity of the wall grid point m at any distance s along the normal direction;
s is gradually increased from 0 until ζ (x s )＝ξ c S at this time c I.e. the distance s between the effective profile and the wall c ；
S3-4, assigning flow field parameters at the effective appearance to corresponding wall grid points;
step S3-5, executing steps 3-1 to 3-4 on wall grid points of all CFD models;
s3-6, converting the flow field parameters of the grid points of the wall surface of the updated CFD model into a grid coordinate system of the piston;
expressed as:
wherein the flow field speed after conversion is
Converted coordinates are
Wherein [ v ] xm ,v ym ,v zm ]For the air flow velocity at the wall grid point m, p m 、ρ m 、ξ m 、a m Respectively the pressure, the gas density, the flow field vorticity and the local sound velocity at the wall grid point m, L pa Representing a conversion matrix of the CFD model coordinate system to the piston grid coordinate system,is the coordinates of the origin of coordinates of the CFD model in the piston grid coordinate system.
Further, the step S4 specifically includes:
interpolation is carried out by a thin plate spline method, and the flow field parameters of the grid points on the wall surface of the CFD model after updatingAnd interpolating to the piston grid, and obtaining a displacement interpolation matrix G from the structural finite element model to the piston grid.
Further, the generalized aerodynamic force influence coefficient matrix expression in the step S5 is as follows:
wherein F represents a aerodynamic force vector,and->Respectively representing a generalized pneumatic damping matrix and a generalized pneumatic stiffness matrix, wherein q represents a structural modal coordinate vector, phi is a modal translational vibration mode matrix at grid points of a wall surface of a structural finite element model, S=diag (S) is an area weighting matrix of a piston grid unit, and N=diag (N) is a diagonal matrix formed by normal directions N of the piston grid units; a is that 0 =diag (ρa) is a diagonal matrix consisting of the product of the local air density ρ and the sound velocity a of each piston grid cell; a is that k ＝A 0 V k K=1, 2,3, wherein +.>Local flow rate for each piston grid cell is x k Component of direction/>Diagonal matrix of components x 1 ,x 2 ,x 3 The directions correspond to the directions of the x axis, the y axis and the z axis respectively.
Further, the step S6 specifically includes:
in the step S6-1, the dynamic equation taking the structural modal coordinate q as the generalized displacement is in the following form in the state space:
wherein M and K are a generalized mass array and a generalized stiffness array respectively;
step S6-2, according to the theory of a linear system, the full necessary condition that the linear time-invariant system does not generate flutter is as follows: state matrix of systemAll eigenvalues λ of r (A) Are all positioned at the left half part of the complex plane, i.e
Re[λ r (A)]＜0,r＝1,2,...,n
Wherein n is the number of eigenvalues of the matrix A;
step S6-3, when calculating multiple working conditions, each speed state corresponds toAnd->And the corresponding state matrix and the characteristic value thereof, thereby obtaining the damping coefficient and the vibration frequency corresponding to each speed state;
and S6-4, drawing a root locus diagram of a real part and an imaginary part of a characteristic value corresponding to the speed from low to high, drawing a VG diagram of a damping coefficient corresponding to the characteristic value, drawing a VF diagram of a vibration frequency corresponding to the characteristic value, and searching a first unstable speed point to serve as a vibration speed boundary.
Compared with the prior art, the invention has the beneficial effects that:
1. compared with the existing CFD/CSD coupling algorithm, the supersonic flutter analysis method based on the viscous local flow piston theory can effectively solve the problems of large calculated amount and low calculation efficiency.
2. Compared with the existing supersonic flutter analysis method based on the piston theory, the supersonic flutter analysis method based on the viscous local flow piston theory provided by the invention considers the viscous effect of real gas, is suitable for high-fidelity three-dimensional pneumatic grids, introduces effective shape automatic correction based on flow field interpolation, and is more convenient for practical engineering application while improving the precision.
3. The supersonic flutter analysis method based on the viscous local flow piston theory is free from the limitations of complex aerodynamic shape, large attack angle flight state and the like, and can analyze the problem of thermoaerodynamic elastic stability by introducing a thermal mode, thus having important significance in the design of supersonic and hypersonic aircrafts.
Drawings
So that the manner in which the above recited embodiments of the present invention and the manner in which the same are attained and can be understood in detail, a more particular description of the invention, briefly summarized above, may be had by reference to the embodiments thereof which are illustrated in the appended drawings, which drawings are intended to be illustrative, and which drawings, however, are not to be construed as limiting the invention in any way, and in which other drawings may be obtained by those skilled in the art without the benefit of the appended claims.
FIG. 1 is a schematic diagram of a piston grid constructed in accordance with the present invention;
FIG. 2 is a schematic illustration of a first order modal translational matrix interpolation of a structure to a piston grid;
FIG. 3 is a schematic illustration of effective profile modification based on CFD model flow field parameters;
FIG. 4 is a V-G diagram reflecting the state matrix eigenvalue damping coefficient with the consequent flow velocity change;
FIG. 5 is a V-F diagram reflecting the change in state matrix eigenvalue vibration frequency with flow velocity.
Detailed Description
In order that the above-recited objects, features and advantages of the present invention will be more clearly understood, a more particular description of the invention will be rendered by reference to the appended drawings and appended detailed description. It should be noted that, without conflict, the embodiments of the present invention and features in the embodiments may be combined with each other.
In the following description, numerous specific details are set forth in order to provide a thorough understanding of the present invention, however, the present invention may be practiced in other ways than those described herein, and therefore the scope of the present invention is not limited to the specific embodiments disclosed below.
Example 1
Taking a certain full-motion control surface as an example, the supersonic flutter analysis method based on the viscous local flow piston theory provided by the invention is adopted for analysis.
The first step: and (5) building a structural finite element model of a certain full-motion control surface and carrying out modal analysis.
The natural frequencies of the model front 20 th order modes are as follows:
| modality | Frequency (Hz) | Modality | Frequency (Hz) |
| 1 | 49.17 | 11 | 2674.12 |
| 2 | 53.90 | 12 | 2763.49 |
| 3 | 78.13 | 13 | 3519.24 |
| 4 | 472.65 | 14 | 3635.11 |
| 5 | 493.93 | 15 | 3905.66 |
| 6 | 713.78 | 16 | 4168.34 |
| 7 | 1223.24 | 17 | 4469.27 |
| 8 | 1580.67 | 18 | 4692.98 |
| 9 | 1750.23 | 19 | 5269.25 |
| 10 | 2085.50 | 20 | 5866.80 |
For facilitating subsequent calculations, the modal patterns normalize the generalized mass.
And a second step of: and establishing a CFD model of a certain full-motion control surface, solving an N-S equation, and carrying out effective shape correction according to flow field parameters.
Incoming flow pressure q ∞ The CFD calculation conditions corresponding to =110 KPa are shown in the following table:
maintaining incoming flow Mach number M in flutter analysis ∞ Density ρ of incoming flow ∞ And wall temperature T w Unchanged, by increasing the incoming flow sound velocity a ∞ Controlling the flowing pressure, and determining the other parameters according to a gas state equation.
After the calculation is completed, substituting the flow field parameters into the following formula to obtain the boundary layer thickness s at each node c
Wherein the gas constant is R= 287.0J/(kg. Times.K), and the reference length is L m ＝0.284m。
And a third step of: and (3) establishing a piston grid, converting the structural finite element modal analysis result and the CFD model wall grid flow field parameters subjected to effective shape correction into a coordinate system where the piston grid is located, and interpolating the coordinate system onto the piston grid.
While pneumatic, CFD model grids and piston grids differ in that the flow field parameters of the CFD model are defined at nodes, while the local velocities, pressures, normal vectors, etc., required by piston theory are defined at the grids. Therefore, the interpolation of the modal displacement is from the wall surface node of the structural finite element model grid to the piston grid node; and interpolation of flow field parameters is from the effective profile nodes of the CFD model to the piston grid centroid.
The built piston grid is divided into an upper part and a lower part as shown in fig. 1, and each side grid is only subjected to aerodynamic force of one side of the outer direction.
Where n is the number of grid points for which the parameters are known; w (w) p Is a parameter component to be interpolated; epsilon is a given constant, called a parameter, for a generally flat function epsilon=10 -2 1, for an odd function, epsilon=10 is preferable -5 ～10 -6 ；D-th dimensional coordinates of the i-th node; r is (r) i 2 ＝||x-x i || 2 Is the distance between the point to be solved and the known point. />For the undetermined coefficients, it is determined by solving the following matrix:
in the middle ofWherein->Is the distance between two known nodes; h is a j The weighting coefficient corresponding to the j-th node is predetermined by the calculator. When all h j When=0, the fitted function passes exactly all node function values; when all h j At → infinity, the approximation function tends to fit to the least squares method.
As shown in FIG. 2, the interpolation result of the first-order mode array type is that the nodes of the structural finite element model grid and the piston grid are not in one-to-one correspondence, and the array type of the structural finite element model is required to be interpolated onto the nodes of the piston grid by using a thin plate interpolation theory. Because the applicable condition of the thin plate interpolation theory is that the known point set and the point set to be interpolated are positioned on the approximately overlapped smooth curved surfaces, the invention adopts a partition interpolation method, namely the upper wall surface and the lower wall surface of the structure finite element model and the upper wall surface and the lower wall surface of the piston grid correspond to each other to interpolate.
The effective profile modification is shown in fig. 3. Firstly, calculating the effective appearance of the CFD model according to the flow field parameters of the CFD model, then carrying out partition interpolation on the distance (namely the thickness of a boundary layer) between the wall surface grid of the CFD model and the effective appearance to the piston grid, and then constructing the corrected three-dimensional piston grid.
Fourth step: and calculating a generalized aerodynamic force influence coefficient matrix and a state matrix according to the piston theory.
The generalized aerodynamic force influence coefficient is calculated as follows:
wherein the method comprises the steps of
Where n is the number of piston grids, ρ l1 ，ρ l2 …ρ ln Local air density for each piston grid cell, a l1 ，a l2 …a ln For each piston grid cell sound velocity, V lxn Representing the local flow velocity component in the x-direction for each piston grid cell, V lyn Representing the component of the local flow velocity of each piston grid cell in the y-direction, V lzn Representing the component of the local flow rate in the z-direction for each piston grid cell.
The mode matrix is simplified because the generalized quality is normalized by the mode matrix type:
in the middle ofA diagonal array consisting of the modal vibration circle frequency (unit is rad/s) of the front 20 th order of the structure, I 20 Is a 20 x 20 identity matrix.
Fifth step: characteristic value analysis is carried out on state matrixes under different flowing pressures, and a V-G, V-F chart is drawn:
let matrix A (q ∞ ) Is lambda j Then the V-G diagram isWith q ∞ The curve of the first 10-order mode is shown in fig. 4. It can be seen from the graph that the damping ratio of the second-order mode increases and the damping ratio of the third-order mode increases and decreases, and the rest modes are basically unchanged. When coming flow dynamic pressure q ∞ When=110 KPa, the second order mode has traversed from below to above the horizontal axis, i.e., the analysis object has fluttered.
The V-F diagram isWith q ∞ The curves varied by the change are shown in fig. 5. Since the amplitude of the change of the higher order modes in the V-G graph is small, the first 3-order modes are focused on in the V-F graph. According to the flutter theory, when the two modes are coupled with each otherWhen the vibration frequency is increased, the vibration frequencies are close to each other, and the incoming flow pressure q of the analysis object can be judged ∞ The reason for chatter when=110 KPa is the coupling of the second and third order modes.
Since the state space equations are rank-increasing, each modality corresponds to two curves on the V-G, V-F plot. When q ∞ When not too large, the two curves correspond to a pair of conjugate complex roots and thus completely overlap on the V-G, V-F plot. The number of curves and the number of modes in fig. 4 and 5 are equal.
In conclusion, whether the analysis object flutters in the calculation range can be judged by observing whether the characteristic value of the real part larger than 0 appears in the root locus diagram; the critical speed of flutter and the crossing mode can be judged by observing the mode branch crossing the transverse axis in the V-G diagram; by observing the modal branches with close frequency in the V-F diagram, the coupling modal order and the flutter frequency causing the flutter of the aircraft can be judged.
In the present invention, unless explicitly specified and limited otherwise, the terms "mounted," "connected," "secured," and the like are to be construed broadly, and may be, for example, fixedly connected, detachably connected, or integrally formed; can be mechanically or electrically connected; can be directly connected or indirectly connected through an intermediate medium, and can be communicated with the inside of two elements or the interaction relationship of the two elements. The specific meaning of the above terms in the present invention can be understood by those of ordinary skill in the art according to the specific circumstances.
In the present invention, unless expressly stated or limited otherwise, a first feature "above" or "below" a second feature may include both the first and second features being in direct contact, as well as the first and second features not being in direct contact but being in contact with each other through additional features therebetween. Moreover, a first feature being "above," "over" and "on" a second feature includes the first feature being directly above and obliquely above the second feature, or simply indicating that the first feature is higher in level than the second feature. The first feature being "under", "below" and "beneath" the second feature includes the first feature being directly under and obliquely below the second feature, or simply means that the first feature is less level than the second feature.
In the present invention, the terms "first," "second," "third," "fourth" are used for descriptive purposes only and are not to be construed as indicating or implying relative importance. The term "plurality" refers to two or more, unless explicitly defined otherwise.
The above description is only of the preferred embodiments of the present invention and is not intended to limit the present invention, but various modifications and variations can be made to the present invention by those skilled in the art. Any modification, equivalent replacement, improvement, etc. made within the spirit and principle of the present invention should be included in the protection scope of the present invention.
Claims (7)
1. The supersonic flutter analysis method based on the viscous local flow piston theory is characterized by comprising the following steps of:
step S1, a structural finite element model is established aiming at a target aircraft, and translational modes, generalized masses M and generalized rigidity K of all nodes on the surface of the aircraft are analyzed and extracted;
s2, aiming at a target aircraft, establishing a CFD model of an external flow field, solving an N-S equation, and deriving near-wall flow field parameters;
s3, interpolating the near-wall flow field to obtain flow field parameters of grid points of each wall, and calculating the effective appearance of the aircraft under the viscous effect according to the vorticity criterion;
s4, establishing an interpolation function according to a thin plate spline method to realize parameter transfer among different grids;
s5, calculating a generalized aerodynamic force influence coefficient matrix according to a local flow piston theory;
and S6, assembling a state space equation, and analyzing the flutter characteristics of the aircraft through the eigenvalues of the state transition matrix.
2. The supersonic flutter analysis method based on the viscous local flow piston theory according to claim 1, wherein the translational modal expression in the step S1 is as follows:
where i denotes the number of the aircraft surface node, j denotes the modal order,representing the translational component of the mode shape of the j-order mode of the surface node, u xij 、u yij And u zij Respectively representing translational components of a jth order modal matrix form of an ith aircraft surface node of a structure finite element model in x-axis, y-axis and z-axis directions of a structure finite element global coordinate system;
the conversion is carried out under a piston grid coordinate system:
wherein L is ps Representing a transformation array of the structural finite element model coordinate system to the piston grid coordinate system.
3. The supersonic flutter analysis method based on the viscous local flow piston theory according to claim 1, wherein the near-wall flow field parameter expression in the step S2 is:
P o ＝[p o ρ o v xo v yo v zo ξ o a o ]
wherein [ v ] xo ,v yo ,v zo ]For the gas flow velocity at flow field node o, p o 、ρ o 、ξ o 、a o The pressure, gas density, flow field vorticity and local sound velocity at the flow field node o, respectively.
4. The supersonic flutter analysis method based on the viscous local flow piston theory according to claim 1, wherein the step S3 specifically comprises:
step S3-1, setting the coordinate of the wall grid point m of the CFD model as x m The normal vector of the wall surface is n m Extracting all flow field nodes above the wall grid point m in the near-wall flow field to be used as a data point set of flow field interpolation;
the extraction mode of the data point set is as follows:
wherein U is m Represents a data point set, d represents a flow field node in the data point set, x d Representing the coordinates of a flow field node d, I 3 Representing a 3 x 3 identity matrix, epsilon representing a preset threshold;
step S3-2, obtaining flow field parameters of the wall grid point m at any distance S along the normal direction by using a common kriking method, wherein the flow field parameters are expressed as follows:
wherein x is s Representing the coordinates at s, lambda d Representing weights, P d Representing a flow field parameter vector at a flow field node d;
weight vector λ= [ λ ] d ]Obtained according to the following formula:
wherein e= [11 … 1]Is a column vector of all 1's, matrixThe covariance matrix between the known nodes is adopted, and the vector r is the covariance between the known nodes and the nodes to be solved;
s3-3, obtaining the effective appearance and the distance S between the wall surfaces according to the vorticity criterion c Expressed as:
wherein Ma ∞ 、c ∞ 、ρ ∞ 、T ∞ 、μ ∞ Mach number, sound velocity, density, temperature and dynamic viscosity, T, of the incoming flow far ahead respectively w T' is Anderson reference temperature, x m Representing the distance from the wall grid point m to the leading edge of the wing, L m Representing the reference length, i.e. the chord length of the wing passing through the wall grid point m, a 1 、a 2 Are effective shape fitting coefficients, xi c Is the vorticity criterion at the wall grid point m, ζ (x s ) The vorticity of the wall grid point m at any distance s along the normal direction;
s is gradually increased from 0 until ζ (x s )＝ξ c S at this time c I.e. the distance s between the effective profile and the wall c ；
S3-4, assigning flow field parameters at the effective appearance to corresponding wall grid points;
step S3-5, executing steps 3-1 to 3-4 on wall grid points of all CFD models;
s3-6, converting the flow field parameters of the grid points of the wall surface of the updated CFD model into a grid coordinate system of the piston;
expressed as:
wherein the flow field speed after conversion is
Converted coordinates are
Wherein [ x ] m y m z m ]Coordinates v representing wall grid points m xm ,v ym ,v zm ]For the air flow velocity at the wall grid point m, p m 、ρ m 、ξ m 、a m Respectively the pressure, the gas density, the flow field vorticity and the local sound velocity at the wall grid point m, L pa Representing a conversion matrix of the CFD model coordinate system to the piston grid coordinate system,is the coordinates of the origin of coordinates of the CFD model in the piston grid coordinate system.
5. The supersonic flutter analysis method based on the viscous local flow piston theory according to claim 1, wherein the step S4 specifically comprises:
interpolation is carried out by a thin plate spline method, and the flow field parameters of the grid points on the wall surface of the CFD model after updatingAnd interpolating to the piston grid, and obtaining a displacement interpolation matrix G from the structural finite element model to the piston grid.
6. The supersonic flutter analysis method based on the viscous local flow piston theory according to claim 1, wherein the generalized aerodynamic force influence coefficient matrix expression in the step S5 is as follows:
wherein F represents a aerodynamic force vector,and->Respectively representing a generalized pneumatic damping matrix and a generalized pneumatic stiffness matrix, wherein q represents a structural modal coordinate vector, phi is a modal translational vibration mode matrix at grid points of a wall surface of a structural finite element model, S=diag (S) is an area weighting matrix of a piston grid unit, and N=diag (N) is a diagonal matrix formed by normal directions N of the piston grid units; a is that 0 =diag (ρa) is a diagonal matrix consisting of the product of the local air density ρ and the sound velocity a of each piston grid cell; a is that k ＝A 0 V k K=1, 2,3, wherein +.>Local flow rate for each piston grid cell is x k Component of direction->Diagonal matrix of components x 1 ,x 2 ,x 3 The directions correspond to the directions of the x axis, the y axis and the z axis respectively.
7. The supersonic flutter analysis method based on the viscous local flow piston theory according to claim 1, wherein the step S6 specifically comprises:
in the step S6-1, the dynamic equation taking the structural modal coordinate q as the generalized displacement is in the following form in the state space:
wherein M and K are a generalized mass array and a generalized stiffness array respectively;
step S6-2, according to the theory of a linear system, the full necessary condition that the linear time-invariant system does not generate flutter is as follows: state matrix of systemAll eigenvalues λ of r (A) Are all positioned at the left half part of the complex plane, i.e
Re[λ r (A)]＜0,r＝1,2,...,n
Wherein n is the number of eigenvalues of the matrix A;
step S6-3, when calculating multiple working conditions, each speed state corresponds toAnd->And the corresponding state matrix and the characteristic value thereof, thereby obtaining the damping coefficient and the vibration frequency corresponding to each speed state;
and S6-4, drawing a root locus diagram of a real part and an imaginary part of a characteristic value corresponding to the speed from low to high, drawing a VG diagram of a damping coefficient corresponding to the characteristic value, drawing a VF diagram of a vibration frequency corresponding to the characteristic value, and searching a first unstable speed point to serve as a vibration speed boundary.
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202311542986.0A CN117494611A (en) | 2023-11-20 | 2023-11-20 | Supersonic flutter analysis method based on viscous local flow piston theory |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202311542986.0A CN117494611A (en) | 2023-11-20 | 2023-11-20 | Supersonic flutter analysis method based on viscous local flow piston theory |
Publications (1)
| Publication Number | Publication Date |
|---|---|
| CN117494611A true CN117494611A (en) | 2024-02-02 |
Family
ID=89678023
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN202311542986.0A Pending CN117494611A (en) | 2023-11-20 | 2023-11-20 | Supersonic flutter analysis method based on viscous local flow piston theory |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN117494611A (en) |
Cited By (1)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN118153210A (en) * | 2024-05-10 | 2024-06-07 | 中国空气动力研究与发展中心高速空气动力研究所 | Nonlinear optimization design method for lifting line of diamond-shaped back missile wing and diamond-shaped back missile wing |
-
2023
- 2023-11-20 CN CN202311542986.0A patent/CN117494611A/en active Pending
Cited By (1)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN118153210A (en) * | 2024-05-10 | 2024-06-07 | 中国空气动力研究与发展中心高速空气动力研究所 | Nonlinear optimization design method for lifting line of diamond-shaped back missile wing and diamond-shaped back missile wing |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN105843073B (en) | A kind of wing structure aeroelastic stability analysis method not knowing depression of order based on aerodynamic force | |
| Carlson | Transonic airfoil analysis and design using Cartesian coordinates | |
| CN108427322B (en) | Modeling method of large flexible aircraft based on online identification | |
| Smith et al. | CFD-based analysis of nonlinear aeroelastic behavior of high-aspect ratio wings | |
| CN117494611A (en) | Supersonic flutter analysis method based on viscous local flow piston theory | |
| CN103646131A (en) | Multi-objective optimization design method considering aeroelasticity restraint and for wing made of composite material | |
| Hu et al. | Adjoint-based aerodynamic optimization of supersonic biplane airfoils | |
| CN109840349B (en) | Fixed-wing aircraft gust response modeling analysis method | |
| CN113569360B (en) | Method for designing flutter-resistant airfoil cluster of wind turbine blade | |
| CN103810341A (en) | Predicating method for wind turbine blade airfoil fluttering | |
| CN110287505B (en) | Aircraft stability analysis method | |
| KR102616901B1 (en) | Aerodynamic layout design method and system for wide-speed-range hypersonic aircraft | |
| CN106096088A (en) | A kind of propeller aeroplane WHIRL FLUTTER ANALYSIS method | |
| CN113868758B (en) | Elastic hypersonic aircraft modeling method based on computational mechanics | |
| CN110287643A (en) | Cavity Flow vibration and stream cause the wind tunnel experiment analogy method of noise coupling characteristic | |
| Hammerton et al. | Optimum distributed wing shaping and control loads for highly flexible aircraft | |
| Asadi et al. | Flutter optimization of a wing–engine system with passive and active control approaches | |
| Nguyen et al. | Flutter analysis of the transonic truss-braced wing aircraft using transonic correction | |
| CN114154434A (en) | Multi-constraint refined pneumatic optimization design method for layout of tailless flying wing | |
| CN112577657A (en) | Method for quickly predicting pulsating load generated by separation shock wave oscillation | |
| CN110705154B (en) | Optimization method for balanced order reduction of open-loop pneumatic servo elastic system model of aircraft | |
| CN115358167B (en) | Flying and launching integrated pneumatic accompanying optimization design method considering engine parameters | |
| CN114036628A (en) | Method for collaborative design of wingspan and control strategy of morphing aircraft | |
| CN115964795A (en) | Deformation control method for morphing aircraft based on disturbance observer | |
| CN114545778A (en) | Helicopter robust tracking control method based on LPV model |
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 |