Numerical Methods and Apparatus for Complex Modal Analysis of Rotor Systems of Rotating Machinery

Abstract

This invention discloses a numerical method, apparatus, equipment, and medium for complex modal analysis of a rotor system of rotating machinery. The method includes: obtaining the differential equation of the free vibration of the rotor-bearing system and determining the corresponding characteristic equation based on the characteristic problem; generating an initial subspace basis matrix when the stiffness matrix in the characteristic equation is a non-singular matrix; constructing a projection problem of the characteristic problem on the subspace by calculating the subspace projection matrix corresponding to the matrix in the characteristic equation; solving the subspace characteristic problem by executing the iterative format of the projection problem; iteratively updating the eigenvalues ​​and eigenvectors of the characteristic equation and normalizing them; and verifying the residuals of the characteristic equation and outputting the solution result of the final characteristic equation if the convergence criterion is satisfied. This enables the direct solution of large-scale asymmetric matrix quadratic generalized characteristic problems, allowing for more efficient complex modal analysis of vibrations in large rotating machinery.

Description

Technical Field

[0001] This invention relates to the field of rotating machinery design, and relates to, but is not limited to, a numerical method and apparatus for complex modal analysis of a rotating machinery rotor system. Background Technology

[0002] Rotating machinery is primarily driven by a rotor system, comprising the rotor, bearings, supports, and other structural components. Steam turbines, generators, centrifuges, and compressors are typical examples of rotating machinery. To ensure the normal operation of rotating machinery, prevent structural resonance, and avoid shaft damage caused by unbalanced mass and eddy currents, it is necessary to calculate the critical speed for structural vibration, stability, and other inherent characteristics during the design process. Therefore, modal analysis of rotating machinery is essential.

[0003] For rotating machinery, the gyroscopic effect caused by the Coriolis force of the rotor must be considered; therefore, the vibration control equations include antisymmetric gyroscopic matrices. The modal space of such systems, which involve damping or asymmetric matrices, is a complex space, and the eigenvalues ​​are usually pairs of conjugate complex numbers. Modal analysis for such systems is called complex modal analysis.

[0004] Typical complex modal analysis methods usually first transform the quadratic generalized characteristic problem into a first-order generalized characteristic problem using state-space methods, and then solve the first-order generalized characteristic problem in state space. Many mature and reliable numerical methods exist for solving the first-order generalized characteristic problem, such as the block Lancsos method, the bioorthogonal Lancsos method, the subspace iteration method, and the synchronous iteration method. However, transforming the quadratic generalized characteristic problem into a first-order generalized characteristic problem doubles the matrix size, thus increasing computer memory usage, computational load, and solution time. Furthermore, the reconstructed matrix loses some of the original matrix's desirable properties, such as decreased matrix sparsity and increased bandwidth, reducing computational and storage efficiency.

[0005] With the general trend of increasing complexity in engineering structures and refinement in mathematical models, the matrix size of the problems to be solved is also growing larger, and traditional solution strategies are gradually becoming insufficient to meet the requirements. For large rotating machinery, in addition to the shaft, dynamic modeling is required for supporting structures such as bearings, seals, and bases, as well as working components such as impellers, blades, and gears on the shaft. Some even involve multi-rotor systems, resulting in complex models, a large workload of numerical calculations, and lengthy computation times.

[0006] Therefore, developing a numerical method for complex modal analysis of large rotating machinery based on a direct solution of the quadratic generalized characteristic problem, which avoids the above-mentioned shortcomings, has important research value and practical significance. Summary of the Invention

[0007] This invention provides a numerical method and apparatus for complex modal analysis of a rotating machinery rotor system, which at least solves the problems of matrix size doubling and matrix property destruction after transforming a quadratic generalized characteristic problem into a first-order generalized characteristic problem in traditional methods.

[0008] The technical solution of this invention is implemented as follows:

[0009] In a first aspect, embodiments of the present invention provide a numerical method for complex modal analysis of a rotor system of rotating machinery, the method comprising:

[0010] The differential equations of the free vibration of the rotor-bearing system are obtained, and the corresponding characteristic equations are determined based on the characteristic problem. The characteristic problem is used to analyze the vibration characteristics of the system. If the stiffness matrix in the characteristic equation is a non-singular matrix, an initial subspace basis matrix is ​​generated. The subspace projection matrix corresponding to the matrix in the characteristic equation is calculated to construct the projection problem of the characteristic problem onto the subspace. The iterative format of the projection problem is executed to solve the subspace characteristic problem. The eigenvalues ​​and eigenvectors of the characteristic equation are iteratively updated and normalized. If the convergence criterion is satisfied, the residuals of the characteristic equation are checked, and the final solution result of the characteristic equation is output.

[0011] In some implementations, determining the corresponding characteristic equation based on the characteristic problem includes: if the characteristic problem is to solve for the m largest eigenvalues, then the characteristic equation is MΦΛ. 2 +CΦΛ+KΦ=0; M is the mass matrix, C is the sum of the damping matrix and the gyroscope matrix, K is the stiffness matrix; Λ is the eigenvalue matrix, which is a diagonal matrix, with each diagonal element representing an eigenvalue to be found; Φ is the eigenvector matrix, with each column of Φ representing a complex eigenvector; if the characteristic problem is to find the minimum m eigenvalues, the characteristic equation is MU+CUT+KUT 2 =0; T is an m-order diagonal matrix, and the diagonal elements t of T s Let λ be the eigenvalue represented by the diagonal elements of Λ. s The reciprocal of λ; U is an n×m matrix, and the s-th column vector of U is the eigenvalue λ. s The corresponding feature vector.

[0012] In some implementations, generating the initial subspace basis matrix includes: generating m orthogonal subspace basis vectors using the Gram-Schmidt method; and combining the m subspace basis vectors column-wise into an n×m matrix U. i , which serves as the basis matrix of the initial subspace.

[0013] In some embodiments, the characteristic equation is MU+CUT+KUT 2When T = 0, the step of constructing the projection problem of the characteristic problem onto the subspace by calculating the subspace projection matrix corresponding to the matrix in the characteristic equation includes: Let T i and U i Substituting the eigenvalue matrix and eigenvector matrix obtained in the i-th iteration into the characteristic equation yields MU. i +CU i T i +KU i T i 2 =R i ; Constructed in subspace U i The orthogonal projection on the surface yields Among them, R i The residual matrix; Let M, C, and K be the projection matrices of M, C, and K, respectively, and Let be the residual matrix of the subspace projection problem and U i The conjugate transpose of .

