Nonlinear blade flutter evaluation method based on accumulated work measurement and fatigue life

Abstract

The invention relates to a nonlinear blade flutter evaluation method based on accumulated work measurement and fatigue life, which comprises the following steps of: obtaining the maximum allowable vibration amplitude of a blade by taking the fatigue life as a target, and enabling the blade to vibrate according to the maximum allowable vibration amplitude on the premise of considering blade nonlinearity caused by overlarge amplitude in simple harmonic vibration of the blade, and solving the unsteady flow by using a harmonic balance method to obtain pneumatic accumulated work, and judging whether the flutter characteristic of the blade meets the safety requirement of the design or not according to the positive and negative accumulated work. According to the method, whether the blade vibration has the risk of amplitude overrun or not is judged on the basis of accumulated work measurement and fatigue life, the problem of conservative chatter linear analysis is effectively solved, and the blade chatter judgment precision is improved.

Description

Technical Field

[0001] The present invention relates to the technical field of impeller machinery design, and particularly to an evaluation method for blade flutter by using accumulated work measurement. Background Art

[0002] A turbomachine is an important part of the core engine of an aeroengine, and the blade is an extremely important part of the turbomachine. Blade flutter is a self-excited vibration, that is, the unsteady aerodynamic force causing the blade oscillation comes from the oscillation of the blade, which is an aeroelastic instability phenomenon. When flutter occurs, the vibration amplitude of the blade initially grows exponentially. If there is not enough mechanical damping, the blade vibration will quickly exceed the limit and cause the blade to break. If there is a certain amount of mechanical damping, the blade vibration may enter a limit cycle state with a constant amplitude and will not exceed the maximum allowable vibration amplitude.

[0003] After flutter occurs, as the vibration amplitude of the blade increases, the blade vibration enters non-linear fluid-structure coupling vibration and may enter limit cycle vibration. However, as long as the vibration amplitude does not exceed the limit, this kind of flutter has no substantial harm.

[0004] At present, the general flutter analysis method is a linear analysis method, which can only judge whether flutter occurs. That is, this flutter analysis method is relatively conservative. It does not pay attention to whether the vibration amplitude will increase to exceed the limit after flutter occurs, so it cannot judge, and even less can it judge in time whether the blade design can meet the safety requirements. Summary of the Invention

[0005] The purpose of the present invention is to avoid the deficiencies of the prior art and provide a non-linear blade flutter evaluation method based on accumulated work measurement and fatigue life, which obtains the maximum allowable vibration amplitude of the blade according to the fatigue life target, and in the numerical analysis of aerodynamic damping, makes the blade vibrate at the maximum allowable vibration amplitude, and then uses the harmonic balance method to solve the unsteady flow, that is, considering the flow non-linearity caused by the excessive blade vibration amplitude, obtains the aerodynamic accumulated work, and finally judges whether the flutter characteristics of the blade meet the requirements according to the positive and negative of the accumulated work.

[0006] To achieve the above purpose, the technical solution adopted by the present invention is: a non-linear blade flutter evaluation method based on accumulated work measurement and fatigue life, including the following steps: Step 1: Establish a computational grid for the fluid domain of the blade, perform a steady calculation of the blade flow field, obtain the characteristic curves of flow rate vs. pressure ratio and flow rate vs. efficiency at different rotational speeds of the blade, so as to determine the calculation conditions for blade flutter; at the same time, based on structural static and dynamic analysis, obtain the modal shapes corresponding to at least the first three vibration modes of the blade vibration and the natural angular frequencies ; Step 2: Take any one of the modal shapes Interpolate it onto the computational grid of the fluid domain of the blade, and calculate the maximum allowable vibration amplitude of the blade according to the maximum allowable alternating stress that satisfies the high-cycle fatigue life of the blade. ; Step 3: Under the conditions of the calculation working condition, make the blade perform simple harmonic vibration with any first-order modal shape , natural angular frequency and maximum allowable vibration amplitude ; On the premise of considering the non-linearity of the blade caused by excessive amplitude in the simple harmonic vibration of the blade, specify the fundamental frequency and harmonic information of the unsteady flow caused by the simple harmonic vibration of the blade, and initialize the unsteady blade flutter calculation based on the steady calculation results of the blade flow field. Thus, use the harmonic balance method to analyze the unsteady flow caused by the blade vibration and obtain the unsteady pressure on the blade surface and the blade vibration velocity ; Finally, integrate the unsteady pressure on the blade surface and the blade vibration velocity within one vibration period to obtain the accumulated work; Step 4: Repeat Step 2 to Step 3, calculate the accumulated work corresponding to each vibration mode obtained, and evaluate the blade flutter. Specifically: If the accumulated work is positive, the work done by the airflow on the blade increases within the period, the amplitude will continue to increase, there is a risk of amplitude exceeding the limit, and the fatigue life requirement of the corresponding vibration mode is not satisfied;

