Method for performing stability analysis and transient calculation on fluid-filled rotor on basis of multi-frequency whirl decomposition

Abstract

Disclosed in the present invention is a method for performing stability analysis and transient calculation on a fluid-filled rotor on the basis of multi-frequency whirl decomposition. The implementation steps of the method comprise: 1) establishing a "dual-rotor" finite element model, acquiring an uncoupled fluid cross-stiffness state space, and calculating feature values corresponding to modes in various orders of a system; 2) on the basis of a multi-frequency whirl hypothesis, calculating, order by order, a fluid cross stiffness under the feature value corresponding to the mode in each order of the system, coupling the fluid cross stiffnesses into a state space equation of the system, and performing iterative calculation to obtain a converged feature value of each order of the system under the influence of fluid-structure interaction; 3) on the basis of the converged feature values of the system, calculating a liquid excitation force of the mode in each order, and decomposing, on the basis of a modal superposition principle, the total displacement response of the system into the modes in the various orders of the system; and 4) in view of an explicit-implicit-segregated solution framework, performing, with respect to the mode in each order, transient dynamic calculation and analysis on a dynamic rotor model obtained after fluid-structure interaction. The method solves the problem of a traditional lumped-mass model in a traditional fluid-structure-interaction rotor stability analysis method only being capable of processing a single whirl frequency.

Description

A liquid-filled rotor stability analysis and transient calculation method based on multi-frequency whirl decomposition TECHNICAL FIELD

[0001] The present application relates to the field of finite element simulation of computer-aided calculation of liquid-solid coupling rotor dynamics, and particularly relates to a liquid-filled rotor stability analysis and transient calculation method based on multi-frequency whirl decomposition. BACKGROUND

[0002] A centrifuge is a typical liquid-filled rotor equipment mainly using liquid as working medium. When the centrifuge is in a partial liquid-filled working condition, the liquid-solid coupling effect of the liquid and the rotor will cause an unstable region in the speed range near the first-order critical speed. In this region, the liquid-filled rotor will have a large amplitude of non-synchronous vibration due to the influence of liquid excitation, which will affect the running stability of the liquid-filled rotor and even cause production accidents, thereby delaying the research and development process of large-scale centrifuges in China.

[0003] The traditional research method for liquid-filled rotors is to consider the entire rotor as a lumped mass model. However, the low-degree-of-freedom model cannot accurately represent the increasingly complex dynamic characteristics of the liquid-filled rotor of the centrifuge, and it also cannot consider the influence of high-order modes of the system. To clarify the instability mechanism of the liquid-filled rotor under liquid-solid coupling and the stability variation law of high-order modes, an accurate model of the liquid-filled rotor of the centrifuge can be established based on the finite element method. However, the existing fluid force model can only consider single-frequency whirl of the rotor, and the introduction of the finite element method means that the contradiction between multi-frequency whirl and the single-frequency whirl fluid force model needs to be solved. Therefore, based on the modal superposition principle, the multi-frequency whirl of the rotor can be considered as the sum of multiple single-frequency whirls, and the displacement can be decomposed into each mode and analyzed and calculated for the stability interval of the liquid-filled rotor. Then, based on the converged system eigenvalues, the linear-nonlinear explicit-implicit separation solution method is used to efficiently solve the transient response of the liquid-filled rotor under liquid-solid coupling. SUMMARY

[0004] The liquid-filled rotor stability analysis and transient calculation method based on multi-frequency whirl decomposition provided by the present application can solve the problems of the traditional analysis of the liquid-solid coupling dynamics of the liquid-filled rotor, such as the difficulty of accurately representing the overall dynamics of the rotor system using the lumped mass model, and the inability to analyze the stability interval variation of the high-order modes of the system. In addition, the present application also solves the problem of the existing single-frequency whirl fluid force model that cannot meet the multi-frequency whirl assumption brought by the finite element method. To improve the transient calculation efficiency of the liquid-filled rotor under multi-modal liquid-solid coupling, the present application also uses the linear-nonlinear explicit-implicit separation transient calculation method based on the G-alpha numerical calculation method, which overcomes the problem of the conventional Newmark-beta method that has no suppression ability for high-frequency mode oscillation response, and finally quickly, stably and accurately solves the nonlinear transient dynamics response of the system.

[0005] To solve the above technical problems, the present application provides the following technical solutions:

[0006] S1, a liquid-filled rotor system finite element model is established, and modal analysis is performed on the liquid-filled rotor system finite element model;