[0014] In some implementations, performing the iterative scheme for the projection problem to solve the subspace feature problem includes: constructing an inner iterative scheme to obtain matrix H. i And satisfy For sequence H i,j Assign an initial value to H. i,j =T i For the j-th iteration, construct a matrix S of order m. i,j system of equations Solve for S i,j Generalized characteristic problem Among them, P i,j Let D be the eigenvector matrix. i,j Let S be the eigenvalue matrix. i,j It is an m×m matrix.

[0015] In some implementations, the step of exiting the outer iteration when the convergence criterion is met, checking the residual of the characteristic equation, and outputting the final solution result of the characteristic equation includes: calculating the difference between the results of two adjacent iterations, and determining that the convergence criterion is met when the difference between the results of two adjacent iterations is less than a specified error, and exiting the outer iteration; substituting the k-th set of eigenvalues ​​and the corresponding eigenvector into the characteristic equation to calculate the residual vector; k is a natural number between 1 and m; if the magnitude of the residual vector is less than the error threshold, determining that the calculation result of the k-th set of eigenvalues ​​and the corresponding eigenvector is correct and outputting it as the final solution result.

[0016] Secondly, embodiments of the present invention provide a numerical apparatus for complex modal analysis of a rotating machinery rotor system, comprising an equation determination module, a subspace generation module, a matrix projection module, a subspace solving module, an iterative update module, and a verification output module, wherein:

[0017] The equation determination module is used to obtain the differential equation of the free vibration of the rotor-bearing system and determine the corresponding characteristic equation based on the characteristic problem; wherein, the characteristic problem is used to analyze the vibration characteristics of the system.

[0018] The subspace generation module is used to generate an initial subspace basis matrix when the stiffness matrix in the characteristic equation is a non-singular matrix.

[0019] The matrix projection module is used to construct the projection problem of the feature problem onto the subspace by calculating the subspace projection matrix corresponding to the matrix in the feature equation.

[0020] The subspace sub-solution module is used to execute the iterative format of the projection problem and solve the subspace feature problem;

[0021] The iterative update module is used to iteratively update the eigenvalues ​​and eigenvectors of the characteristic equation and normalize them;

[0022] The verification output module is used to verify the residuals of the characteristic equation and output the final solution result of the characteristic equation, provided that the convergence criterion is met.

[0023] Thirdly, embodiments of the present invention provide an electronic device, including a memory and a processor, wherein the memory stores a computer program that can run on the processor, and the processor executes the program to implement the steps in the above-described method.

[0024] Fourthly, embodiments of the present invention provide a computer-readable storage medium having a computer program stored thereon, which, when executed by a processor, implements the steps in the above-described method.

[0025] The beneficial effects of the technical solutions provided in the embodiments of the present invention include at least the following:

[0026] In this embodiment of the invention, for the problem of solving the characteristic equation of a rotor-bearing system, an inner-layer iterative scheme is added to the traditional subspace iterative method. This extends the projection method, which is applicable to first-order generalized problems, to directly solve second-order characteristic problems. This enables the direct solution of large-scale asymmetric matrix second-order generalized characteristic problems, avoiding the disadvantages of traditional methods such as doubling the matrix size and destroying matrix properties. This allows for more efficient complex modal analysis of vibration of large rotating machinery. Attached Figure Description

[0027] To more clearly illustrate the technical solutions in the embodiments of the present invention, the accompanying drawings used in the description of the embodiments will be briefly introduced below. Obviously, the drawings described below are only some embodiments of the present invention. For those skilled in the art, other drawings can be obtained based on these drawings without creative effort, wherein:

[0028] Figure 1 A flowchart illustrating the numerical method for complex modal analysis of a rotating machinery rotor system provided in an embodiment of the present invention;

[0029] Figure 2 This is a schematic diagram of the rotor-bearing system provided in an embodiment of the present invention;

[0030] Figure 3 A flowchart illustrating the numerical method for complex modal analysis of a rotating machinery rotor system provided in an embodiment of the present invention;

[0031] Figure 4 The first three modal frequency convergence curves of the rotor-bearing system provided in this embodiment of the invention;

[0032] Figure 5 A schematic diagram comparing the convergence curves of the traditional subspace iteration method and the method of the present invention, provided in an embodiment of the present invention.

[0033] Figure 6 A schematic diagram of the composition structure of a numerical device for complex modal analysis of a rotor system of rotating machinery provided in an embodiment of the present invention;

[0034] Figure 7 This is a schematic diagram of the hardware entity of an electronic device provided in an embodiment of the present invention. Detailed Implementation

[0035] To make the objectives, technical solutions, and advantages of the embodiments of the present invention clearer, the technical solutions of the embodiments of the present invention will be clearly and completely described below with reference to the accompanying drawings. Obviously, the described embodiments are only some embodiments of the present invention, not all embodiments. The following embodiments are used to illustrate the present invention, but are not intended to limit the scope of the present invention. Based on the embodiments of the present invention, all other embodiments obtained by those skilled in the art without creative effort are within the scope of protection of the present invention.

[0036] In the following description, references are made to “some embodiments,” which describe a subset of all possible embodiments. However, it is understood that “some embodiments” may be the same subset or different subsets of all possible embodiments and may be combined with each other without conflict.

[0037] It should be noted that the terms "first\second\third" involved in the embodiments of the present invention are only used to distinguish similar objects, and do not represent a specific ordering of the objects. It can be understood that "first\second\third" may exchange the specific order or precedence where permitted, so that the embodiments of the present invention described herein can be implemented in an order other than those illustrated or described herein.

[0038] Those skilled in the art can understand that unless otherwise defined, all terms (including technical terms and scientific terms) used herein have the same meaning as commonly understood by those of ordinary skill in the art to which the embodiments of the present invention pertain. It should also be understood that terms such as those defined in general dictionaries should be understood to have meanings consistent with their meanings in the context of the prior art, and unless specifically defined as herein, they will not be interpreted with idealized or overly formal meanings.