[0007] Furthermore, in the steady calculation of the blade flow field, solve the FANS equation for determining the calculation working condition of blade flutter, and express the FANS equation in the cylindrical coordinate system as:

[0008] In the formula: is the physical time, is the conserved variable per unit volume of the flow field, The expression of is: where , , are the velocity components of the fluid in the axial, circumferential, and radial directions, is the total energy of the fluid per unit mass; are the coordinate positions in the axial, circumferential, and radial directions respectively, are the velocity components of the computational grid of the fluid domain moving in the axial, circumferential, and radial directions, respectively, , , represent the axial, circumferential, and radial fluxes generated due to the movement of the computational grid of the fluid domain, respectively; are the convective fluxes of the fluid in the axial, circumferential, and radial directions, respectively. The convective flux represents the physical quantity transported per unit area per unit time, The expression formula of is: is the total enthalpy of the fluid per unit mass, is the static pressure of the fluid; are the viscous fluxes of the fluid in the axial, circumferential, and radial directions, respectively. The viscous flux is used to describe the diffusive transport of momentum and energy caused by viscous stress and heat conduction, The expression of

[0009]

[0010] is: In the formula, is the normal stress between fluids, is the shear stress between fluids, is the heat flux between fluids, , , are the axial, circumferential, and radial components of the heat flux between fluids; is the source term, and the source term is used to describe the influence of the internal or external of the flow field system on the conserved quantity, The expression of is: After solving the FANS equation, the characteristic curves of flow rate versus pressure ratio and flow rate versus efficiency at different rotational speeds of the blade are obtained, and the calculation conditions are selected according to the characteristic curves for flutter analysis; Among them, the calculation conditions include the design condition, the blockage condition, and the near-stall condition.

[0011] Furthermore, obtaining the modal shapes and the natural angular frequencies corresponding to at least the first three vibration modes of the blade vibration described in step one, specifically: Based on the structural dynamics equation of the blade vibration: , In the formula: is the mass matrix, is the damping matrix, is the stiffness matrix, f is the aerodynamic force acting on the blade surface, is the vibration displacement of the blade, is the first derivative of, that is, the vibration velocity of the blade; is the acceleration of the blade; In modal analysis, the aerodynamic force f and the damping matrix C are not considered, and the following formula one is obtained: , Let be formula two: , In the formula: is the natural angular frequency of the blade, t is the physical time. Substituting formula two into formula one, formula three is obtained: , Formula three is the characteristic equation. Solving formula three can obtain the modal shapes corresponding to the vibration modes of each order of the blade and the natural angular frequency , as well as the modal stress corresponding to the modal shape and the modal displacement .

[0012] Furthermore, the specific content of step two is as follows: The maximum allowable alternating stress satisfying the high-cycle fatigue life of the blade is obtained from the Goodman diagram, and further, the maximum allowable vibration amplitude L of the blade is obtained as: , In the formula, is the modal stress, is the modal displacement.

[0013] Furthermore, in step three: Specify the fundamental frequency of the unsteady flow caused by the blade vibration, which is the natural frequency of the blade: , In the formula, is the natural angular frequency; At the same time, specify the inter-blade phase angle , and the inter-blade phase angle That is, the phase difference when the blades vibrate in a traveling wave mode, which is denoted by the calculation formula as: , In the formula: is the vibration nodal diameter of the blade, is the number of blades; When is an even number, takes values in [- / 2, / 2), where - / 2 and / 2 differ by , and only one of them is taken during calculation; When is an odd number, takes values in [-( -1) / 2, ( -1) / 2]; In addition, to analyze the unsteady flow caused by blade vibration, the specified harmonic information is: , where is the harmonic number, and is a positive integer greater than 1.

