CN113254863B - Adaptive orthogonal decomposition method for model reduction of multi-disk rotor system - Google Patents
Adaptive orthogonal decomposition method for model reduction of multi-disk rotor system Download PDFInfo
- Publication number
- CN113254863B CN113254863B CN202110492921.4A CN202110492921A CN113254863B CN 113254863 B CN113254863 B CN 113254863B CN 202110492921 A CN202110492921 A CN 202110492921A CN 113254863 B CN113254863 B CN 113254863B
- Authority
- CN
- China
- Prior art keywords
- matrix
- equation
- expressed
- order
- pom
- 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.)
- Active
Links
- 238000000034 method Methods 0.000 title claims abstract description 76
- 238000000354 decomposition reaction Methods 0.000 title claims abstract description 17
- 230000003044 adaptive effect Effects 0.000 title claims description 8
- 230000009467 reduction Effects 0.000 title claims description 8
- 239000011159 matrix material Substances 0.000 claims abstract description 62
- 230000004044 response Effects 0.000 claims description 38
- 238000001907 polarising light microscopy Methods 0.000 claims description 33
- 229920006324 polyoxymethylene Polymers 0.000 claims description 33
- 238000006073 displacement reaction Methods 0.000 claims description 26
- 230000033001 locomotion Effects 0.000 claims description 12
- 239000013598 vector Substances 0.000 claims description 12
- 230000009466 transformation Effects 0.000 claims description 11
- 238000005070 sampling Methods 0.000 claims description 10
- 230000001052 transient effect Effects 0.000 claims description 8
- 239000012530 fluid Substances 0.000 claims description 5
- 238000013016 damping Methods 0.000 claims description 4
- 230000005484 gravity Effects 0.000 claims description 4
- 238000012886 linear function Methods 0.000 claims description 3
- 239000000314 lubricant Substances 0.000 claims description 3
- 238000005461 lubrication Methods 0.000 claims description 3
- 238000004364 calculation method Methods 0.000 abstract description 11
- 238000010606 normalization Methods 0.000 abstract description 3
- 238000010586 diagram Methods 0.000 description 3
- 230000010354 integration Effects 0.000 description 2
- 238000004458 analytical method Methods 0.000 description 1
- 230000009286 beneficial effect Effects 0.000 description 1
- 230000000739 chaotic effect Effects 0.000 description 1
- 230000007423 decrease Effects 0.000 description 1
- 230000007547 defect Effects 0.000 description 1
Classifications
-
- 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
- G06F17/13—Differential equations
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01M—TESTING STATIC OR DYNAMIC BALANCE OF MACHINES OR STRUCTURES; TESTING OF STRUCTURES OR APPARATUS, NOT OTHERWISE PROVIDED FOR
- G01M13/00—Testing of machine parts
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01M—TESTING STATIC OR DYNAMIC BALANCE OF MACHINES OR STRUCTURES; TESTING OF STRUCTURES OR APPARATUS, NOT OTHERWISE PROVIDED FOR
- G01M13/00—Testing of machine parts
- G01M13/04—Bearings
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01M—TESTING STATIC OR DYNAMIC BALANCE OF MACHINES OR STRUCTURES; TESTING OF STRUCTURES OR APPARATUS, NOT OTHERWISE PROVIDED FOR
- G01M13/00—Testing of machine parts
- G01M13/04—Bearings
- G01M13/045—Acoustic or vibration analysis
-
- 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
-
- Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
- Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
- Y02T—CLIMATE CHANGE MITIGATION TECHNOLOGIES RELATED TO TRANSPORTATION
- Y02T90/00—Enabling technologies or technologies with a potential or indirect contribution to GHG emissions mitigation
Landscapes
- Physics & Mathematics (AREA)
- General Physics & Mathematics (AREA)
- Engineering & Computer Science (AREA)
- Mathematical Physics (AREA)
- Computational Mathematics (AREA)
- Mathematical Analysis (AREA)
- Mathematical Optimization (AREA)
- Pure & Applied Mathematics (AREA)
- Data Mining & Analysis (AREA)
- Theoretical Computer Science (AREA)
- Software Systems (AREA)
- Databases & Information Systems (AREA)
- Algebra (AREA)
- General Engineering & Computer Science (AREA)
- Operations Research (AREA)
- Acoustics & Sound (AREA)
- Computing Systems (AREA)
- Complex Calculations (AREA)
Abstract
The invention provides a multi-disk rotor system model reduced order self-adaptive orthogonal decomposition method, which is characterized in that a solution of a system dynamics control equation under different modeling parameters is obtained in advance through a fourth-order Dragon library tower method, then an autocorrelation matrix is obtained according to a snapshot matrix formed by the solution, further a reduced order POM under new parameters is obtained through calculating characteristic values of the autocorrelation matrix, and finally the complex original system dynamics control equation under the new model parameters is simplified through the reduced order POM. The method can obtain the optimal reduced order model of the complex rotor system, can maintain the orthogonal normalization characteristic of the reduced order POM, eliminates the local characteristic of the ITGM method, is used in a wide frequency parameter range, and has accurate calculation, higher calculation efficiency and strong robustness.
Description
Technical Field
The invention relates to the field of dynamics, in particular to a model reduced-order orthogonal decomposition method.
Background
The rotor bearing system is one of the most important components in mechanical systems such as aircraft engines, gas turbines, large generators, and the like. Typically, rotor bearing systems include multiple rotating disks, which are very complex and highly dimensional, and there is a need to reduce rotor bearing system dimensions using suitable system dimension reduction methods.
The POD method obtains the most important components of an infinite-dimension continuous system or a finite-dimension DOFs (Degrees of Freedom) system through an inherent orthogonal mode (POMs). The POD method can greatly reduce the degree of freedom of the complex system, and can improve the calculation efficiency of the complex system on the premise of keeping higher calculation precision, so that the POD method is widely applied to the related field of fluid mechanics. However, the Reduced Order Model (ROM) obtained by the POD method is often lacking in robustness when system parameters are changed. Compared with the IGM method, the ITGM method has more reduced modes required by the approximation of the original system and has poorer approximation precision.
Disclosure of Invention
In order to overcome the defects of the prior art, the invention provides a multi-disk rotor system model reduced-order self-adaptive orthogonal decomposition method, namely an IGM method. The method comprises the steps of obtaining solutions of system dynamics control equations under different modeling parameters in advance through a fourth-order Dragon library tower method, obtaining an autocorrelation matrix according to a snapshot matrix formed by the solutions, obtaining a reduced order POM under new parameters through calculating characteristic values of the autocorrelation matrix, and simplifying the complex original system dynamics control equations under the new model parameters through the reduced order POM.
In order to solve the above problems in the POD method, the present invention proposes an adaptive POD method (IGM method) based on the glasman manifold and the lagrangian interpolation. Firstly, a fourth-order Dragon-Gregory tower method is adopted to obtain solutions of system dynamics control equations under different modeling parameters in advance, then an autocorrelation matrix is obtained according to a snapshot matrix of the solutions, a reduced order POM under new parameters can be obtained through the autocorrelation matrix, and finally the complex original system dynamics control equations under the new model parameters are simplified through the reduced order POM. The result of analysis of the model response can show that the IGM method not only obtains the optimal ROM (Reduced-ordered Models) of the complex physical system, maintains the orthogonal normalization characteristic of the Reduced POM, but also eliminates the limitation of the ITGM method. The IGM method is accurate in calculation in a wider parameter domain, high in calculation efficiency and high in robustness.
The technical scheme adopted by the invention for solving the technical problems comprises the following steps:
step one, obtaining a motion differential equation of a high-dimensional nonlinear rotor system;
the rotor bearing system is a multi-disc rotor bearing system and is supported by two oil lubrication short shaft journal bearings at two ends, a base loosening fault occurs at a left support part, and 14 discs are arranged on a shaft;
the shafts are elastic, the mass of each section of shaft is concentrated at the mass center of the corresponding disc, the gyro moment of the disc is not considered, and the equivalent rigidity of the shafts between the turntables is k i I=1..15, the horizontal and vertical displacement of the geometric center of the disk and the center of mass of the support end is denoted as X i ,Y i I=1..16, nonlinear oil film forces in horizontal and vertical directions are denoted as F x ,F y Vibration of loosening failure occurs only in the vertical direction, the corresponding displacement is Y 17 Base k l And base c) l Is a piecewise linear function, expressed as:
the differential equation of motion of the rotor bearing system is expressed as:
wherein X= [ X ] 1 ,...X 16 ,Y 1 ...Y 16 ,Y 17 ]M, C, K are mass, damping and stiffness matrices, respectively, F is the resultant of nonlinear fluid film forces, imbalance forces and gravity of each disk, expressed as:
M＝diag(m 1 ,…,m 16 ,m 1 ,…m 16 ,m 17 )
C＝diag(c 1 ,…,c 16 ,c 1 ,…c 16 ,c 17 )
carrying out dimensionless transformation on a dynamics control equation of the system, wherein the formula (3) is dimensionless by the following formula:
wherein f x And f y Is a dimensionless nonlinear oil film force,is the Sommerfeld number, μ is the lubricant viscosity, c is the oil film thickness, L is the bearing length, R is the bearing radius, ω is the angular frequency, τ is the dimensionless time, and W is half the weight of the rotor system;
thereby establishing a high-dimensional nonlinear rotor system model;
2. inputting the response signals into the snapshot matrix and calculating an autocorrelation matrix;
manifold G (M, N) s ) Is provided with a plurality of sampling time sequences N s Covering M-dimensional column vectors of corresponding solutions;
n pre-computed snapshot signal sets of M-dimensional system are set toω i ∈[ω α ,ω β ]Wherein each parameter omega i Sampling length t at s Same, snapshot set->Seen as manifold G (M, N s ) N points->Selecting a point S of the manifold 0 As a reference point, by Lagrange interpolationThe method obtains a snapshot matrix y with corresponding parameter omega on manifold * (ω,t s ) Expressed as:
manifold G (M, N) s ) The interpolation points above are expressed as:
wherein Σ is i ∈R M×M Is a matrix of singular values that are used,V i ∈R M×M is an orthogonal matrix, then equation (7) is rewritten as:
let y (omega) i ,t s ) The solution representing the original system with parameter value ω also belongs to manifold G (M, N s ) By thin singular value decomposition, y (ω) i ,t s ) Expressed as:
y(ω,t s )＝UΣV T (9)
according to Lagrange interpolation method, the truncation error epsilon N (ω,t s ) The definition is as follows:
wherein the superscript (n+1) denotes the (n+1) -th derivative of the parameter ω, denoted asWhen N → infinity, the truncation error ε N (ω,t s )→0；y(ω i ,t s ) The rewriting is as follows:
in connection with equation (9), equation (11) is expressed as:
diagonal matrix Σ=diag [ σ ] 1 ,σ 2 ,...,σ M ]Is y (omega) i ,t s ) Singular values of (1), wherein sigma 1 ≥σ 2 ≥,...,≥σ M 0 or more, the singular value energy ratio ε is expressed as:
UΣV T is divided into main partsAnd a truncated error fraction epsilon m (ω,t s ) Substituting the main part and the error part into formula (12) to obtain:
wherein ε m (ω,t s ) And epsilon N (ω,t s ) Are minor and omitted, and therefore, equation (14) is rewritten as:
autocorrelation matrix S c The definition is as follows:
the dimensionless control equation of the original system is solved by a fourth-order Runge-Kutta (RK 4) method;
calculating displacement response by inputting different modeling parameters and storing response signals in a snapshot matrixCalculating an autocorrelation matrix S according to equation (16) c The method comprises the following steps:
3. by means of an autocorrelation matrix S c Obtaining the reduced order POM with new parameters
The reduced order POM of the parameter ω is calculated,expressed as:
then equation (17) is rewritten as:
by calculating a correlation matrix S c The feature vector of (2) obtaining reduced order POMs;
4. reduced POM for dynamic differential equation reduction
The reduced POM is defined by the eigenvector of the autocorrelation matrix, sampled at a rotational speed ω' with a sampling length t s Order-makingFeature vectors representing the autocorrelation matrix, the feature vectors being arranged in descending order of the corresponding feature values, provided thatRepresenting a linear subspace having n reduced order POMs;
transient response of complex high-dimensional space is projected onto low-dimensional subspace through coordinate transformation:
wherein,u(t)＝[u 1 (t),u 2 (t),...,u n (t)] T />the dimensionless displacement response, the coordinate transformation matrix and the POD modal coordinates of the original system are respectively, and the ordinary differential equation (10) is simplified into:
for the IGM method with the same ROM degree of freedom and k points, singular value decomposition is only required to be executed once, so that the computing efficiency of the IGM method is far higher than that of the ITGM method.
The singular value energy ratio epsilon is epsilon more than or equal to 95 percent.
When the dimensionless control equation is solved by a fourth-order Runge-Kutta (RK 4) method, the initial displacement and the velocity are respectively X i ＝Y i ＝Y 17 ＝0.1，The time step is 0.001 pi.
The invention has the beneficial effects that:
(1) The invention not only can obtain the optimal reduced order model of the complex rotor system, but also can maintain the orthogonal normalization characteristic of the reduced order POM, and eliminates the local characteristic of the ITGM method.
(2) The method can be used in a wide frequency parameter range, and has the advantages of accurate calculation, high calculation efficiency and strong robustness. The method is not only suitable for rotor bearing models, but also suitable for other dynamic models.
Drawings
FIG. 1 is a rotor bearing system of a 600MW turbine engine.
Figure 2 is a schematic view of a multi-disc rotor bearing system,
FIG. 3 is a graph of adaptive interpolation based on the Grassman manifold, and FIG. 3 (a) is a four-point interpolation of the ITGM method; fig. 3 (b) is a four-point interpolation of IGM method.
Fig. 4 is a comparison of IGM and ITGM methods of two-point interpolation at rotational speed ω=550 rad/s: FIG. 4 (a) shows the displacement response, the horizontal axis shows the number of rotations, and the vertical axis shows the longitudinal displacement; fig. 4 (b) is an axial trace diagram, and the abscissa is the axial displacement.
Fig. 5 is a comparison of IGM and ITGM methods of three-point interpolation at rotational speed ω=550 rad/s, with displacement response of fig. 5 (a), rotation times on the horizontal axis and longitudinal displacement on the vertical axis; fig. 5 (b) is an axial trace diagram, and the abscissa represents the axial displacement.
Fig. 6 is a comparison of IGM and ITGM methods of five-point interpolation at rotational speed ω=550 rad/s: fig. 6 (a) shows displacement response, the horizontal axis shows the number of rotations, and the vertical axis shows longitudinal displacement; fig. 6 (b) is an axial trace diagram, and the abscissa is the axial displacement.
Fig. 7 is a graph showing the convergence of the ROM response of the IGM method at a rotational speed ω=550 rad/s, fig. 7 (a) is a graph showing the convergence of the 2-point ROM response, fig. 7 (b) is a graph showing the convergence of the 3-point ROM response, fig. 7 (c) is a graph showing the convergence of the 5-point ROM response, and the horizontal axis shows the number of rotations and the vertical axis shows the longitudinal displacement.
Detailed Description
The invention will be further described with reference to the drawings and examples.
The method comprises the following specific implementation steps:
firstly, a high-dimensional nonlinear rotor system model is established, the transient response of a rotor is calculated through a fourth-order Dragon-grid-tower method under different modeling parameter values, and response signals are input into a snapshot matrixIs a kind of medium.
Subsequently, an autocorrelation matrix S is calculated c And pass through an autocorrelation matrix S c To obtain the reduced-order POM after the model parameters are changed.
And finally, deriving a coordinate transformation matrix by using the reduced order POM, and simplifying a complex high-dimensional nonlinear rotor system differential equation through coordinate transformation to achieve the purpose of reducing the order.
Step one, obtaining a motion differential equation of a high-dimensional nonlinear rotor system
FIG. 1 is a rotor bearing system of a 600MW turbine engine, which is a multi-disc rotor bearing system as shown in FIG. 2, with a base loosening failure at the left support of FIG. 2, with 14 discs on the shaft supported by two oil lubricated stub journal bearings at both ends. Is supported by two oil lubrication short journal bearings at two ends, a base loosening fault occurs at the left support, and 14 discs are arranged on the shaft.
Assuming that the axes are elastic and that the mass of each segment of the axis is concentrated at the center of mass of the corresponding disk, the gyroscopic moment of the disk is not considered because the distribution of the disk is symmetrical about the midpoint of the axis with less influence on the system response. O (O) 1 ,O 16 Representing the geometric centers of the left and right bearings, respectively. The geometric center of each disc is O i (i=2..15), centroid O (i=2..15), eccentricity e i (i=2..15). The equivalent mass of each disc and the corresponding shaft mass are m i (i=1..16) indicates that the equivalent damping of each disc is c i (i=1..16). In addition, the equivalent stiffness of the shaft between the turntables is k i (i=1..15) horizontal and vertical displacement of the geometric center of the disk and the center of mass of the support end is denoted as X i ,Y i (i=1..16), the nonlinear oil film force in the horizontal and vertical directions is denoted as F x ,F y . Assuming vibration in which a loosening failure occurs only in the vertical direction, the corresponding displacement is Y 17 . Base k l And base c) l Is typically a piecewise linear function, expressed as:
considering the nonlinear fluid film forces, the unbalance forces, and the overall system gravity, the differential equation of motion of the rotor bearing system is expressed as:
wherein X= [ X ] 1 ,...X 16 ,Y 1 ...Y 16 ,Y 17 ]M, C, K are mass, damping and stiffness matrices, respectively, F is the resultant of nonlinear fluid film forces, imbalance forces and gravity of each disk, expressed as:
M＝diag(m 1 ,…,m 16 ,m 1 ,…m 16 ,m 17 )
C＝diag(c 1 ,…,c 16 ,c 1 ,…c 16 ,c 17 )
in order to facilitate solving the transient response of a multi-degree of freedom system, the system dynamics control equation is usually subjected to dimensionless transformation, and the formula (3) is dimensionless by the following formula:
wherein f x And f y Is non-dimensional nonlinear oil film force，Is the Sommerfeld number, μ is the lubricant viscosity, c is the oil film thickness, L is the bearing length, R is the bearing radius, ω is the angular frequency, τ is the dimensionless time, and W is half the weight of the rotor system.
Obtaining the high-dimensional nonlinear rotor system model established in the first step
2. Inputting the response signals into the snapshot matrix and calculating the autocorrelation matrix
Based on the theory of the dynamics system, the entire solution set of the second order ordinary differential equation is micro-manifold, and the order of the micro-manifold is larger than 2. Thus, the snapshot matrix of responses under different parameters can be considered as one point on the glasman manifold. Manifold G (M, N) s ) Is provided with a plurality of sampling time sequences N s The M-dimensional column vectors of the corresponding solutions are overlaid.
N pre-computed snapshot signal sets of M-dimensional system are set toω i ∈[ω α ,ω β ]Wherein each parameter omega i Sampling length t at s The same applies. Then snapshot set +.>Can be regarded as manifold G (M, N s ) N points->Selecting a point S of the manifold 0 As a reference point, a snapshot matrix y with a corresponding parameter omega on a manifold can be obtained by a Lagrange interpolation method * (ω,t s ) As shown in fig. 3. Expressed as:
according to Grassmann manifold correlation theory in differential geometryPrinciple of outlier decomposition, manifold G (M, N s ) The interpolation points above are expressed as:
y(ω i ,t s )＝U i Σ i V i t ,i＝0,...,N-1 (28)
wherein Σ is i ∈R M×M Is a matrix of singular values that are used,V i ∈R M×M is an orthogonal matrix. Then equation (7) may be rewritten as:
let y (omega) i ,t s ) Representing a solution of the original system with parameter value ω, which also belongs to manifold G (M, N s ). Using thin singular value decomposition method, y (ω) i ,t s ) Expressed as:
y(ω,t s )＝UΣV T (30)
according to Lagrange interpolation method, the truncation error epsilon N (ω,t s ) The definition is as follows:
wherein the superscript (n+1) denotes the (n+1) -th derivative of the parameter ω, denoted asWhen N → infinity, the truncation error ε N (ω,t s ) The more interpolation points, the closer to the solution of the original system at the parameter ω, is indicated by 0. y (omega) i ,t s ) Rewritten as
In connection with equation (9), equation (11) is expressed as:
diagonal matrix Σ=diag [ σ ] 1 ,σ 2 ,...,σ M ]Is y (omega) i ,t s ) Singular values of (1), wherein sigma 1 ≥σ 2 ≥,...,≥σ M And is more than or equal to 0. The singular value energy ratio ε is expressed as:
the energy contribution criterion is typically epsilon 95% or more, even epsilon 99% or more, and the value of m is much smaller than the dimension of the original system (this condition is typically met). Therefore, U ΣV T Is divided into main partsAnd a truncated error fraction epsilon m (ω,t s ) Substituting the main part and the error part into formula (12) to obtain:
wherein ε m (ω,t s ) And epsilon N (ω,t s ) All minor parts may be omitted, so equation (14) is rewritten as:
the dimension is high, resulting in a snapshot matrix y * (ω,t s ) The singular value decomposition of (c) can take a significant amount of time and memory. Autocorrelation matrix S c The definition is as follows:
the dimensionless control equation of the original system is solved by a fourth-order Runge-Kutta (RK 4) method, which is a direct numerical integration method, and the initial displacement and the velocity are respectively X i ＝Y i ＝Y 17 ＝0.1，The time step is 0.001 pi. Calculating displacement response by inputting different modeling parameters and storing response signals in a snapshot matrixThe snapshot matrix has N s The snapshot matrix sampled from the transient response in row 33 columns may reduce the model order better. Because the transient response of the system contains free vibration, forced vibration, and more kinetic information in the original system. The present invention therefore chooses to obtain snapshot matrices with different frequency parameter values from the transient response. The autocorrelation matrix S can be calculated according to the formula (16) c ：
3. By means of an autocorrelation matrix S c Obtaining the reduced order POM with new parameters
The reduced order POM of the parameter ω is calculated,expressed as:
then equation (17) is rewritten as:
as can be seen from the above, by calculating the correlation matrix S c The feature vectors of (2) obtain reduced order POMs. Although snapshot matrixNot the solution of the original system under the parameter ω, but it is located in the vicinity of the solution of the original system. Furthermore, due to the snapshot matrix->The corresponding interpolation points in the wide parameter area contain various motion states (such as quasi-periodic motion, chaotic motion), thus, the snapshot matrix is +.>The obtained reduced-order POM has more abundant modal information than the case of a normal snapshot matrix. Thus, the adaptive POD method not only can obtain the optimal POM of a complex physical system, but also can maintain the orthonormal of the reduced POM. Meanwhile, the method can be used in a wide parameter range due to the fact that the limitation of local space properties of point tangent lines on the manifold is eliminated.
4. Reduced POM for dynamic differential equation reduction
The reduced order POM is defined by the eigenvectors of the autocorrelation matrix. Let us assume that the sampling is performed at a rotational speed ω' with a sampling length t s . Order theFeature vectors representing the autocorrelation matrix, the feature vectors being arranged in descending order of the corresponding feature values. Is provided with->Representing a linear subspace with n reduced order POMs.
Transient response of a complex high-dimensional space can be projected onto a low-dimensional subspace through coordinate transformation:
wherein,u(t)＝[u 1 (t),u 2 (t),...,u n (t)] T />the dimensionless displacement response, the coordinate transformation matrix and the POD modal coordinates of the original system are respectively. The ordinary differential equation (10) is reduced to:
5. calculation results of IGM method
For the simple motion case, as shown in fig. 4 to 6, the displacement response and the axis locus of the ROM are obtained in 2-point, 3-point, 5-point interpolation by IGM and ITGM methods, and compared with the numerical results at ω=500 rad/s rotation speed. In fig. 4, 17 reduced-order POMs are obtained by IGM method during two-point interpolation, and then the calculated responses are consistent with the response of the original system. In fig. 5, in the case of three-point interpolation, the IGM method is used, and only 12 reduced-order POMs are used to calculate the response consistent with the original system. In fig. 6, only 4 reduced-order POMs are needed to perform 5-point interpolation to perfectly match the response of the original system.
Converging response curves of ROMs with different degrees of freedom at ω=550 rad/s rotational speed were obtained by IGM methods of 2, 3 and 5 points, respectively, as shown in fig. 7. It can be derived that the more degrees of freedom the ROM, the more accurate the response is obtained using the IGM method at 3 different interpolation points. In addition, as the number of interpolation points increases, the DOF number of the ROM decreases.
During complex locomotion, ITGM and IGM methods can be found to be equivalent in the local parameter region. However, since the ITGM method cannot compute dense interpolation points over a wider parameter area, it is not suitable for complex motion in the case of a wider parameter area.
For the same ROM degree of freedom, the k-point ITGM method needs to execute the singular value decomposition of k+1, while the k-point IGM method only needs to execute the singular value decomposition once, and meanwhile, the calculation time of the IGM five-point method is faster than the numerical integration speed of an original system by more than 20 times, so that the calculation efficiency of the IGM method is far higher than that of the ITGM method.
Claims (3)
1. The adaptive orthogonal decomposition method for reducing the model order of the multi-disk rotor system is characterized by comprising the following steps of:
step one, obtaining a motion differential equation of a high-dimensional nonlinear rotor system;
the rotor bearing system is a multi-disc rotor bearing system and is supported by two oil lubrication short shaft journal bearings at two ends, a base loosening fault occurs at a left support part, and 14 discs are arranged on a shaft;
the shafts are elastic, the mass of each section of shaft is concentrated at the mass center of the corresponding disc, the gyro moment of the disc is not considered, and the equivalent rigidity of the shafts between the turntables is k i I=1..15, the horizontal and vertical displacement of the geometric center of the disk and the center of mass of the support end is denoted as X i ,Y i I=1..16, nonlinear oil film forces in horizontal and vertical directions are denoted as F x ,F y Vibration of loosening failure occurs only in the vertical direction, the corresponding displacement is Y 17 Base k l And base c) l Is a piecewise linear function, expressed as:
the differential equation of motion of the rotor bearing system is expressed as:
wherein X= [ X ] 1 ,...X 16 ,Y 1 ...Y 16 ,Y 17 ]M, C, K are mass, damping and stiffness matrices, respectively, F is the resultant of nonlinear fluid film forces, imbalance forces and gravity of each disk, expressed as:
M＝diag(m 1 ,…,m 16 ,m 1 ,…m 16 ,m 17 )
C＝diag(c 1 ,…,c 16 ,c 1 ,…c 16 ,c 17 )
carrying out dimensionless transformation on a dynamics control equation of the system, wherein the formula (3) is dimensionless by the following formula:
wherein f x And f y Is a dimensionless nonlinear oil film force,is the Sommerfeld number, μ is the lubricant viscosity, c is the oil film thickness, L is the bearing length, R is the bearing radius, ω is the angular frequency, τ is the dimensionless time, and W is half the weight of the rotor system;
thereby establishing a high-dimensional nonlinear rotor system model;
2. inputting the response signals into the snapshot matrix and calculating an autocorrelation matrix;
manifold G (M, N) s ) Is provided withMultiple sampling time sequences N s Covering M-dimensional column vectors of corresponding solutions;
n pre-computed snapshot signal sets of M-dimensional system are set toWherein each parameter omega i Sampling length t at s Same, snapshot set->Seen as manifold G (M, N s ) N points onSelecting a point S of the manifold 0 As a reference point, obtaining a snapshot matrix y with a corresponding parameter omega on the manifold by a Lagrange interpolation method * (ω,t s ) Expressed as:
manifold G (M, N) s ) The interpolation points above are expressed as:
y(ω i ,t s )＝U i Σ i V i t ,i＝0,...,N-1 (7)
wherein Σ is i ∈R M×M Is a matrix of singular values that are used,V i ∈R M×M is an orthogonal matrix, then equation (7) is rewritten as:
let y (omega) i ,t s ) The solution representing the original system with parameter value ω also belongs to manifold G (M, N s ) By thin singularValue decomposition method, y (omega) i ,t s ) Expressed as:
y(ω,t s )＝UΣV T (9)
according to Lagrange interpolation method, the truncation error epsilon N (ω,t s ) The definition is as follows:
wherein the superscript (n+1) denotes the (n+1) -th derivative of the parameter ω, denoted asWhen N → infinity, the truncation error ε N (ω,t s )→0；y(ω i ,t s ) The rewriting is as follows:
in connection with equation (9), equation (11) is expressed as:
diagonal matrix Σ=diag [ σ ] 1 ,σ 2 ,...,σ M ]Is y (omega) i ,t s ) Singular values of (1), wherein sigma 1 ≥σ 2 ≥,...,≥σ M 0 or more, the singular value energy ratio ε is expressed as:
UΣV T is divided into main partsAnd a truncated error fraction epsilon m (ω,t s ) Substituting the main part and the error part into formula (12) to obtain:
wherein ε m (ω,t s ) And epsilon N (ω,t s ) Are minor and omitted, and therefore, equation (14) is rewritten as:
autocorrelation matrix S c The definition is as follows:
the dimensionless control equation of the original system is solved by a fourth-order Runge-Kutta (RK 4) method;
calculating displacement response by inputting different modeling parameters and storing response signals in a snapshot matrixCalculating an autocorrelation matrix S according to equation (16) c The method comprises the following steps:
3. by means of an autocorrelation matrix S c Obtaining the reduced order POM with new parameters
The reduced order POM of the parameter ω is calculated,expressed as:
then equation (17) is rewritten as:
by calculating a correlation matrix S c The feature vector of (2) obtaining reduced order POMs;
4. reduced POM for dynamic differential equation reduction
The reduced POM is defined by the eigenvector of the autocorrelation matrix, sampled at a rotational speed ω' with a sampling length t s Order-makingFeature vectors representing the autocorrelation matrix, the feature vectors being arranged in descending order of the corresponding feature values, provided thatRepresenting a linear subspace having n reduced order POMs;
transient response of complex high-dimensional space is projected onto low-dimensional subspace through coordinate transformation:
wherein, the dimensionless displacement response, the coordinate transformation matrix and the POD modal coordinates of the original system are respectively, and the ordinary differential equation (10) is simplified into:
for the IGM method with the same ROM degree of freedom and k points, singular value decomposition is only required to be executed once, so that the computing efficiency of the IGM method is far higher than that of the ITGM method.
2. The adaptive orthogonal decomposition method for model reduction of a multi-disk rotor system according to claim 1, wherein:
the singular value energy ratio epsilon is epsilon more than or equal to 95 percent.
3. The adaptive orthogonal decomposition method for model reduction of a multi-disk rotor system according to claim 1, wherein:
when the dimensionless control equation is solved by a fourth-order Runge-Kutta method, the initial displacement and the speed are respectively X i ＝Y i ＝Y 17 ＝0.1，The time step is 0.001 pi.
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202110492921.4A CN113254863B (en) | 2021-05-07 | 2021-05-07 | Adaptive orthogonal decomposition method for model reduction of multi-disk rotor system |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202110492921.4A CN113254863B (en) | 2021-05-07 | 2021-05-07 | Adaptive orthogonal decomposition method for model reduction of multi-disk rotor system |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN113254863A CN113254863A (en) | 2021-08-13 |
| CN113254863B true CN113254863B (en) | 2023-11-07 |
Family
ID=77223829
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN202110492921.4A Active CN113254863B (en) | 2021-05-07 | 2021-05-07 | Adaptive orthogonal decomposition method for model reduction of multi-disk rotor system |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN113254863B (en) |
Citations (5)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| DE10137909A1 (en) * | 2001-01-18 | 2002-08-08 | Siemens Ag | Process for the simulation of mechatronic systems |
| DE102013007621A1 (en) * | 2013-05-02 | 2014-11-06 | André Schleicher | Bearing reluctance motor with identical slot pitch in stand and rotor |
| CN111339706A (en) * | 2020-03-09 | 2020-06-26 | 西北工业大学 | POD-based rotor-bearing system model secondary order reduction method |
| CN111368447A (en) * | 2020-03-09 | 2020-07-03 | 西北工业大学 | Improved transient time sequence-based nonlinear POD dimension reduction method |
| CN111563323A (en) * | 2020-04-27 | 2020-08-21 | 西北工业大学 | Transient POD method based on minimum error of bifurcation parameter |
Family Cites Families (1)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US7512525B2 (en) * | 2005-01-05 | 2009-03-31 | Chang Gung University | Multi-point model reductions of VLSI interconnects using the rational Arnoldi method with adaptive orders |
-
2021
- 2021-05-07 CN CN202110492921.4A patent/CN113254863B/en active Active
Patent Citations (5)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| DE10137909A1 (en) * | 2001-01-18 | 2002-08-08 | Siemens Ag | Process for the simulation of mechatronic systems |
| DE102013007621A1 (en) * | 2013-05-02 | 2014-11-06 | André Schleicher | Bearing reluctance motor with identical slot pitch in stand and rotor |
| CN111339706A (en) * | 2020-03-09 | 2020-06-26 | 西北工业大学 | POD-based rotor-bearing system model secondary order reduction method |
| CN111368447A (en) * | 2020-03-09 | 2020-07-03 | 西北工业大学 | Improved transient time sequence-based nonlinear POD dimension reduction method |
| CN111563323A (en) * | 2020-04-27 | 2020-08-21 | 西北工业大学 | Transient POD method based on minimum error of bifurcation parameter |
Non-Patent Citations (2)
| Title |
|---|
| 基于BPOD的气动弹性降阶及其在主动控制中的应用;徐敏;安效民;曾宪昂;陈士橹;;中国科技论文在线(10);全文 * |
| 基于特征正交分解的非定常气动力建模技术;姚伟刚;徐敏;叶茂;;力学学报(04);全文 * |
Also Published As
| Publication number | Publication date |
|---|---|
| CN113254863A (en) | 2021-08-13 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| Jin et al. | Nonlinear dynamic analysis of a complex dual rotor-bearing system based on a novel model reduction method | |
| Harvey III et al. | Steady and unsteady computation of impeller‐stirred reactors | |
| CN110532625A (en) | Aero-engine turbine disk-twin the modeling method of rotor-bearing system number | |
| Hashish et al. | Finite element and modal analyses of rotor-bearing systems under stochastic loading conditions | |
| CN107341288B (en) | Optimization method for controlling vibration of combined cycle unit by adjusting elevation of bearing | |
| CN111339706A (en) | POD-based rotor-bearing system model secondary order reduction method | |
| CN110941184A (en) | Sliding mode vibration active control method for electromagnetic bearing flexible rotor different-position system | |
| Hassan et al. | A new modal-based approach for modelling the bump foil structure in the simultaneous solution of foil-air bearing rotor dynamic problems | |
| Yang et al. | Vibration reduction optimum design of a steam-turbine rotor-bearing system using a hybrid genetic algorithm | |
| Aini et al. | Vibration modeling of rotating spindles supported by lubricated bearings | |
| CN112965512A (en) | Unmanned aerial vehicle wind-resistant control method based on propeller model | |
| CN113254863B (en) | Adaptive orthogonal decomposition method for model reduction of multi-disk rotor system | |
| Chen et al. | Effects of journal static eccentricity on dynamic responses of a rotor system under base motions using FDM inertia model | |
| CN113434983B (en) | Rapid calculation method for nonlinear dynamic characteristics of sliding bearing rotor system | |
| Lu et al. | Nonlinear Dynamic Behavior Analysis of Dual-Rotor-Bearing Systems with Looseness and Rub–Impact Faults | |
| CN111368447A (en) | Improved transient time sequence-based nonlinear POD dimension reduction method | |
| Zhang et al. | Modeling of Robot's Low-Speed Motion Nonlinear Dynamics Based on Phase Space Reconstruction Neural Network | |
| Kim et al. | Steady-state analysis of a nonlinear rotor-housing system | |
| CN111563323A (en) | Transient POD method based on minimum error of bifurcation parameter | |
| Ghalayini et al. | The Application of the Arbitrary-Order Galerkin Reduction Method to the Dynamic Analysis of a Rotor With a Preloaded Single-Pad Foil-Air Bearing | |
| Li et al. | Dynamic characteristics of opposed-conical gas-dynamic bearings | |
| Nandy | Design of nonlinear controllers and speed and angular position estimator for an electrically controlled cycloidal propeller | |
| CN104361145A (en) | Rotor dynamics modeling method based on axle closely attached coordinate system | |
| CN112896555B (en) | Self-balancing control method for rotating speed of attitude control flywheel | |
| CN116090135B (en) | Damper-rotor system response analysis method |
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 |