[0007] S2, the cross stiffness of the fluid acting on the liquid-solid coupling system in each modal is obtained, and is iterated with the modal analysis to obtain the converged system characteristic value;

[0008] S3, based on the modal superposition principle, the external force and total displacement of the rotor are decomposed into each modal, and the fluid excitation force is calculated step by step, and the linear-nonlinear explicit-implicit separation method is combined to solve the converged modal displacement response of the liquid-filled rotor system step by step.

[0009] S1 includes the following steps:

[0010] S1.1 measures the actual size of the rotor, and divides the nodes, wherein the centrifuge sleeve is regarded as an outer rotor with a rotational speed of 0, the shaft is regarded as an inner rotor, and a rotor finite element model is constructed based on Timoshenko beam theory and finite element method;

[0011] S1.2 converts the rotor dynamics equation into a state space form, and performs modal analysis to calculate the modal characteristic values and characteristic vectors under the action of fluid-structure coupling;

[0012] S2 includes the following steps:

[0013] S2.1, according to the fluid cross stiffness formula under the single frequency vortex, the fluid cross stiffness corresponding to each pair of conjugate characteristic values of each modal is calculated one by one, the cross stiffness is coupled into the state space in S1, and each modal characteristic value is iteratively solved until stable convergence;

[0014] S2.2, for each pair of conjugate characteristic values of each modal, the characteristic value with the largest real part is selected as the stability criterion of the liquid-solid coupling liquid-filled rotor system under the rotational speed at this moment.

[0015] S2.3, the stability of each modal is analyzed, and the rotational speed range with the characteristic value real part greater than 0 is the instability interval of the liquid-solid coupling liquid-filled rotor.

[0016] S3 includes the following steps:

[0017] S3.1, the external excitation and displacement response of the rotor are decomposed into each modal, wherein the initial time displacement is simplified as the displacement response without fluid excitation;

[0018] S3.2, based on the converged system characteristic value, the fluid excitation force corresponding to each modal is calculated step by step;

[0019] S3.2 Based on the transient calculation method of apparent-implicit separation, the linear and nonlinear nodes of the rotor are separated, and the G-alpha numerical calculation method is used to obtain the integral formula of the rotor motion equation;

[0020] S3.3 The Newton-Raphson iteration formula is used to solve the nonlinear node response under each order mode step by step, and the displacement, velocity and acceleration responses of all nodes under each order mode are inversely deduced.

[0021] S3.4 Based on the modal superposition principle, the total response of the liquid-filled rotor under the liquid-solid coupling action is obtained by superimposing the displacement of each order mode.

[0022] The present application has the following advantages:

[0023] The liquid-filled rotor stability analysis and transient calculation method based on multi-frequency vortex decomposition provided by the present application regards the multi-frequency vortex of the liquid-filled rotor as the sum of multiple single-frequency vortexes by introducing the finite element modeling method and combining the modal superposition principle. The accurate modeling of the centrifuge rotor under complex elastic support conditions is realized, so that compared with the traditional lumped mass model, the transient simulation result of the liquid-filled rotor under the liquid-solid coupling action is closer to the actual working condition. The present application introduces the linear-nonlinear apparent-implicit separation algorithm to extract the nonlinear nodes of the rotor, realizes the iteration calculation dimension of the high-dimensional rotor finite element model, and greatly reduces the time of calculating the transient response. By combining the G-alpha numerical calculation method with the explicit trapezoidal method to calculate the transient response of the system step by step, the problem of lack of suppression ability of the traditional Newmark-beta numerical calculation method for high-order modes of the system is avoided. BRIEF DESCRIPTION OF DRAWINGS

[0024] Fig. 1 is a liquid-filled rotor finite element model established in the embodiment of the present application.

[0025] Fig. 2 is a stability analysis method flowchart under the action of multi-frequency vortex in the embodiment of the present application.

[0026] Fig. 3 is a curve diagram of the first four order stability of the liquid-filled rotor system under different fluid viscosity conditions in the embodiment of the present application.

[0027] Fig. 4 is a liquid-filled rotor stability analysis and dynamics calculation device based on multi-frequency vortex decomposition in the embodiment of the present application.

[0028] Fig. 5 is the decomposition result of the unbalance response of each order rotor of the liquid-filled system in the embodiment of the present application.

[0029] Fig. 6 is a simulation transient response spectrum waterfall diagram of the liquid-filled rotor system in the embodiment of the present application.