[0014] Furthermore, in step three, the steady calculation results of the blade flow field are used to initialize the unsteady blade flutter calculation, and then the harmonic balance method is used to analyze the unsteady flow caused by blade vibration to obtain the unsteady pressure acting on the blade surface corresponding to the unsteady flow and the blade vibration velocity . Among them, the formula of the harmonic balance method is: , In the formula: is the variable of the unsteady flow field around the blade, is the time-averaged quantity of the unsteady flow field, is the physical time, n is the harmonic number used to analyze the unsteady flow field, i is the order corresponding to different harmonics, and are the Fourier coefficients corresponding to the i th harmonic, is the i th harmonic corresponding angular frequency of the unsteady flow field; the harmonic balance method is used to represent the flow field variable of the unsteady flow as the sum of the time-averaged quantity and the pulsating quantity, and approximate the pulsating quantity with a Fourier series; The expression formula for representing the unsteady pressure acting on the blade surface using the harmonic balance method is: , where: is the time-averaged value of the unsteady pressure on the blade surface; Similarly, the vibration velocity of the blade is expressed as: , where: is the time-averaged value of the vibration velocity of the blade.

[0015] Furthermore, the accumulated work described in Step 3 is: , where: is the vibration velocity of the blade, is the fluid domain grid normal vector on the blade surface, is the differential of the blade surface area, is the blade vibration period, dt is the differential of physical time, and S is the total blade surface area.

[0016] Furthermore, based on the structural statics and dynamics analysis, at least the first three modal vibration modes and the natural angular frequencies are obtained.

[0017] The beneficial effects of the present invention are as follows: The evaluation method proposed by the present invention obtains the maximum allowable vibration amplitude of the blade according to the fatigue life target, performs simple harmonic vibration with the maximum allowable amplitude, uses the harmonic balance method to solve the unsteady flow, and calculates the aerodynamic accumulated work, which takes into account the flow nonlinearity caused by excessive blade vibration amplitude. When the flow nonlinearity caused by blade vibration is strong, high-order harmonics are introduced to handle the unsteady flow, and then the positive and negative of the accumulated work are used to judge whether the flutter characteristics of the blade meet the safety requirements of the design, and to judge whether there is a risk of amplitude overrun in blade vibration, effectively solving the problem that the linear flutter analysis is too conservative. BRIEF DESCRIPTION OF THE DRAWINGS

[0018] Figure 1 is a schematic flow chart of the method of the present invention; Figure 2 is the computational domain of the blade simulation in the specific example of the present invention; Figure 3 is a schematic diagram of the modal displacement of the blade in the specific example of the present invention; Figure 4 is a schematic diagram of the modal stress of the blade in the specific example of the present invention; Figure 5 is a schematic diagram of the accumulated work density distribution on the blade surface when the blade vibrates in the first mode and the vibration node diameter is 0 in the specific example of the present invention. Detailed implementation manners

[0019] The principles and features of the present invention will be described below in conjunction with the accompanying drawings. The examples given are only used to explain the present invention and are not intended to limit the scope of the present invention.

[0020] To achieve the above object, the present invention provides the following specific implementation manners: Example 1: As Figure 1 shown, a non-linear blade flutter evaluation method based on accumulated work measurement and fatigue life includes the following steps: S01. Establish a computational grid for the fluid domain of the blade and perform a steady calculation of the blade flow field. In the steady calculation of the blade flow field, solve the FANS equation for determining the calculation conditions of blade flutter, and express the FANS equation in the cylindrical coordinate system as:

[0021] In the formula: is the physical time, is the conserved variable per unit volume in the flow field, The expression of is: , is the fluid density, , , are the velocity components of the fluid in the axial, circumferential, and radial directions, is the total energy of the fluid per unit mass; are the coordinate positions in the axial, circumferential, and radial directions respectively, are the velocity components of the computational grid of the fluid domain moving in the axial, circumferential, and radial directions respectively, , , respectively represent the axial, circumferential, and radial fluxes generated due to the movement of the computational grid of the fluid domain; are the convective fluxes of the fluid in the axial, circumferential, and radial directions respectively. The convective flux represents the physical quantity transported through a unit area per unit time, The expression formula of is: , is the total enthalpy of the fluid per unit mass, is the static pressure of the fluid; are the viscous fluxes of the fluid in the axial, circumferential, and radial directions respectively. The viscous flux is used to describe the diffusive transport of momentum and energy caused by viscous stress and heat conduction, The expression of is:

[0022]

[0023] , wherein, is the normal stress between fluids, is the shear stress between fluids, is the heat flux between fluids, , , are respectively the components of the heat flux between fluids in the axial, circumferential, and radial directions; is the source term, which is used to describe the influence of the inside or outside of the flow field system on the conserved quantity, The expression of , After solving the FANS equation, the characteristic curves of flow rate vs. pressure ratio and flow rate vs. efficiency at different rotational speeds of the blade are obtained, thereby determining the calculation conditions for blade flutter. Then, the calculation conditions can be selected according to the characteristic curves for the subsequent flutter analysis; Among them, the calculation conditions include the design condition, the choking condition, and the near-stall condition.

[0024] S02. Based on the structural statics and dynamics analysis, obtain the modal shapes and the natural angular frequencies corresponding to at least the first three vibration modes of the blade vibration, specifically: Based on the structural dynamics equation of the blade vibration: , wherein: is the mass matrix, is the damping matrix, is the stiffness matrix, f is the aerodynamic force acting on the blade surface, is the vibration displacement of the blade, is the first derivative of which is the vibration velocity of the blade; is the acceleration of the blade; f In modal analysis, the aerodynamic force and the damping matrix C are not considered, and the following is obtained, Equation 1: Let be Equation 2: , wherein: is the natural angular frequency of the blade, tis the physical time. Substituting Equation (2) into Equation (1) gives Equation (3): , Equation (3) is the characteristic equation. Solving Equation (3) yields the mode shapes corresponding to the vibration modes of the blade and the natural angular frequencies , as well as the modal stress corresponding to the mode shape and the modal displacement .

[0025] S03. Interpolate any one-order mode shape onto the fluid domain computational grid of the blade. The maximum allowable alternating stress satisfying the high-cycle fatigue life of the blade is obtained from the Goodman diagram. Furthermore, the maximum allowable vibration amplitude L of the blade is obtained as follows: , wherein is the modal stress and is the modal displacement.

[0026] S04. Under the computational condition, let the blade perform simple harmonic vibration with any one-order mode shape , natural angular frequency and maximum allowable vibration amplitude . On the premise of considering the non-linearity of the blade caused by excessive amplitude in the simple harmonic vibration of the blade, specify the fundamental frequency of the unsteady flow caused by the blade vibration, which is the natural frequency of the blade: , wherein is the natural angular frequency; Meanwhile, specify the inter-blade phase angle . The inter-blade phase angle is the phase difference when the blades vibrate in a traveling wave mode, and is recorded by the calculation formula as follows: , wherein: is the vibration nodal diameter of the blade and is the number of blades; When is an even number, takes values in [- / 2, / 2), where - / 2 and / 2 differ by . Only one of them is taken during calculation; When is an odd number, The value range is [-( -1) / 2, ( -1) / 2]; In addition, to analyze the unsteady flow caused by blade vibration, the specified harmonic information is as follows: , wherein, is the harmonic number, and is a positive integer greater than 1.

[0027] S05. Initialize the unsteady blade flutter calculation based on the steady calculation results of the blade flow field, and then use the harmonic balance method to analyze the unsteady flow of the blade caused by blade vibration, and obtain the unsteady pressure acting on the blade surface corresponding to the unsteady flow of the blade and the blade vibration velocity . Among them, the formula of the harmonic balance method is: In the formula: is the variable of the unsteady flow field around the blade, is the time-averaged quantity of the unsteady flow field, is the physical time, n is the harmonic number used to analyze the unsteady flow field, i is the order corresponding to different harmonics, and are the Fourier coefficients corresponding to the i th harmonic, is the i th harmonic corresponding angular frequency of the unsteady flow field; the harmonic balance method is used to represent the flow field variable of the unsteady flow as the sum of the time-averaged quantity and the pulsating quantity, and approximate the pulsating quantity with a Fourier series; The expression formula for representing the unsteady pressure acting on the blade surface by the harmonic balance method is: . In the formula: is the time-averaged quantity of the unsteady pressure acting on the blade surface; Similarly, the expression formula for the vibration velocity of the blade obtained is: . In the formula: is the time-averaged quantity of the vibration velocity of the blade.