[0039] Modal analysis is the basis of dynamic analysis. It transforms the vibration differential equation from physical coordinates to modal coordinates, and obtains a set of decoupled equations described by modal parameters, thereby identifying the natural vibration characteristics of the structure such as frequency, damping ratio and mode shape, which provides an important basis for vibration characteristic analysis, fault diagnosis and optimization design of engineering structures.

[0040] The key technology of modal analysis is to solve the eigenvalues and eigenvectors of the characteristic equation, that is, the solution of the generalized eigenvalue problem. In practical engineering problems, there are various non-conservative systems, that is, systems with energy dissipation during vibration, such as damped systems.

[0041] The key technology of complex modal analysis is the solution of the quadratic generalized eigenvalue problem, that is, in the form of MΦΛ 2 +CΦΛ+KΦ=0, which is the solution of the characteristic equation. Wherein, M is the mass matrix, C is used herein to represent the sum of the damping matrix and the antisymmetric gyro matrix, and K is the stiffness matrix. Λ is the eigenvalue matrix, which is a diagonal matrix, and each diagonal element represents a complex eigenvalue; Φ is the eigenvector matrix, and each column of Φ represents a complex eigenvector.

[0042] The subspace iteration method is a method that combines iteration method and Rayleigh-Ritz method and uses them alternately. It not only uses Rayleigh-Ritz method to reduce degrees of freedom, but also uses iteration method to make the mode shapes gradually approach the required accuracy in calculation. In the subspace iteration method, n trial vectors are selected first (n<N, N is the total degrees of freedom of the system), and iterations are performed on these n vectors at the same time. Generally, the degrees of freedom of a structure are in thousands or tens of thousands, while only dozens of mode shapes need to be solved. The subspace iteration method does not require global solution, but stops when obtaining the required result. The subspace iteration method uses the vectors obtained by the iteration method as Rayleigh-Ritz method vectors, and then uses the Rayleigh-Ritz method to solve n approximate eigenpairs, which is reduced to solving the eigenpair problem in the degenerate subspace.

[0043] This invention provides a numerical method for complex modal analysis of a rotating machinery rotor system, applicable to electronic devices. These electronic devices include, but are not limited to, mobile phones, laptops, tablets, handheld internet devices, multimedia devices, streaming media devices, mobile internet devices, wearable devices, or other types of electronic devices. The functions implemented by this method can be achieved by a processor in the electronic device calling program code. The program code can be stored in a computer storage medium; therefore, the electronic device includes at least a processor and a storage medium. The processor can be used to process the complex modal analysis numerical process of the rotating machinery rotor system, and the memory can be used to store the data required and generated during the complex modal analysis numerical process of the rotating machinery rotor system.

[0044] Figure 1 A flowchart illustrating a numerical method for complex modal analysis of a rotor system of rotating machinery provided in an embodiment of the present invention is shown below. Figure 1 As shown, the method includes at least the following steps:

[0045] Step S110: Obtain the differential equation of the free vibration of the rotor-bearing system, and determine the corresponding characteristic equation based on the characteristic problem.

[0046] Here, the rotor-bearing system refers to a system in large rotating machinery such as marine propulsion shafts, compressors, pumps, etc. The differential equations of the free vibration of this system can be obtained through methods such as the finite element method or the transfer matrix method.

[0047] The characteristic problem is used to analyze the vibration characteristics of the system, such as identifying the system's natural frequencies, damping ratios, and configuration. Since the bearing stiffness matrix is ​​asymmetric and the gyroscope matrix is ​​antisymmetric, the differential equation for the free vibration of the rotor-bearing system is a linear ordinary differential equation for a steady system. The corresponding characteristic problem is a quadratic generalized characteristic problem, and the modal analysis of this rotor-bearing system belongs to complex modal analysis. In some embodiments, the characteristic problem is to solve for several higher-order modes or several lower-order modes of the rotor-bearing system.

[0048] In some implementations, if the characteristic problem is to find the largest m eigenvalues, the characteristic equation is:

[0049] MΦΛ 2 +CΦΛ+KΦ=0 formula (1);

[0050] Where M is the mass matrix, C is the sum of the damping matrix and the gyroscope matrix, K is the stiffness matrix; Λ is the eigenvalue matrix, which is a diagonal matrix, and each diagonal element represents an eigenvalue to be determined; Φ is the eigenvector matrix, and each column of Φ represents a complex eigenvector.

[0051] In some implementations, if the characteristic problem is to find the minimum m eigenvalues, the characteristic equation is:

[0052] MU+CUT+KUT 2 =0 formula (2);

[0053] Where T is an m-order diagonal matrix, and the diagonal elements t of T s Let λ be the eigenvalue represented by the diagonal elements of Λ. s The reciprocal of λ; U is an n×m matrix, and the s-th column vector of U is the eigenvalue λ. s The corresponding feature vector.

[0054] Step S120: If the stiffness matrix in the characteristic equation is a non-singular matrix, generate the initial subspace basis matrix.

[0055] Here, for practical engineering problems, constraints are typically introduced to ensure that the stiffness matrix is ​​non-singular. A vector space can contain many subspaces. Subspace basis vectors define a certain m-dimensional subspace of an n-dimensional complex space. A total of m subspace basis vectors are combined column-wise to form an n×m matrix U. i A set of orthogonal basis vectors can be generated using the Gram-Schmidt method.

[0056] In some implementations, the initial subspace basis matrix is ​​generated through the following process: generating m orthogonal subspace basis vectors using the Gram-Schmidt method; and combining the m subspace basis vectors column-wise into an n×m matrix U. i , which serves as the basis matrix of the initial subspace.

[0057] Step S130: By calculating the subspace projection matrix corresponding to the matrix in the characteristic equation, the projection problem of the characteristic problem on the subspace is constructed.

[0058] It should be noted that, in practical engineering, it is usually necessary to solve for several low-order modes. Subsequent solutions focus on solving for the smallest m eigenvalues, i.e., the characteristic equation is MU + CUT + KUT. 2 The calculation process is described using 0 as an example.

[0059] In some embodiments, the characteristic equation is MU+CUT+KUT 2 When T = 0, the above step S130 is further implemented through the following process: Let T i and U i Substituting the eigenvalue matrix and eigenvector matrix obtained in the i-th iteration into the characteristic equation yields:

[0060] MUi +CU i T i +KU i T i 2 =R i Formula (3);

[0061] Constructed in subspace U i The orthogonal projection on the surface yields:

[0062]

[0063] Among them, R i The residual matrix; Let M, C, and K be the projection matrices of M, C, and K, respectively, and Let be the residual matrix of the subspace projection problem and U i The conjugate transpose of .

[0064] because Both are m*m square matrices, U i It is a matrix composed of column-wise combinations of m subspace basis vectors. As the iteration proceeds, these m subspace basis vectors gradually converge to the first m eigenvectors.

[0065] Step S140: Execute the iterative format of the projection problem to solve the subspace feature problem.

[0066] Here, the quadratic generalized eigenvalues ​​and eigenvectors of the projection problem are solved through the following procedure:

[0067] Construct the inner iterative scheme to obtain matrix H i And satisfy:

[0068]

[0069] For sequence H i,j Assign an initial value to H. i,j =T i ;

[0070] For the j-th iteration, construct the m-order matrix S. i,j The system of equations:

[0071]

[0072] Solve for S i,j The problem of generalized characteristics:

[0073]

[0074] Among them, P i,jLet D be the eigenvector matrix. i,j Let S be the eigenvalue matrix. i,j It is an m×m matrix.

[0075] It is worth noting that S i,j Since it is an m×m matrix, the solution process belongs to the generalized characteristic problem of small-scale matrices, which will not take up a long computation time or a lot of memory.

[0076] Because of D i,j Yes An approximation for determining H i,j+1 To approximate the result, we need to perform square root extraction:

[0077]

[0078] In the formula, d kk It is D i,j The kth diagonal element, the sequence {n k The values ​​of each term in} are either 0 or 1, used to determine the sign of the square root. The sign can be obtained through trial and error. and Substituting each into the characteristic equation of the projection problem, we obtain the residuals r. + and r - :

[0079]

[0080]

[0081] If r + <r - Then take the positive sign; let n k = 0; otherwise, take the negative sign and let n k =1. This process involves only a few small matrix multiplications, so it does not consume much computational resources.

[0082] Update H i,j+1 The estimated value H i,j+1 =G i,j If the result of the Nth iteration satisfies the convergence criterion, then let H i =H i,N P i =P i,N .

[0083] The above outlines the main steps of the inner layer iteration. As the outer layer iteration proceeds, the eigenvalues ​​and eigenvectors gradually converge, reducing the number of iterations required for the inner layer until it is no longer necessary.

[0084] Step S150: Iteratively update the eigenvalues ​​and eigenvectors of the characteristic equation and normalize them.