[0030] Fig. 7 is an experimental vibration response spectrum waterfall diagram of the liquid-filled rotor system in the embodiment of the present application. DETAILED DESCRIPTION

[0031] As shown in Figure 1, a kind of based on flexible basis transfer function's whole rotor dynamics transient calculation method and device of embodiment of the application:

[0032] The present example is based on the principle of modal decomposition, and the multi-frequency whirl of the liquid-filled rotor is regarded as the sum of multiple single-frequency whirls. The stability of each modal of the liquid-filled rotor system is analyzed step by step, the stability interval of each modal is analyzed, and the linear-nonlinear explicit-implicit separation algorithm is used for rotor transient response analysis. The specific steps are as follows:

[0033] The scope of protection of the present application is not limited to the description of the present embodiment.

[0034] First, a finite element model of the liquid-filled system is established, and modal analysis is performed on the liquid-filled system.

[0035] 1) The liquid-filled rotor system is divided into beam elements, disc elements and bearing elements. The rotor shaft is modeled based on the Timoshenko beam theory, the disc elements are regarded as concentrated mass elements, and the bearing elements are regarded as support elements to provide support reaction for the rotor system. Finally, the dynamic equation of the rotor system is established as follows:

[0036] In the formula, M, C, G and K are the total mass, damping, gyroscopic effect and stiffness matrix of the assembled liquid-filled rotor system, And q are the acceleration, velocity and displacement vectors of the liquid-filled rotor system, Ω is the rotational speed of the liquid-filled rotor system, F is the external excitation force on the liquid-filled rotor system, including unbalance force F u , support force F b And fluid excitation force F l , that is: F=F u +F b +F l

[0037] 2) Convert the rotor system dynamic equation into state space form:

[0038] In the formula, A and B are coefficient matrices in the process of converting the rotor dynamic equation into state space, z is the state variable, and Q is the state space force vector:

[0039] Solve the eigenvalues of the state space equation, obtain the characteristic values λ m And modal matrix ψ of each modal of the liquid-filled system, wherein: λ mn , λ m(n+1) = σ m ±jw m , m=1, 2, 3...., n=1, 3

[0040] The subscript m represents the modal order corresponding to the solution, and n represents the eigenvalue sequence number corresponding to the modal order. Since each modal order corresponds to two sets of conjugate roots, i.e., four eigenvalues, n e {1, 2, 3, 4}. The real part σ m is the stability criterion of the mth order modal of the liquid-filled rotor system, and the imaginary part w m is the mth order natural frequency of the liquid-filled rotor system, and j is the imaginary unit.

[0041] Secondly, the cross stiffness of the fluid acting on the liquid-solid coupling system in each modal order is obtained, and the modal analysis iteration cycle is obtained to obtain the converged system eigenvalue;

[0042] 1) The fluid dynamics equation is established based on the fluid N-S equation, and the finite difference method is used to solve the differential equation to derive the fluid cross stiffness of the system in the mth order modal under the fluid-structure interaction:

[0043] In the formula, S m is the equivalent cross stiffness of the mth order of the liquid-filled rotor system, k XX , k XY , k YX and k YY are the equivalent stiffness coefficients of the fluid under the liquid-solid coupling action, and the form is:

[0044] m l is the mass of the liquid in the rotor drum, v is the dynamic viscosity of the fluid, r o is the inner wall radius of the rotor drum, and are auxiliary variables constructed when solving the N-S fluid dynamics equation, and the purpose is to decouple the equation in two perturbation directions. The superscripts "'" and " " represent the third derivative and the second derivative of the corresponding function, respectively.

[0045] 2) The mth order modal fluid cross stiffness is coupled with the rotor system state space expression to calculate the nth eigenvalue λ mn of the mth order modal of the new liquid-filled rotor system under the liquid-solid coupling action. When λ mn satisfies the iterative convergence condition (the two-norm of the calculated eigenvalues before and after satisfies the convergence accuracy 1e-4), it is considered that the nth eigenvalue of the mth order modal of the rotor system is calculated, otherwise the fluid cross stiffness in 1) is recalculated and solved until the convergence requirement is met. Figure 2 is a flowchart of the stability analysis method under the action of multi-frequency whirl.

[0046] 3) The largest eigenvalue in the two sets of eigenvalues of the mth order is selected as the mth order stability criterion, and if the largest real part σ mIf ≥0, the rotor system is unstable at the mth order, otherwise, the system is stable at the mth order. Fig. 3 is a diagram of the first four orders of stability curves of the liquid-filled rotor system under different fluid viscosity conditions.