[0028] S06. Integrate the unsteady pressure acting on the blade surface and the blade vibration velocity within one vibration period to obtain the accumulated work as: , In the formula, is the vibration velocity of the blade, is the normal vector of the fluid domain grid on the blade surface, is the differential of the blade surface area, is the vibration period of the blade, dt is the differential of physical time, and S is the total surface area of the blade.

[0029] S07. Repeat S03 to S06 to calculate the accumulated work corresponding to each vibration mode obtained, and evaluate the blade flutter. Specifically: If the accumulated work is positive, the work done by the air flow on the blade increases within the period, and the amplitude will continue to increase. Then there is a risk of amplitude exceeding the limit, and the fatigue life requirement for the corresponding vibration mode is not met; If the accumulated work is negative and does not exceed the maximum allowable vibration amplitude, there is no risk of amplitude exceeding the limit, and the fatigue life requirement for the corresponding vibration mode is met.

[0030] Example 2: The same as Example 1, except that in S02, based on the structural statics and dynamics analysis, at least the first three modal shapes and the natural angular frequencies are obtained.

[0031] To further illustrate the solution and beneficial effects of the present invention, a specific example with NASARotor 67 is provided: NASA Rotor 67 is a certain stage of the fan rotor in a multi-stage compressor designed and tested by NASA in the late 1980s. The number of blades in one circumference is 22, and the designed rotational speed is 16,043 revolutions per minute.

[0032] As Figures 2 - 4 shown, the non-linear blade flutter evaluation method based on the measurement of the aerodynamic accumulated work of the blade includes the following steps: Step 1. Perform a steady calculation of the blade flow field to obtain the flow rate-pressure ratio and flow rate-efficiency characteristic curves at different rotational speeds, and determine the calculation conditions for blade flutter. Specifically: Generate the computational grid of the fluid domain of the blade, as Figure 2 shown. Solving the FANS equation is required for the flow field calculation. The FANS equation in the cylindrical coordinate system is as follows:

[0033] In the formula: is the conserved variable, is the source term, is the physical time, are the convective fluxes in the axial, circumferential, and radial directions respectively, The viscous fluxes in the axial, circumferential, and radial directions respectively; The velocity components of the computational grid of the fluid domain moving in the axial, circumferential, and radial directions respectively, , , respectively represent the fluxes in the axial, circumferential, and radial directions generated due to the movement of the computational grid of the fluid domain; After the flow field calculation is completed, the characteristic curve of Rotor 67 is obtained; the design rotational speed of Rotor 67 is 16,043 revolutions per minute, its inlet total pressure is the standard atmospheric pressure, the inlet total temperature is 288.15 K, the design point back pressure is 122,000 Pa, the near stall point back pressure is 126,000 Pa, and the near choke point back pressure is 90,000 Pa. Flutter analysis is carried out for each working condition.

[0034] Step 2: Conduct structural statics and dynamics analyses to obtain the modal vibration modes and the natural angular frequencies ; specifically: Let be the mass matrix, be the damping matrix, be the stiffness matrix, f be the aerodynamic force, be the displacement, then based on the structural dynamics equation of blade vibration, there is: , In the formula: is the first derivative of the displacement, which is the velocity; is the second derivative of the displacement, which is the acceleration; In modal analysis, the aerodynamic force f and the damping matrix C are not considered, and we get, Formula 1: , Let be Formula 2: , In the formula: is the natural angular frequency of the blade, is the time. Substituting Formula 2 into Formula 1, we get Formula 3: , Formula 3 is the characteristic equation. Solving Formula 3 can obtain the modal vibration mode and the natural angular frequency , as well as the modal stress corresponding to the modal vibration mode and the modal displacement . After calculation, the first vibration mode of the blade is obtained, and its corresponding natural frequency f is 537.9 Hz; Step 3. Interpolate the calculated modal vibration modes onto the CFD grid on the blade surface. For the vibration modes that are potentially prone to flutter, from the maximum allowable alternating stress that meets the high-cycle fatigue life, the maximum allowable vibration amplitude can be calculated; specifically: , For the given design life, take the maximum value of the alternating stress from the Goodman diagram of Rotor 67 as , and the maximum allowable vibration amplitude is the maximum value of the alternating stress divided by the modal stress multiplied by the modal displacement . The modal displacement and the modal stress are respectively as shown in Figure 3 , Figure 4 . This mode is the first-order bending mode. From the figure, the maximum modal displacement is 7.04356, and the maximum modal stress is . Thus, the maximum allowable vibration amplitude of the blade is .