[0085] Here, the new estimate T of the updated eigenvalues ​​is... i+1 =H i , solve for KU i+1 ′=-(MU i P i +CU i T i+1 P i The new estimates of the eigenvectors are updated, and finally the eigenvectors are normalized.

[0086] Step S160: If the convergence criterion is met, check the residuals of the characteristic equation and output the final solution result of the characteristic equation.

[0087] Here, the convergence criterion indicates that the difference between the results of two adjacent iterations is within a specified error range. The solved eigenvalues ​​and their corresponding eigenvectors are fed into the original characteristic equation to calculate the residual of the characteristic equation, i.e., the residual vector, and then checked whether it is less than the error threshold.

[0088] In some implementations, step S160 is performed by the following process: calculating the difference between the results of two adjacent iterations, and determining that the convergence criterion is met when the difference between the results of two adjacent iterations is less than a specified error, and exiting the outer iteration; substituting the k-th set of eigenvalues ​​and the corresponding eigenvector into the characteristic equation to calculate the residual vector; k is a natural number between 1 and m; if the magnitude of the residual vector is less than the error threshold, determining that the calculation result of the k-th set of eigenvalues ​​and the corresponding eigenvector is correct and outputting it as the final solution result.

[0089] In this embodiment of the invention, for the problem of solving the characteristic equation of a rotor-bearing system, an inner-layer iterative scheme is added to the traditional subspace iterative method. This extends the projection method, which is applicable to first-order generalized problems, to directly solve second-order characteristic problems. This enables the direct solution of large-scale asymmetric matrix second-order generalized characteristic problems, avoiding the disadvantages of traditional methods such as doubling the matrix size and destroying matrix properties. This allows for more efficient complex modal analysis of vibration of large rotating machinery.

[0090] The numerical method for complex modal analysis of the above-mentioned rotating machinery rotor system is described below with specific examples. However, it is worth noting that the specific embodiments are only for better illustrating the present invention and do not constitute an improper limitation of the present invention.

[0091] This invention provides a numerical method for complex modal analysis of large rotating machinery, aiming to perform vibration complex modal analysis of large rotating machinery more efficiently. This method eliminates the need for reassembling the matrix using the state-space method. By adding inner iteration steps to the traditional subspace iterative method, it achieves direct solution to the quadratic generalized eigenvalue problem of large-scale asymmetric matrices, and the number of inner iterations gradually decreases until it is no longer necessary. This avoids the drawbacks of traditional methods, such as doubling the matrix size and violating matrix properties.

[0092] Example 1 illustrates this by taking the calculation of the first three natural frequencies of a rotor-bearing system as an example. The rotor-bearing system model is as follows: Figure 2 As shown, since the bearing stiffness matrix is ​​asymmetric and the gyroscope matrix is ​​antisymmetric, the modal analysis of this rotor-bearing system belongs to complex modal analysis.

[0093] The differential equations for the free vibration of this system can be obtained through methods such as the finite element method or the transfer matrix method. Where M is the mass matrix, C is the sum of the damping matrix and the gyroscope matrix, and K is the stiffness matrix. All of these matrices are positive definite, sparse, 1000-order square matrices. For practical engineering problems, constraints can usually be introduced to ensure that the stiffness matrix K is a non-singular matrix.

[0094] This is a linear ordinary differential equation for a time-invariant system, and its general solution can be obtained by solving the following characteristic equation: MΦΛ 2 +CΦΛ+KΦ=0; where Λ is the eigenvalue matrix, which is a 6th order diagonal matrix, and each diagonal element represents an eigenvalue to be found; Φ is the eigenvector matrix with 1000×6 rows and columns, and each column represents an eigenvector.

[0095] The three smallest third-order natural frequencies essentially correspond to the imaginary parts of the first three pairs of conjugate eigenvalues, i.e., the first six eigenvalues ​​with the smallest imaginary parts. For practical engineering structures, the real part is usually much smaller than the imaginary part, and the six eigenvalues ​​with the smallest modulus can generally be solved.

[0096] First, perform a substitution on the characteristic equation. Let T be a 6th-order diagonal matrix, whose diagonal elements are the reciprocals of the diagonal elements of Λ. Let U be a 1000×6 matrix, whose s-th column vector is the eigenvalue λ. s The corresponding eigenvectors. The characteristic equation after substitution is MU + CUT + KUT. 2 =0.

[0097] This is a typical quadratic generalized eigenvalue problem. The complex modal analysis numerical method provided in this invention does not perform order reduction processing on this equation, but directly solves this quadratic generalized eigenvalue problem.

[0098] like Figure 3The diagram shown is a logical flowchart of the numerical method for complex modal analysis of a rotating machinery rotor system provided in an embodiment of the present invention, which includes the following steps:

[0099] S1, input stiffness matrix K, mass matrix M and damping matrix C.

[0100] S2, determine whether K is non-singular.

[0101] Here, if the determination result is yes, continue to execute step S3; otherwise, if the determination result is no, end the process.

[0102] S3 generates the initial subspace basis vectors.

[0103] Generate six orthogonal subspace basis vectors using the Gram-Schmidt method. Specify or randomly generate six linearly independent 1000-row complex column vectors u1 to u6. Execute the outer loop, letting k = 1, 2, ..., 5. In the k-th step, perform the following procedure: Among them, u k+1 Let u be the basis vector of the (k+1)th subspace (to be orthogonalized); s Let be the basis vector of the s-th subspace (which has been orthogonalized), where s = 1, 2, ..., k.

[0104] Then, the vectors are normalized. The six vectors are arranged column-wise into a 1000×6 matrix U1 as the basis matrix of the initial subspace, which defines a 6-dimensional subspace in the 1000-dimensional complex space.

[0105] S4, construct the projection problem of the original feature problem onto the subspace.

[0106] Remember T i and U i Substituting the eigenvalue matrix and eigenvector matrix obtained in the i-th iteration into the characteristic equation yields... Where R i Let be the residual matrix.

[0107] Calculate the subspace projection matrix of the original M, C, and K matrices, and construct it in subspace U. i The orthogonal projection on the surface yields in,

[0108] S5, solve for the quadratic generalized eigenvalues ​​and eigenvectors of the projection problem.

[0109] Construct the inner iterative scheme to obtain matrix H i It satisfies: For sequence H i,j Assign an initial value to H. i,j =T iFor the j-th iteration, construct a matrix S of order 6. i,j The system of equations (m=6) is Then solve for S. i,j Generalized characteristic problem Where P i,j Let D be the eigenvector matrix. i,j Let S be the eigenvalue matrix. Note that S... i,j It is a 6×6 matrix.

[0110] S6, the square root sum determination symbol.

[0111] Because of D i,j Yes An approximation for determining H i,j+1 To approximate the result, we need to perform square root extraction: Where m = 6, d kk It is D i,j The kth diagonal element, the sequence {n k The values ​​of each term in} are 0 or 1, which are used to determine the sign of the square root.

[0112] Symbols can be obtained through trial and error. and Substitute each value into the characteristic equation of the projection problem and choose the one with the smaller residual. This process involves only a few multiplications of small matrices, and therefore does not consume many computational resources.

[0113] S7 determines that the residual of the projection problem is sufficiently small.

[0114] Using formula H i,j+1 =G i,j Update H i,j+1 The estimated value of this iterative scheme can have a relatively relaxed convergence criterion, and ||H| can be used. i,j+1 -H i,j | m1 ≤1×10 -2 ||H i,1 || m1 If the result of the Nth iteration satisfies the convergence criterion, i.e., the result is yes, then let H... i =H i,N P i =P i,N If not, proceed to step S8; otherwise, it indicates that the residual of the projection problem is not sufficiently small, the judgment result is negative, and it is necessary to return to step S5.

[0115] Steps S5 and S7 above are the main steps of the inner layer iteration. As the outer layer iteration proceeds, the eigenvalues ​​and eigenvectors gradually converge, and the number of iterations required for the inner layer iteration will decrease until no further inner layer iteration is needed.

[0116] S8 updates the eigenvalues ​​and eigenvectors and normalizes them.

[0117] First, use formula T i+1 =H i Update the new estimates of the eigenvalues; then solve the following equation to update the new estimates of the eigenvectors: KU i+1 ′=-(MU i P i +CU i T i+1 P i Finally, the eigenvectors are normalized using the following formula:

[0118] S9 determines whether the convergence has occurred.

[0119] Calculate the difference between the results of two adjacent iterations. The difference between the results of calculating the k-th eigenvalue is expressed as... The difference between the calculated results of the k-th eigenvalue vector is expressed as:

[0120] If the difference between the results of two adjacent iterations is less than the specified error, the iteration is considered to have converged, and the iteration is exited and step S10 is executed; otherwise, the iteration process is returned to step S2 and the iteration process continues.

[0121] S10, check whether the residual of the characteristic equation is sufficiently small.

[0122] Since the characteristic equation involves substitution, the eigenvalues ​​need to be taken in reverse order, while the eigenvectors remain unchanged. Substituting the k-th set of eigenvalues ​​and their corresponding eigenvectors into the characteristic equation yields the residual vector R. (k) =Mu (k) (t (k) ) 2 +Cu (k) t (k) +Ku (k) Determine ||R (k) The value of ||2 is less than the specified error to check whether the residual of the characteristic equation is sufficiently small.

[0123] S11 outputs the eigenvalues ​​and their corresponding eigenvectors.

[0124] If ||R (k) If ||2 is less than the specified error, then the calculation results of the k-th eigenvalue and the corresponding eigenvector can be considered accurate and can be output as the solution result.

[0125] The calculations were performed using both the traditional subspace iteration method and the numerical method for complex modal analysis proposed in this invention, and the results are shown in Table 1. The convergence curves of the first six modal frequencies (imaginary eigenvalues) are shown below. Figure 4 As shown.

[0126] Table 1 Comparison of calculation results from the two methods

[0127]

[0128] The numerical method for complex modal analysis provided by this invention converges to sufficient accuracy after 11 iterations. Comparative verification shows that the calculation results of the two algorithms are in agreement. This example demonstrates that the numerical method for complex modal analysis provided by this invention can accurately solve typical asymmetric quadratic generalized eigenvalue problems and is applicable to typical modal analysis problems of rotating machinery.

[0129] Example 2 calculates the first six modal frequencies of a dual-rotor system of a certain type of aero-engine. The computational model in this example is more complex than that in Example 1, and the algorithm of this invention is compared with the traditional subspace iterative algorithm. A quadratic characteristic problem solver was developed based on the complex modal analysis numerical method provided in this embodiment of the invention, and this solver was integrated into a self-developed rotor dynamics analysis software. The computational performance of the solver was verified through a series of examples. The computation time of the subspace iterative method and the method of this invention was compared, demonstrating the superiority of the complex modal analysis numerical method proposed in this invention. The calculation results are shown in Table 2, and the convergence curves of the subspace iterative method and the complex modal analysis numerical method provided by this invention are compared. Figure 5 As shown, the comparison of calculation results is shown in Table 3.

[0130] Table 2 Comparison of Calculation Results from Two Methods

[0131]

[0132] Table 3 Comparison of Calculation Results of Two Methods

[0133]

[0134] The above results demonstrate that the calculation results from the two methods are in good agreement. The method of this invention, namely the numerical method for complex modal analysis, not only yields accurate results but also boasts faster calculation speed and higher efficiency.

[0135] This invention provides a numerical method for directly solving the quadratic generalized eigenvalue problem of large-scale asymmetric matrices. By adding an inner-layer iterative step to the traditional subspace iterative method, the direct solution to the quadratic generalized eigenvalue problem of large-scale asymmetric matrices is achieved. Compared with traditional methods, this method has the following advantages: 1) lower computational cost and shorter time required to solve the quadratic generalized eigenvalue problem; 2) smaller matrix size and lower memory usage; 3) preservation of the original matrix structure, retaining its sparsity, narrow bandwidth, and other excellent properties; 4) applicability to eigenvalue problems of asymmetric matrices, requiring only that the matrix be nonsingular.