[0047] Step 3: Based on the modal superposition principle, the external force and total displacement of the rotor are decomposed into each order of modal, and the fluid excitation force is calculated step by step, and the linear-nonlinear explicit-implicit separation method is combined to solve each order of modal displacement response of the liquid-filled rotor system after iteration and convergence.

[0048] 1) In order to calculate the system response at each order of modal, the unbalance force, support force and other forces acting on the rotor system need to be decomposed into the corresponding modal. According to the rotor dynamics equation, the physical space coordinates are transformed into modal space coordinates through modal matrix: z = Ψg

[0049] Wherein, z is the state variable of the liquid-filled rotor system in the physical space, g is the state variable of the liquid-filled rotor system in the modal space, Ψ is the complex modal matrix of the liquid-filled rotor, ψ i is the modal column vector of the i-th column of the liquid-filled rotor system: Ψ = [ψ1ψ2…ψ 2n ]

[0050] It is brought into the state space expression:

[0051] Then the mth order n-th modal coordinate of the force acting in the modal space is:

[0052] In the formula, is the mth order n-th component in the modal space. Since there are two groups of characteristic values in the mth order modal, the external force of the system in the mth order modal in the physical space is:

[0053] In the formula, is the rest of the external force in the physical space except the fluid force, Ψ 4(m-1)+n is the eigenvector corresponding to the mth order n-th eigenvalue in the complex modal matrix of the liquid-filled rotor system. Since it is a non-square matrix, the generalized inverse should be used when inverting it. Fig. 4 is the decomposition result of the rest of the external force except the fluid force.

[0054] 2) According to the formula in the preceding text, the mth order fluid force of the system is:

[0055] In the formula, is the mth order liquid excitation force of the liquid-filled rotor system, and the mth order displacement q m of the liquid-filled rotor system is:

[0056] Figure 5 is the decomposition result of the unbalance response of each order of the liquid-filled rotor system. Then the external excitation force and F m are:

[0057] 3) Linear-nonlinear explicit-implicit separation transient dynamics solution is performed on the system, the linear and nonlinear nodes of the rotor system are separated based on the explicit-implicit separation transient calculation method, since the rotor dynamics modeling in this patent is based on the finite element method, there are high-order modes of the system, and the traditional Newmark-β method has no suppression ability for high-frequency mode oscillation response, and the G-α method can introduce numerical damping to suppress the low-damping numerical oscillation effect caused by high-order modes. Therefore, the integral formula of the rotor motion equation is obtained by using the G-α calculation method.

[0058] After separating and reordering the linear and nonlinear nodes of the rotor system, the motion equation can be written as:

[0059] In the formula, m represents the mth mode, subscript i represents the linear node, and s represents the nonlinear node, M ii , M is , M si and M ss are the mass sub-matrices of the liquid-filled rotor system after linear-nonlinear separation and reordering, and the damping matrix C and the stiffness matrix K are the same. and are the mth order unbalance force and other forces of the linear and nonlinear nodes of the liquid-filled rotor system, is the mth order nonlinear fluid force of the liquid-filled rotor system, and are the acceleration response, velocity response and displacement response of the linear nodes of the liquid-filled rotor system, and are the mth order acceleration response, velocity response and displacement response of the nonlinear nodes of the liquid-filled rotor system, w m is the mth order natural frequency of the liquid-filled rotor system.

[0060] Based on the recursive formula of the G-α numerical integration method, the velocity and acceleration formulas of the liquid-filled rotor system can be obtained by the following formula:

[0061] In the formula, and are the mth order acceleration response, velocity response and displacement response of the corresponding nodes of the liquid-filled rotor system at t+Δt, The r-2 order derivative of the displacement response of the corresponding node of the liquid-filled rotor system at time t, when r is 2, 3 and 4 respectively, represents the displacement response, the velocity response and the acceleration response of the liquid-filled rotor system respectively. r r are the rth coefficients of the G-alpha numerical algorithm used to solve the liquid-filled rotor system, which are in the following form:

[0062] In the formula, α m , α f and γ are parameters introduced in the G-alpha numerical calculation process, ρ i is the numerical damping introduced in the G-alpha numerical calculation method, and the greater ρ i , the greater the numerical damping of the algorithm, and the more obvious the suppression effect on the high-frequency modal oscillation response. Further derivation of the integral relationship satisfied by the nonlinear node is as follows:

[0063] In the formula, H and f0 m are the coefficient matrices derived in the process of solving the displacement response of the mth order nonlinear node of the liquid-filled rotor system, and their specific forms are as follows: Π pq = α1M pq + β1C pq + K pq , p, q = i, s

[0064] 4) To solve the nonlinear node displacement, a function is constructed, and the essence of calculating the displacement of the nonlinear node by using the Newton-Raphson method is to solve the solution of , that is, the displacement of the nonlinear node at time t+Δt. The iterative form is as follows:

[0065] In the formula, n is the iteration step number, and are the displacement responses of the mth order iteration process of the nonlinear node of the liquid-filled rotor system at the n+1th and nth steps respectively, is the function to be solved corresponding to the constructed Jacobian matrix, and the expression of the Jacobian matrix J is as follows:

[0066] After several iterations, the displacement response of the nonlinear node of the liquid-filled rotor system can be obtained, and the mth order displacement response of the linear node of the liquid-filled rotor system at time t+Δt is solved by iteration ​​

[0067] wherein, is the force vector of the mth order of the liquid-filled rotor system under the action of linear nodes and nonlinear nodes at time t, and its specific form is:

[0068] The displacement of the mth order nonlinear node of the liquid-filled rotor system is The linear node displacement of the liquid-filled rotor system is calculated by substituting the above formula into the formula: The linear node velocity of the liquid-filled rotor system is solved by substituting it into the G-α recursive formula: And the acceleration After completing the rotor dynamics transient calculation of all modes step by step, the superposition is performed, and the transient response calculation of the liquid-filled rotor system under the influence of multi-frequency whirl is completed:

[0069] wherein, q k is the total displacement response of the liquid-filled rotor system under the influence of each order mode of linear or nonlinear nodes, is the mth order displacement response of the liquid-filled rotor system, and x is the modal order required in the process of solving the transient response of the liquid-filled rotor system

[0070] Fig. 6 is a liquid-filled rotor system transient response displacement spectrum waterfall diagram based on the multi-frequency whirl decomposition principle.

[0071] The application also discloses a liquid-filled rotor stability analysis and dynamics calculation device based on multi-frequency whirl decomposition, and specific components are referred to, comprising:

[0072] The decomposition module is configured to perform system response and excitation force decomposition based on modal superposition theory;

[0073] The analysis module is configured to perform modal analysis of the “double rotor” model. The modal characteristics of the rotor model are embodied as the characteristic values and characteristic vectors of the system, and the system stability is reflected by the characteristic values. The system modal information is provided to the fluid module, and the iteration is performed with the fluid module until the calculation result converges;

[0074] The fluid module is configured to perform fluid dynamics equation solving method based on finite difference method to obtain the fluid excitation force part on the rotor wall, and feed back the fluid cross stiffness under fluid-structure coupling action to the analysis module, and the iteration is performed with the analysis module until the calculation result converges;

[0075] The coupling module is configured to perform the coupling of the fluid equivalent cross stiffness and the state space equation of the “double rotor” model, and update the rotor dynamics model under the liquid-solid coupling action of the system;

[0076] Computing module: configured to perform a linear-nonlinear decoupling calculation method based on G-alpha numerical solution method to solve the coupled dynamic model step by step.

[0077] To verify the rationality of the algorithm, the transient simulation of the liquid-filled rotor system is carried out based on the method of the patent, and compared with the experimental phenomenon:

[0078] In this example, the rotor finite element model is shown in Figure 1, and the intermediate bearing is at nodes 11 and 19, and the electromagnetic bearing is at node 27, and the electromagnetic excitation force can be applied to control the rotor vibration, and the unbalance amount of 3.5g·mm∠0° is applied at node 25, and the fluid force is applied at node 36 (at the liquid disc). In order to reflect the vibration response characteristics of the entire liquid-filled rotor system under fluid-structure coupling, the rotor is set to accelerate from 0rpm to 3000rpm with an acceleration of 4πrad / s 2 , and the frequency spectrum waterfall plot of the vibration signal is drawn, as shown in Figure 6.

[0079] Figure 7 is the experimental result of the liquid-filled rotor system at 1200rpm, and it can be seen that the non-synchronous vibration frequency of fluid excitation appears near the main frequency, and in the simulation result of Figure 6, the vibration component under the action of fluid excitation also appears obviously near the main frequency. This verifies the effectiveness of the liquid-filled rotor transient solution method based on multi-frequency vortex decomposition proposed in the patent, solves the problem that the traditional lumped mass model in the traditional liquid-solid coupling rotor stability analysis method can only handle single vortex frequency, and realizes efficient solution of multi-order modal of complex finite element model.