[0035] Step 4. Let the blade perform simple harmonic vibration with the modal vibration mode , the natural angular frequency and the maximum allowable vibration amplitude . Use the harmonic balance method to calculate the unsteady flow induced by the blade vibration. In the unsteady flow calculation, it is necessary to consider the flow nonlinearity caused by the relatively large vibration amplitude of the blade, that is, in the unsteady flow calculation, it is necessary to analyze the multiple frequencies of the blade vibration frequency. Specifically: Use the harmonic balance method to analyze the unsteady flow induced by the blade vibration. The harmonic balance method can be expressed as: , where: is the variable of the unsteady flow field around the blade, is the time-averaged quantity of the unsteady flow field, is the time, is the number of harmonics used to analyze the unsteady flow field, i is the order corresponding to different harmonics, and are the Fourier coefficients corresponding to the i th harmonic, is the iThe angular frequency of the unsteady flow field corresponding to the subharmonic; the harmonic balance method is used to represent the flow field variables of the unsteady flow as the time-averaged quantity plus the pulsating quantity, and the pulsating quantity is approximated by a Fourier series;

[0036] The fundamental frequency of the unsteady flow caused by the blade vibration is the natural frequency of the blade: , In addition, the inter-blade phase angle also needs to be specified , the inter-blade phase angle characterizes the phase difference when each blade vibrates in a traveling wave mode, that is: , In the formula: is the vibration node diameter, is the number of blades. When is even, takes values in [- / 2, / 2), where - / 2 and / 2 differ by , and only one of them needs to be considered in the calculation; when is odd, takes values in [-( -1) / 2, ( -1) / 2].

[0037] The total number of blades in the full annulus of this example is 22, so its vibration node diameter ranges from [-11, 11). The nonlinearity caused by the blade vibration is relatively strong, and higher-order harmonics also need to be considered in the calculation. If =1 harmonic is used to analyze the unsteady flow caused by the vibration and the 0 node diameter is considered, that is, when the inter-blade phase angle is 0, the harmonic information input is , 0°.

[0038] When performing unsteady calculations, if the nonlinearity of the flow field is relatively strong, it will cause large grid deformation, and the harmonic number independence verification also needs to be carried out.

[0039] Step Five: Based on the unsteady calculation results, calculate the work done by the unsteady aerodynamic force caused by the blade vibration within one vibration period, that is, the accumulated work; Judge whether the blade meets the fatigue life requirements of this mode according to the sign of the accumulated work: the accumulated work being positive means that the work done by the air flow on the blade increases within the period, the amplitude will continue to increase, and there is a risk of amplitude exceeding the limit; the accumulated work being negative means there is no risk of amplitude exceeding the limit; specifically: Accumulated work It refers to the work done by the unsteady aerodynamic force caused by the blade vibration within one vibration period, which is: , When the vibration node diameter is 0, the calculated accumulated work value is -0.2093035 J, and the distribution of the accumulated work density on the blade surface is as Figure 5 shown. The aerodynamic accumulated work on the blade surface is not positive. Therefore, when the blade vibrates in the first-order mode with a vibration node diameter of 0, the blade amplitude will not continue to increase, meeting the fatigue life requirements.

[0040] Step Six: Traverse all possible vibration modes and interblade phase angles to determine whether the blade meets the safety requirements.

[0041] The vibration node diameter range of Rotor 67 is [-11, 11). Traverse all possible vibration modes and interblade phase angles, and determine whether the blade meets the safety requirements according to the positive or negative of the accumulated work.

[0042] From the above steps, the maximum allowable vibration amplitude of the blade can be obtained according to the fatigue life target. Based on this, non-linear analysis can be carried out to determine whether there is a risk of excessive amplitude in blade vibration, solving the problem that the general linear flutter analysis is too conservative.

[0043] Thus, in the numerical analysis of aerodynamic damping by the method of the present invention, the blade performs simple harmonic vibration with the maximum allowable vibration amplitude. The harmonic balance method is used to solve the unsteady flow caused by the blade vibration, and the aerodynamic accumulated work is calculated, considering the flow non-linearity caused by the excessive blade vibration amplitude. According to the positive or negative of the aerodynamic accumulated work, it is determined whether the flutter characteristics of the blade meet the design requirements, completing the non-linear analysis of blade flutter.