[0136] Based on the foregoing embodiments, this invention further provides a numerical device for complex modal analysis of a rotating machinery rotor system. The device includes all the modules included, which can be implemented by a processor in an electronic device; of course, it can also be implemented by specific logic circuits. In the implementation process, the processor can be a central processing unit (CPU), a microprocessor (MPU), a digital signal processor (DSP), or a field programmable gate array (FPGA), etc.

[0137] Figure 6 This is a schematic diagram of the composition structure of a numerical device for complex modal analysis of a rotor system of rotating machinery provided in an embodiment of the present invention, as shown below. Figure 6 As shown, the device 600 includes an equation determination module 610, a subspace generation module 620, a matrix projection module 630, a subspace solving module 640, an iterative update module 650, and a verification output module 660, wherein:

[0138] The equation determination module 610 is used to obtain the differential equation of the free vibration of the rotor-bearing system and determine the corresponding characteristic equation based on the characteristic problem; wherein, the characteristic problem is used to analyze the vibration characteristics of the system.

[0139] The subspace generation module 620 is used to generate an initial subspace basis matrix when the stiffness matrix in the characteristic equation is a non-singular matrix.

[0140] The matrix projection module 630 is used to construct the projection problem of the feature problem on the subspace by calculating the subspace projection matrix corresponding to the matrix in the feature equation.

[0141] The subspace sub-solution module 640 is used to execute the iterative format of the projection problem and solve the subspace feature problem;

[0142] The iterative update module 650 is used to iteratively update the eigenvalues ​​and eigenvectors of the characteristic equation and normalize them.

[0143] The verification output module 660 is used to verify the residual of the characteristic equation and output the final solution result of the characteristic equation when the convergence criterion is met.

[0144] In some possible embodiments, the equation determination module 610 is specifically used to determine the characteristic equation if the characteristic problem is to find the largest m eigenvalues. 2+CΦΛ+KΦ=0; M is the mass matrix, C is the sum of the damping matrix and the gyroscope matrix, K is the stiffness matrix; Λ is the eigenvalue matrix, which is a diagonal matrix, with each diagonal element representing an eigenvalue to be found; Φ is the eigenvector matrix, with each column of Φ representing a complex eigenvector; if the characteristic problem is to find the minimum m eigenvalues, the characteristic equation is MU+CUT+KUT 2 =0; T is an m-order diagonal matrix, and the diagonal elements t of T s Let λ be the eigenvalue represented by the diagonal elements of Λ. s The reciprocal of λ; U is an n×m matrix, and the s-th column vector of U is the eigenvalue λ. s The corresponding feature vector.

[0145] In some possible embodiments, the subspace generation module 620 is specifically used to generate m orthogonal subspace basis vectors using the Gram-Schmidt method; and to combine the m subspace basis vectors column-wise into an n×m matrix U. i , which serves as the basis matrix of the initial subspace.

[0146] In some possible embodiments, the characteristic equation is MU+CUT+KUT 2 When T = 0, the matrix projection module 630 is specifically used to record T. i and U i Substituting the eigenvalue matrix and eigenvector matrix obtained in the i-th iteration into the characteristic equation yields MU. i +CU i T i +KU i T i 2 =R i ; Constructed in subspace U i The orthogonal projection on the surface yields Among them, R i The residual matrix; Let M, C, and K be the projection matrices of M, C, and K, respectively, and Let be the residual matrix of the subspace projection problem and U i The conjugate transpose of .

[0147] In some possible embodiments, the subspace sub-solver module 640 constructs an inner-layer iterative scheme to obtain matrix H. i And satisfy For sequence H i,j Assign an initial value to H. i,j =T i For the j-th iteration, construct a matrix S of order m. i,jsystem of equations Solve for S i,j Generalized characteristic problem Among them, P i,j Let D be the eigenvector matrix. i,j Let S be the eigenvalue matrix. i,j It is an m×m matrix.

[0148] In some possible embodiments, the verification output module 660 is specifically used to calculate the difference between the results of two adjacent iterations, and when the difference between the results of two adjacent iterations is less than a specified error, it is determined that the convergence criterion is met and the outer iteration is skipped; the k-th set of eigenvalues ​​and the corresponding eigenvectors are substituted into the characteristic equation to calculate the residual vector; k is a natural number between 1 and m; when the magnitude of the residual vector is less than the error threshold, the calculation result of the k-th set of eigenvalues ​​and the corresponding eigenvectors is determined to be correct and is output as the final solution result.