[0080] The above is a specific embodiment of the present application, but the protection scope of the present application is not limited thereto, and those skilled in the art can make appropriate adjustments without departing from the principles of the present application, and these adjustments shall be included in the protection scope of the present application.

Claims

1. A method for stability analysis and transient calculation of a liquid-filled rotor based on multi-frequency whirl decomposition, characterized in that, The implementation steps of the method include the following: S1, a finite element model of the liquid-filled rotor system is established, and modal analysis is performed on the finite element model of the liquid-filled rotor system; S2, the cross stiffness of the fluid acting on the liquid-solid coupling system in each mode is obtained, and is iterated with the modal analysis to obtain the converged system characteristic value; S3, based on the modal superposition principle, the external force and total displacement of the rotor are decomposed into each mode, and the fluid excitation force is calculated step by step, and the linear-nonlinear explicit-implicit separation method is combined to solve the liquid-filled rotor system step by step to obtain the converged modal displacement response of each order; S2 includes the following steps: S2.1 According to the fluid cross stiffness formula under single frequency whirl, the fluid cross stiffness corresponding to each pair of conjugate characteristic values of each mode is calculated one by one, the cross stiffness is coupled into the state space in S1, and each modal characteristic value is iteratively solved until stable convergence is achieved; S2.2 For each pair of conjugate characteristic values of each mode, the characteristic value with the largest real part is selected as the stability criterion of the liquid-solid coupling liquid-filled rotor system at the current speed; S2.3 Stability analysis is performed on each mode, and the speed range in which the real part of the characteristic value is greater than 0 is the instability interval of the liquid-solid coupling liquid-filled rotor.

2. The method of claim 1, wherein, In S1: S1.1 Measure the actual size of the rotor and divide the nodes, regard the centrifuge sleeve as an outer rotor with a speed of 0, regard the shaft as an inner rotor, construct a rotor finite element model based on Timoshenko beam theory and finite element method, and establish a liquid-filled rotor dynamics equation; S1.2 Convert the rotor dynamics equation into a state space form, and perform modal analysis to calculate the modal characteristic values and characteristic vectors under the action of fluid-structure coupling.

3. The liquid-filled rotor stability analysis and transient calculation method based on multi-frequency whirl decomposition according to claim 1, S3 includes the following steps: S3.1 Decompose the external excitation and displacement response of the rotor into each mode, wherein the initial time displacement is simplified as the displacement response without fluid excitation; S3.2 Calculate the fluid excitation force corresponding to each mode based on the converged system characteristic value; S3.2 Based on the explicit-implicit separation transient calculation method, separate the linear and nonlinear nodes of the rotor, use the G-alpha numerical calculation method to obtain the integral formula of the rotor motion equation; S3.3 Solve the nonlinear node response of each mode step by step using the Newton-Raphson iteration formula, and then back-propagate the displacement, velocity and acceleration response of all nodes of each mode; S3.4 Based on the modal superposition principle, superimpose the displacement of each mode to obtain the total response of the liquid-filled rotor under the action of liquid-solid coupling.

4. The method of claim 1, wherein the system for implementing the method comprises: The decomposition module, the analysis module, the fluid module, the coupling module and the calculation module are connected; The analysis module inputs the system modal information to the fluid module, the fluid module is connected with the coupling module, the coupling module is connected with the analysis module and the calculation module respectively, and the calculation module is connected with the decomposition module through the next order calculation; The decomposition module is configured to perform system response and excitation force decomposition based on modal superposition theory; The analysis module is configured to perform modal analysis of the "dual rotor" model, wherein modal characteristics of the rotor model are embodied as eigenvalues and eigenvectors of each order of the system, and system stability is reflected by the eigenvalues, and the fluid module is provided with system modal information, and iteration is performed with the fluid module until the calculation result converges; The fluid module is configured to perform a finite difference method-based fluid dynamics equation solving method to obtain fluid excitation force acting on the rotor wall surface, and feedback fluid cross stiffness under fluid-structure coupling to the analysis module, and iteration is performed with the analysis module until the calculation result converges; The coupling module is configured to perform coupling of the fluid equivalent cross stiffness and the "dual rotor" model state space equation, and update the rotor dynamics model under the fluid-structure coupling of the system; The calculation module is configured to perform a linear-nonlinear decoupling calculation method based on the G-alpha numerical solving method to solve the coupled dynamics model step by step.

Images