[0044] The above are only the preferred embodiments of the present invention and are not intended to limit the present invention. Any modifications, equivalent replacements, improvements, etc. made within the spirit and principle of the present invention shall be included in the protection scope of the present invention.

Claims

1. A non-linear blade flutter evaluation method based on accumulated work measurement and fatigue life, characterized in that It includes the following steps: Step 1: Establish a computational grid for the fluid domain of the blade, perform a steady calculation of the blade flow field, obtain the characteristic curves of flow rate versus pressure ratio and flow rate versus efficiency at different rotational speeds of the blade, so as to determine the calculation conditions for blade flutter; meanwhile, based on structural static and dynamic analyses, obtain the modal shapes corresponding to at least the first three vibration modes of blade vibration and the natural angular frequencies ; Step 2: Interpolate any first-order modal vibration mode onto the fluid domain calculation grid of the blade, and calculate the maximum allowable vibration amplitude of the blade based on the maximum allowable alternating stress that satisfies the high-cycle fatigue life of the blade ; Step 3: Under the computational working conditions, let the blade perform simple harmonic vibration with any first-order modal shape , natural angular frequency and maximum allowable vibration amplitude ; On the premise of considering the nonlinearity of the blade caused by excessive amplitude in the simple harmonic vibration of the blade, the fundamental frequency and harmonic information of the unsteady flow caused by the specified simple harmonic vibration of the blade are specified, and the unsteady blade flutter calculation is initialized based on the steady calculation results of the blade flow field. Thus, the unsteady flow caused by the blade vibration is analyzed by the harmonic balance method to obtain the unsteady pressure received on the blade surface and the blade vibration velocity ; Finally, integrate the unsteady pressure on the blade surface and the blade vibration velocity within one vibration period to obtain the accumulated work; and the blade vibration velocity within one vibration period to obtain the accumulated work; Step 4: Repeat Step 2 to Step 3, calculate the accumulated work corresponding to each order of vibration mode obtained, and evaluate the blade flutter, specifically: If the accumulated work is positive, the work done by the airflow on the blade increases within the period, the amplitude will continue to increase, there is a risk of amplitude exceeding the limit, and the fatigue life requirement of the corresponding vibration mode is not met; If the accumulated work is negative and does not exceed the maximum allowable vibration amplitude, there is no risk of amplitude exceeding the limit, and the fatigue life requirement of the corresponding vibration mode is met.

2. The non-linear blade flutter evaluation method based on accumulated work measurement and fatigue life as described in claim 1, wherein In the steady calculation of the blade flow field, solve the FANS equation for determining the calculation conditions of blade flutter, and express the FANS equation in the cylindrical coordinate system as: In the formula: is the physical time, is the conserved variable per unit volume of the flow field, The expression of is: ， is the fluid density, , , are the velocity components of the fluid in the axial, circumferential, and radial directions, is the total energy of the fluid per unit mass; are the coordinate positions in the axial, circumferential, and radial directions respectively, are the velocity components of the computational grid of the fluid domain moving in the axial, circumferential, and radial directions respectively, , , represent the axial, circumferential, and radial fluxes generated due to the movement of the computational grid of the fluid domain respectively; are the convective fluxes of the fluid in the axial, circumferential, and radial directions, respectively. The convective flux represents the physical quantity transported per unit area per unit time. The expression formula is: ， is the total enthalpy of the fluid per unit mass, is the static pressure of the fluid; are the viscous fluxes of the fluid in the axial, circumferential, and radial directions respectively. The viscous flux is used to describe the diffusive transport of momentum and energy caused by viscous stress and heat conduction, The expression of ， In the formula, is the normal stress between fluids, is the shear stress between fluids, is the heat flux between fluids, , , are respectively the components of the heat flux between fluids in the axial, circumferential, and radial directions; is the source term, which is used to describe the influence on the conserved quantity inside or outside the flow field system. The expression of ， After solving the FANS equation, the characteristic curves of flow rate vs. pressure ratio and flow rate vs. efficiency at different rotational speeds of the blade are obtained, and the calculation conditions are selected according to the characteristic curves for flutter analysis; Among them, the calculation conditions include the design condition, the choking condition, and the near-stall condition.