[0149] It should be noted that the descriptions of the above device embodiments are similar to those of the above method embodiments, and have similar beneficial effects. For technical details not disclosed in the device embodiments of the present invention, please refer to the descriptions of the method embodiments of the present invention for understanding.

[0150] It should be noted that, in the embodiments of the present invention, if the above-mentioned numerical method for complex modal analysis of the rotating machinery rotor system is implemented in the form of a software functional module and sold or used as an independent product, it can also be stored in a computer-readable storage medium. Based on this understanding, the technical solution of the embodiments of the present invention, or the part that contributes to related technologies, can be embodied in the form of a software product. This computer software product is stored in a storage medium and includes several instructions to cause an electronic device to execute all or part of the methods described in the various embodiments of the present invention. The aforementioned storage medium includes various media capable of storing program code, such as USB flash drives, portable hard drives, read-only memory (ROM), magnetic disks, or optical disks. Thus, the embodiments of the present invention are not limited to any specific hardware and software combination.

[0151] Correspondingly, embodiments of the present invention provide a computer-readable storage medium storing a computer program thereon. When executed by a processor, the computer program implements the steps in the numerical method for complex modal analysis of any of the rotating machinery rotor systems described in the above embodiments. Correspondingly, embodiments of the present invention also provide a computer program product, which, when executed by a processor of an electronic device, is used to implement the steps in the numerical method for complex modal analysis of any of the rotating machinery rotor systems described in the above embodiments.

[0152] Based on the same technical concept, embodiments of the present invention provide an electronic device for implementing the numerical method for complex modal analysis of a rotating machinery rotor system described in the above method embodiments. Figure 7 This is a hardware entity diagram of an electronic device provided in an embodiment of the present invention, such as... Figure 7 As shown, the electronic device 700 includes a memory 710 and a processor 720. The memory 710 stores a computer program that can run on the processor 720. When the processor 720 executes the program, it implements the steps in the numerical method for complex modal analysis of the rotating machinery rotor system according to any embodiment of the present invention.

[0153] The memory 710 is configured to store instructions and applications executable by the processor 720, and can also cache data to be processed or already processed by the processor 720 and various modules in the electronic device (e.g., image data, audio data, voice communication data and video communication data), which can be implemented by flash memory or random access memory (RAM).

[0154] When the processor 720 executes the program, it implements the steps of the numerical method for complex modal analysis of the rotating machinery rotor system described above. The processor 720 typically controls the overall operation of the electronic equipment 700.

[0155] The aforementioned processor can be at least one of the following: Application Specific Integrated Circuit (ASIC), Digital Signal Processor (DSP), Digital Signal Processing Device (DSPD), Programmable Logic Device (PLD), Field Programmable Gate Array (FPGA), Central Processing Unit (CPU), Controller, Microcontroller, and Microprocessor. It is understood that other electronic devices can also implement the functions of the aforementioned processor, and this embodiment of the invention does not specifically limit their use.

[0156] The aforementioned computer storage media / memory can be read-only memory (ROM), programmable read-only memory (PROM), erasable programmable read-only memory (EPROM), electrically erasable programmable read-only memory (EEPROM), magnetic random access memory (FRAM), flash memory, magnetic surface memory, optical disc, or compact disc read-only memory (CD-ROM), etc.; or it can be various electronic devices that include one or any combination of the above-mentioned memories, such as mobile phones, computers, tablet devices, personal digital assistants, etc.

[0157] It should be noted that the descriptions of the storage medium and device embodiments above are similar to those of the method embodiments above, and have similar beneficial effects. For technical details not disclosed in the storage medium and device embodiments of the present invention, please refer to the descriptions of the method embodiments of the present invention for understanding.

[0158] It should be understood that the phrase "one embodiment" or "an embodiment" throughout the specification means that a specific feature, structure, or characteristic related to the embodiment is included in at least one embodiment of the invention. Therefore, "in one embodiment" or "in an embodiment" appearing throughout the specification does not necessarily refer to the same embodiment. Furthermore, these specific features, structures, or characteristics can be combined in any suitable manner in one or more embodiments. It should be understood that in the various embodiments of the invention, the sequence numbers of the above-described processes do not imply a sequential order of execution; the execution order of each process should be determined by its function and internal logic, and should not constitute any limitation on the implementation process of the embodiments of the invention. The sequence numbers of the above-described embodiments of the invention are merely descriptive and do not represent the superiority or inferiority of the embodiments.

[0159] It should be noted that, in this document, the terms "comprising," "including," or any other variations thereof are intended to cover non-exclusive inclusion, such that a process, method, article, or apparatus that comprises a list of elements includes not only those elements but also other elements not expressly listed, or elements inherent to such a process, method, article, or apparatus. Unless otherwise specified, an element defined by the phrase "comprising one..." does not exclude the presence of other identical elements in the process, method, article, or apparatus that includes that element.

[0160] In the several embodiments provided by this invention, it should be understood that the disclosed devices and methods can be implemented in other ways. The device embodiments described above are merely illustrative. For example, the division of units is only a logical functional division, and in actual implementation, there may be other division methods, such as: multiple units or components can be combined, or integrated into another system, or some features can be ignored or not executed. In addition, the coupling, direct coupling, or communication connection between the various components shown or discussed can be through some interfaces, and the indirect coupling or communication connection between devices or units can be electrical, mechanical, or other forms.

[0161] The units described above as separate components may or may not be physically separate. The components shown as units may or may not be physical units. They may be located in one place or distributed across multiple network units. Some or all of the units may be selected to achieve the purpose of the embodiments of the present invention, depending on actual needs.

[0162] Furthermore, in the various embodiments of the present invention, all functional units can be integrated into one processing unit, or each unit can be a separate unit, or two or more units can be integrated into one unit. The integrated unit can be implemented in hardware or in a combination of hardware and software functional units. Alternatively, if the integrated unit of the present invention is implemented as a software functional module and sold or used as an independent product, it can be stored in a computer-readable storage medium. Based on this understanding, the technical solution of the embodiments of the present invention, or the part that contributes to related technologies, can be embodied in the form of a software product. This computer software product is stored in a storage medium and includes several instructions to cause the device automatic test line to execute all or part of the methods described in the various embodiments of the present invention. The aforementioned storage medium includes various media capable of storing program code, such as mobile storage devices, ROMs, magnetic disks, or optical disks.

[0163] The methods disclosed in the several method embodiments provided by this invention can be arbitrarily combined without conflict to obtain new method embodiments. The features disclosed in the several method or device embodiments provided by this invention can be arbitrarily combined without conflict to obtain new method or device embodiments.

[0164] The above description is merely an embodiment of the present invention, but the scope of protection of the present invention is not limited thereto. Any variations or substitutions that can be easily conceived by those skilled in the art within the technical scope disclosed in the present invention should be included within the scope of protection of the present invention. Therefore, the scope of protection of the present invention should be determined by the scope of the claims.

Claims

1. A numerical method for complex modal analysis of a rotating machine rotor system, characterized by, The method includes: The differential equations of the free vibration of the rotor-bearing system are obtained, and the corresponding characteristic equations are determined based on the characteristic problem; wherein, the characteristic problem is used to analyze the vibration characteristics of the system. When the stiffness matrix in the characteristic equation is a non-singular matrix, an initial subspace basis matrix is ​​generated; the initial subspace basis matrix defines an m-dimensional subspace in an n-dimensional complex space; m is the number of eigenvalues ​​to be solved. By calculating the subspace projection matrix of the matrix in the characteristic equation onto the subspace defined by the initial subspace basis matrix, the projection problem of the characteristic problem onto the subspace is constructed. Execute the iterative format of the projection problem to solve the subspace feature problem; Using the solution results of the subspace feature problem, the eigenvalues ​​and eigenvectors of the feature equation are iteratively updated and normalized; If the convergence criterion is met, check the residuals of the characteristic equation and output the final solution of the characteristic equation; The step of determining the corresponding characteristic equation based on the characteristic problem includes: If the characteristic problem is to find the m largest eigenvalues, the characteristic equation is: ; It is a quality matrix. It is the sum of the damping matrix and the gyroscope matrix. It is the stiffness matrix; It is a diagonal matrix, where each diagonal element represents an eigenvalue to be found. It is the eigenvector matrix of the aforementioned feature problem. Each column represents a complex eigenvector; If the characteristic problem is to find the minimum m eigenvalues, then the characteristic equation is: ; It is an m-order diagonal matrix. diagonal element t s for Eigenvalues ​​represented by diagonal elements The reciprocal of U; U is The s-th column vector of the 1 / 2 matrix U is the eigenvalue. The corresponding feature vector; The order of is equal to the number of eigenvalues ​​sought; The generation of the initial subspace basis matrix includes: generating m orthogonal subspace basis vectors using the Gram-Schmidt method; the number of subspace basis vectors being equal to the number of eigenvalues ​​being solved; and combining the m subspace basis vectors column-wise into... 1-th order matrix , which serves as the basis matrix of the initial subspace.

2. The method according to claim 1, characterized in that, The characteristic equation is In the case of constructing the projection problem of the characteristic problem onto the subspace defined by the initial subspace basis matrix by calculating the subspace projection matrix corresponding to the matrix in the characteristic equation, the steps include: remember and Substituting the estimated eigenvalue matrix and the subspace eigenvector matrix obtained in the i-th iteration into the characteristic equation yields... i is greater than or equal to 1; Constructed in The orthogonal projection on the defined subspace yields ; in, The residual matrix; , , They are respectively , , The projection matrix of the matrix, and , , ; Let be the residual matrix of the subspace projection problem and ; express The conjugate transpose of .

3. The method according to claim 2, characterized in that, The iterative format for performing the projection problem, and solving the subspace feature problem, includes: Constructing inner iterative schemes to obtain matrices And satisfy ; For sequence Assign initial values, let ; For the j-th iteration, construct the m-th order matrix. system of equations ; Solve Generalized characteristic problem ; in, For the eigenvector matrix of the generalized feature problem, For the eigenvalue matrix of the generalized characteristic problem, yes 1-order matrix.

4. The method according to claim 2, characterized in that, The process of checking the residuals of the characteristic equation and outputting the final solution result of the characteristic equation, under the condition of satisfying the convergence criterion, includes: Calculate the difference between the results of two adjacent iterations, and if the difference between the results of two adjacent iterations is less than a specified error, determine that the convergence criterion is met and exit the outer iteration. Substitute the kth set of eigenvalues ​​and the corresponding eigenvectors into the characteristic equation to calculate the residual vector; k is a natural number between 1 and m. If the magnitude of the residual vector is less than the error threshold, the calculation result of the k-th set of eigenvalues ​​and the corresponding eigenvector is determined to be correct and is output as the final solution result.

5. A numerical apparatus for complex modal analysis of a rotating machinery rotor system, applied to the method described in any one of claims 1 to 4, characterized in that, The device includes an equation determination module, a subspace generation module, a matrix projection module, a subspace solving module, an iterative update module, and a verification output module, wherein: The equation determination module is used to obtain the differential equation of the free vibration of the rotor-bearing system and determine the corresponding characteristic equation based on the characteristic problem; wherein, the characteristic problem is used to analyze the vibration characteristics of the system. The subspace generation module is used to generate an initial subspace basis matrix when the stiffness matrix in the characteristic equation is a non-singular matrix; the initial subspace basis matrix defines an m-dimensional subspace in an n-dimensional complex space; m is the number of eigenvalues ​​to be solved; The matrix projection module is used to construct the projection problem of the characteristic problem in the subspace by calculating the subspace projection matrix of the matrix in the characteristic equation on the subspace defined by the initial subspace basis matrix; The subspace solving module is used to execute the iterative format of the projection problem and solve the subspace feature problem; The iterative update module is used to iteratively update the eigenvalues ​​and eigenvectors of the characteristic equation and normalize them using the solution results of the subspace characteristic problem. The verification output module is used to verify the residuals of the characteristic equation and output the final solution result of the characteristic equation, provided that the convergence criterion is met.

6. An electronic device comprising a memory and a processor, the memory storing a computer program executable on the processor, characterized in that, When the processor executes the program, it implements the steps of the method according to any one of claims 1 to 4.

7. A computer-readable storage medium having a computer program stored thereon, characterized in that, When executed by a processor, the computer program implements the steps of the method described in any one of claims 1 to 4.

Images