3. The non-linear blade flutter evaluation method based on accumulated work measurement and fatigue life as claimed in claim 1, wherein Obtaining the modal shapes corresponding to at least the first three vibration modes of the blade vibration and the natural angular frequencies , specifically: Based on the structural dynamics equation of blade vibration: ， In the formula: is the mass matrix, is the damping matrix, is the stiffness matrix, f is the aerodynamic force acting on the blade surface, is the vibration displacement of the blade, is the first derivative of, that is, the vibration velocity of the blade; is the acceleration of the blade; In modal analysis, aerodynamic forces are not considered f and the damping matrix C, to obtain Equation (1): ， Let be Formula 2: ， In the formula: is the natural angular frequency of the blade, t is the physical time. Substituting Formula 2 into Formula 1 gives Formula 3: ， Equation 3 is the characteristic equation. Solving Equation 3 yields the modal shapes corresponding to the various vibration modes of the blade and the natural angular frequencies , as well as the modal stresses corresponding to the modal shapes and the modal displacements .

4. The non-linear blade flutter evaluation method based on accumulated work measurement and fatigue life as claimed in claim 1, wherein The specific content of Step 2 is: Obtain the maximum allowable alternating stress that meets the high-cycle fatigue life of the blade from the Goodman diagram , and then obtain the maximum allowable vibration amplitude of the blade L as follows: ， In the formula, is the modal stress, is the modal displacement.

5. The non-linear blade flutter evaluation method based on accumulated work measurement and fatigue life according to claim 1, characterized in that, In Step 3: Specify the fundamental frequency of the unsteady flow caused by blade vibration, which is the natural frequency of the blade : ， In the formula, is the natural angular frequency; Meanwhile, specify the inter-blade phase angle , the inter-blade phase angle is the phase difference when the blades vibrate in a traveling wave mode, and is denoted by the calculation formula as follows: ， In the formula: is the vibration node diameter of the blade, is the number of blades; When is even, takes values in [- / 2, / 2), where - / 2 and / 2 differ by , and only one of them is taken during calculation; When is odd, the value is in the range of [-( - 1) / 2, ( - 1) / 2]; In addition, to analyze the unsteady flow caused by blade vibration, the specified harmonic information is: ， Among them, is the harmonic number, and is a positive integer greater than 1.

6. The non-linear blade flutter evaluation method based on accumulated work measurement and fatigue life according to claim 1, characterized in that, In step 3, the unsteady blade flutter calculation is initialized based on the steady calculation results of the blade flow field, and then the harmonic balance method is used to analyze the unsteady blade flow caused by blade vibration, so as to obtain the unsteady pressure on the blade surface corresponding to the unsteady blade flow and the blade vibration velocity , where the formula of the harmonic balance method is as follows: ， Wherein: is the variable of the unsteady flow field around the blade, is the time-averaged quantity of the unsteady flow field, is the physical time, n is the harmonic number used to analyze the unsteady flow field, i is the order corresponding to different harmonics, and is the i th Fourier coefficient corresponding to the harmonic, is the i th angular frequency of the unsteady flow field corresponding to the harmonic; the harmonic balance method is used to represent the flow field variables of the unsteady flow as the sum of the time-averaged quantity and the pulsating quantity, and approximate the pulsating quantity with a Fourier series; The unsteady pressure on the blade surface is expressed by the harmonic balance method The expression formula is as follows: ， In the formula: is the time-averaged value of the unsteady pressure on the blade surface; Similarly, the vibration velocity of the obtained blade The expression is as follows: ， In the formula: is the time-averaged value of the vibration velocity of the blade.

7. The non-linear blade flutter evaluation method based on accumulated work measurement and fatigue life as claimed in claim 1, wherein The accumulated work described in Step 3 is as follows: ， In the formula, is the vibration velocity of the blade, is the fluid domain grid normal vector on the blade surface, is the differential of the blade surface area, is the blade vibration period, dt is the differential of physical time, and S is the total blade surface area.

8. The non-linear blade flutter evaluation method based on accumulated work measurement and fatigue life according to any one of claims 1-7, characterized in that Based on the structural static and dynamic analyses, at least the first three modal vibration modes prone to flutter are obtained and the natural angular frequencies .

Images