Fatigue life analysis method of cylindrical roller bearings based on multi-factor coupling

Abstract

A method for analyzing the fatigue life of cylindrical roller bearings based on multi-factor coupling, S1, establishing a thermal-fluid-solid coupling dynamic model of the bearing; S2, establishing a contact subsurface stress calculation model; S3, calculating the fatigue life of the bearing based on the I-H theory. The beneficial effects of the present invention are: different from the prior art that uses a single factor to calculate the fatigue life of the bearing, the present invention is based on the actual working conditions of the bearing, and considers the effects of lubrication, heat and flexible deformation of the ring on the bearing life under actual working conditions, so that the fatigue life obtained by the present invention is more in line with the actual working conditions, providing sufficient theoretical support for the design and production end.

Description

Technical Field

[0001] The invention relates to a bearing fatigue life analysis method, in particular to a cylindrical roller bearing fatigue life analysis method based on multi-factor coupling. Background Art

[0002] Cylindrical thin-walled roller bearings in modern high-end equipment generally include rolling elements, outer rings, inner rings and other components, and generally work under high-speed, heavy-load, high-temperature and harsh working conditions such as aerospace, wind power generation, and industrial gearboxes. Specifically for the above fields, fatigue life analysis of cylindrical thin-walled roller bearings must consider the coupling effects of lubrication, heat and flexible deformation of thin-walled rings under complex and harsh working conditions.

[0003] At present, there are many basic theories for analyzing bearing fatigue life. In the bearing fatigue life analysis theory, the LP fatigue life analysis method has provided important development for the analysis of bearing fatigue life and has been widely used. However, the LP theory has certain limitations, such as the inability to deal with infinite life problems. Taking into account the limitations of the LP theory, researchers have proposed an IH theoretical model that can consider the stress limit value based on the LP theory. Based on the IH theory, there are currently many literatures that provide analysis methods for studying bearing fatigue life. For example, some scholars have studied the influence of friction effects on bearing life, and some scholars have proposed a fatigue life analysis method based on damage accumulation theory.

[0004] However, under complex and harsh working conditions, the rheological properties of the lubricating oil inside the bearing, the heat generation and deformation of the bearing, and the flexible deformation of the thin-walled bearing ring are coupled with each other. All current analysis methods cannot analyze the influence of the coupling characteristics of lubrication-heat-flexible deformation on the fatigue life of the bearing. In order to more reasonably analyze the fatigue life of the bearing, a bearing fatigue life analysis method is needed that can comprehensively consider the complex coupling relationship between lubrication, heat and flexible deformation. Summary of the invention

[0005] In view of the deficiencies in the prior art, the present invention provides a method for analyzing the fatigue life of a cylindrical roller bearing taking into account lubrication-heat-flexible deformation coupling.

[0006] To achieve the above purpose, the technical solution of the present invention is as follows: a method for analyzing the fatigue life of cylindrical roller bearings based on multi-factor coupling,

[0007] S1, establish the bearing thermal-fluid-solid coupling dynamics model. The specific steps are as follows:

[0008] S101, establishing a ferrule flexible body model;

[0009] S102, establishing a bearing lubrication model considering nonlinear viscoelastic characteristics;

[0010] S103, establishing a bearing thermal network model taking into account lubricating oil temperature;

[0011] S104, establishing a bearing dynamics model based on the interaction relationship between various bearing components;

[0012] S105, coupling the four models in S101-S104 is realized by parameter coupling, and a thermal-fluid-solid coupling dynamic model of the bearing is established;

[0013] S2, establish a contact subsurface stress calculation model;

[0014] S3, calculates the fatigue life of the bearing based on the IH theory.

[0015] Furthermore, in step S101, the ferrule flexible body model is established using curved beam elements based on the finite element method.

[0016] Furthermore, in step S102, the constitutive model of the nonlinear viscoelastic fluid characteristics is expressed as follows:

[0017]

[0018] In the above formula, τ is the shear stress, (U s is the sliding speed between the rolling element and the ring, h is the center oil film thickness), G is the shear modulus, τ0 is the yield shear stress, η is the viscosity of the lubricating oil,

[0019] The calculation formula of the center oil film thickness h is as follows:

[0020] h=R y 3.05U 0.69 g 0.56 W -0.1

[0021] Where U is the dimensionless velocity parameter, U = η0u / (E eq R y );

[0022] W is the dimensionless load, W = Q / (E eq R y l);

[0023] g is a dimensionless material parameter, g = λE eq (In the above expression, Q is the contact load, η0 is the reference viscosity, u is the rolling speed, E eq is the equivalent elastic modulus, l is the rolling element length, λ is the viscosity-pressure coefficient of the lubricating oil, R y is the equivalent radius in the rolling direction).

[0024] Furthermore, in step S102, in view of the influence of pressure p and temperature T on shear modulus G, yield shear stress τ0 and lubricating oil viscosity η, the modified calculation formulas of shear modulus G, yield shear stress τ0 and lubricating oil viscosity η are as follows:

[0025]

[0026] In the above formula, η0, are the reference viscosity, reference shear modulus and reference yield shear stress respectively; α i , β j , γ k , α G , β G , α τ and β τ is a specific constant; T and T0 are the current temperature and reference temperature respectively; p is the current pressure.

[0027] Further, according to step S102, the friction traction curve is obtained, and the specific steps are as follows:

[0028] S1021, the constitutive model of nonlinear viscoelastic fluid characteristics uses the Runge-Kutta method to obtain the shear stress τ at the contact surface between the rolling element and the ring;

[0029] S1022, integrating the shear stress τ over the contact area to obtain the friction traction force;

[0030] S1023, comparing the friction traction force and the contact load Q to obtain a friction traction coefficient;

[0031] S1024, by changing the temperature T and sliding speed U s Different values ​​of T and U s The friction traction curve below.

[0032] Furthermore, in the step S1024, a BP neural network is used to perform function fitting on the friction traction curve so that the friction traction curve can be directly used in the dynamic model later.

[0033] Furthermore, in step S103, the calculation formula for the thermal deformation of each bearing component is obtained based on the temperature T of each bearing component as follows:

[0034] u t =ΓL(TT α )

[0035] In the above formula, Γ is the thermal expansion coefficient, T α is the ambient temperature, and L is the characteristic length of each bearing element.

[0036] Furthermore, in step S103, a temperature node is set at each bearing element, and the following calculation formula is listed based on the law of conservation of energy:

[0037] Node 1: (T1-T 10 ) / R 1cv10 +(T1-T2) / R 2cc1 =0

[0038] Node 2: (T2-T1) / R 1cc2 +(T2-T3) / R 3cc2 =0

[0039] Node 3: (T3-T2) / R 2cc3 +(T3-T4) / R 4cc3 +(T3-T0) / R 0cv3 =Q1

[0040] Node 4: (T4-T5) / R 5cc4 +(T4-T3) / R 4cc3 =0

[0041] Node 5: (T5-T4) / R 5cc4 +(T5-T6) / R 6cc5 +(T5-T0) / R 0cv5 =Q2

[0042] Node 6: (T6-T5) / R 6cc5 +(T6-T7) / R 7cc6 =0

[0043] Node 7: (T7-T6) / R 7cc6 +(T7-T8) / R 8cc7 +(T7-T0) / R 7cv0 =Q3

[0044] Node 8: (T8-T7) / R 8cc7 +(T8-T9) / R 9cc8 =0

[0045] Node 9: (T9-T8) / R 9cc8 +(T9-T 10 ) / R 10cv9 =0

[0046] In the above formula, T1 is the outer ring temperature node, T2 is the outer raceway temperature node, T3 is the outer raceway lubricant temperature node, T4 is the rolling element-outer raceway outer node, T5 is the rolling element temperature node, T6 is the rolling element-inner raceway inner node, T7 is the inner raceway lubricant temperature node, T8 is the inner raceway temperature node, T9 is the inner ring temperature node, T 10is the air temperature node, T0 is the lubricant temperature node; Q1, Q2, Q3 are the friction heat generation rate at the contact between the rolling element and the outer raceway, the friction heat generation rate between the rolling element and the lubricating oil in the bearing cavity, and the friction heat generation rate at the contact between the rolling element and the inner raceway, respectively; R is the thermal resistance value (in the subscript, cv represents heat convection, cc represents heat conduction, and Arabic numerals represent temperature nodes).

[0047] Furthermore, in the step S104, a GUPTA bearing model is used to establish a dynamic model of the bearing.

[0048] Furthermore, the step S105 establishes a bearing thermal-fluid-solid coupling dynamics model according to the following steps:

[0049] S1051, initializing the temperature, position and speed of each bearing component;

[0050] S1052, calculating the current friction traction coefficient based on the bearing lubrication model in step S102;

[0051] S1053, calculating the geometric approach and contact load between the bearing components and the friction heat generation rate based on the GUPTA dynamic model in step S104;

[0052] S1054, substituting the friction heat generation rate obtained in step S1053 into the bearing thermal network model in step S103 to update the temperature of each bearing component, and obtain the thermal deformation of each bearing component;

[0053] S1055, updating the geometric dimensions of each bearing component based on the elastic deformation of each bearing component calculated by the ring flexible body model in S101 and the thermal deformation of each bearing component calculated in step S1054;

[0054] S1056, substituting the updated temperature obtained in step S1054 into the bearing lubrication model in step S102 to update the viscosity and traction coefficient of the lubricating oil;

[0055] S1057, substituting the geometric dimensions and traction coefficients updated in the above steps S1055 and S1056 into the GUPTA model to update the position and speed of the bearing element;

[0056] S1058, repeatedly iterate the above steps S1051-S1057 until the termination condition is met and then exit, and the termination condition may be that the calculation time is reached.

[0057] Furthermore, in the step S2, discrete convolution (DC) and fast Fourier transform (FFT) methods are used to calculate the sub-surface stress of the contact area.

[0058] Furthermore, the step S2 obtains the sub-surface stress of the contact area according to the following steps:

[0059] S201, calculate the contact half-width b, and discretize the entire contact area into N along the rolling element length direction (x direction), contact half-width direction (y direction) and depth direction (z direction) x 、N y and N z points;

[0060] S202, calculate the influence coefficient sequence {D} N (symbol{} N represents a sequence of length N), and through zero padding and wrapping, the sequence {D} N Expands to { D} 2N (symbol{} 2N represents a sequence of length 2N);

[0061] S203, for the sequence { D} 2N Perform fast Fourier transform and get

[0062] S204, calculate the pressure sequence {P} N ; and by padding with zeros, {P} N Expands to { P} 2N ;

[0063] S205, for the sequence { P} 2N Using Fourier transform, we get the sequence

[0064] S206, the sequence in step S203 and step S205 Multiply each unit to get the frequency

[0065] S207, the frequency sequence in step S206 Using inverse Fourier transform, we get the sequence {V} 2N ;

[0066] S208, the sequence {V} obtained from step S206 2N Extract the nth (n∈[0,N-1]) element and obtain the stress distribution.

[0067] Furthermore, the calculation formula of the influence coefficient in step S202 is as follows:

[0068]

[0069] In the above formula, x + =x m -x ξ +Δ1 / 2,x - =x m -x ξ -Δ1 / 2,y + =y n -y η +Δ2 / 2,y - =y n -y η -Δ2 / 2, ξ and η are the local coordinates of the discrete unit, m, n, k are the element numbers of the unit in the x, y, z directions respectively, q, r = x, y or z (i.e. qr is a combination of x, y, z in pairs).

[0070] Furthermore, for the contact in the elastic half-space, the stress field adopts the function T Nqr The calculation formula is as follows:

[0071]

[0072] T Nyy (x,y,z)=T Nxx (y,x,z)

[0073]

[0074] In the above formula, When z = 0, (sign() is the sign function).

[0075] Furthermore, the step S3 calculates the fatigue life of the bearing according to the following steps:

[0076] S301, calculating the LP life of the rolling element-rotating ring raceway contact and the rolling element-stationary ring raceway contact respectively based on the LP;

[0077] S302, according to the data obtained in step S301, the IH life of each rolling element at the contact point with the rotating ring and the stationary ring is obtained;

[0078] S303, obtaining the fatigue life of the bearing according to the data obtained in step S302.

[0079] Furthermore, in the step S301, for the rotating ring, the calculation formula of the LP life at the rolling element-raceway contact is as follows:

[0080]

[0081] For the stationary ring, the LP life at the rolling element-raceway contact is calculated as follows:

[0082]

[0083] In the above formula, q μj and q vj is the contact load between the rotating ring and the stationary ring (the subscripts μ and v represent the rotating ring and the stationary ring respectively, and j represents the jth rolling element); q cμj and q cvj is the rated dynamic load of the rotating ring and the stationary ring; where q cμj and q cvj The calculation formula is as follows:

[0084]

[0085] In the above formula, D is the rolling element diameter, d m is the bearing pitch diameter, α is the contact angle, l is the length of the rolling element, and Z is the number of rolling elements.

[0086] Further, according to the formula in step S301, the IH life L of the contact point between the j-th rolling element and the rotating ring and the stationary ring in step S302 is obtained. IHμj and L IHvj , and its calculation formula is as follows:

[0087]

[0088] In the above formula, σ VMμj,max and σ VMvj,max are the maximum Von-Mises stress at the contact point between the jth rolling element and the rotating ring and the stationary ring, respectively, VM,lim is the limiting value of Von-Mises stress. and These are fixed values ​​and vary only with the bearing type.

[0089] Furthermore, the calculation formula in step S303 is as follows:

[0090]

[0091] The above formula is used to superimpose the IH life of each rolling element obtained in step S302 at the contact point with the rotating ring and the stationary ring to obtain the fatigue life of the entire bearing.

[0092] The beneficial effects of the present invention are as follows: Different from the prior art that uses a single factor to calculate the fatigue life of the bearing, the present invention is based on the actual working conditions of the bearing, and takes into account the influence of lubrication, heat and flexible deformation of the ring on the bearing life under the actual working conditions, so that the fatigue life obtained by the present invention is more in line with the actual working conditions, providing sufficient theoretical support for the design and production end;

[0093] Secondly, coupling lubrication, heat and flexible deformation of the ring through parameters can not only reduce the difficulty of calculation, but also reduce the computing power cost required for calculation, thereby improving the calculation efficiency; and by comprehensively considering the impact of multiple factors on the bearing life, and considering the interaction between multiple factors as a correction for each independent factor, the fatigue life data obtained by the present invention is further closer to the actual working conditions. BRIEF DESCRIPTION OF THE DRAWINGS

[0094] Figure 1 A schematic diagram of a flow chart of an embodiment of the present invention;

[0095] Figure 2 It is a schematic diagram of the overall composition of the thermal-fluid-solid coupling dynamics model of a bearing according to an embodiment of the present invention;

[0096] Figure 3 is a friction traction curve obtained based on step S102 in an embodiment of the present invention;

[0097] Figure 4 A BP neural network model for fitting the friction traction curve in an embodiment of the present invention;

[0098] Figure 5 is the contact load between the rolling element and the inner ring (rotating ring) at different temperatures in the embodiment of the present invention;

[0099] Figure 6 is the contact load between the rolling element and the outer ring (stationary ring) at different temperatures in the embodiment of the present invention;

[0100] Figure 7 is the bearing life when the temperature is T=300K in the embodiment of the present invention;

[0101] Figure 8 is the bearing life when the temperature is T=310K in the embodiment of the present invention. DETAILED DESCRIPTION

[0102] The present invention provides a method for analyzing the fatigue life of a cylindrical roller bearing based on multi-factor coupling. Figure 1-8 As shown: In this embodiment, the bearing is a cylindrical roller bearing, and its inner ring is a rotating ring, and its outer ring is a stationary ring. The following describes the calculation of the fatigue life of the cylindrical roller bearing according to the analysis method in this embodiment.

[0103] S1, establish the bearing thermal-fluid-solid coupling dynamic model, and its fixing steps are as follows:

[0104] S101, establishing a flexible body model of the ferrule based on the finite element method;

[0105] The flexible ring is modeled using two-dimensional curved beam elements.

[0106] S102, establishing a bearing lubrication model that can take into account nonlinear viscoelastic characteristics;

[0107] The constitutive model of the nonlinear viscoelastic fluid used can be expressed as:

[0108]

[0109] In the above formula, τ is the shear stress, (U s is the sliding speed between the rolling element and the ring, h is the center oil film thickness), G is the shear modulus, τ0 is the yield shear stress, η is the viscosity of the lubricating oil,

[0110] The calculation formula of the center oil film thickness h is as follows:

[0111] h=R y 3.05U 0.69 g 0.56 W -0.1

[0112] Where U is the dimensionless velocity parameter, U = η0u / (E eq R y );

[0113] W is the dimensionless load, W = Q / (E eq R y l);

[0114] g is a dimensionless material parameter, g = λE eq (In the above expression, Q is the contact load, η0 is the reference viscosity, u is the rolling speed, E eq is the equivalent elastic modulus, l is the rolling element length, λ is the viscosity-pressure coefficient of the lubricating oil, R y is the equivalent radius in the rolling direction). Here, the influence of pressure and temperature on fluid parameters is more obvious, so it is necessary to correct the shear modulus G, yield shear stress τ0 and lubricating oil viscosity η according to the pressure p and temperature T. The calculation formula is as follows:

[0115]

[0116] In the above formula, η0, are the reference viscosity, reference shear modulus and reference yield shear stress respectively; αi , β j , γ k , α G , β G , α τ and β τ is a specific constant; T and T0 are the current temperature and reference temperature respectively; p is the current pressure.

[0117] According to step S102, a friction traction curve can be obtained, and the specific steps are as follows:

[0118] S1021, the constitutive model of nonlinear viscoelastic fluid characteristics is numerically integrated using the Runge-Kutta method, and the shear stress τ at the contact surface between the rolling element and the ring is obtained after integration;

[0119] S1022, integrating the shear stress τ over the contact area to obtain the friction traction force;

[0120] S1023, comparing the friction traction force and the contact load Q to obtain a friction traction coefficient;

[0121] S1024, by changing the temperature T and sliding speed U s Different values ​​of T and U s The friction traction curve below;

[0122] S1025, after using BP neural network to perform function fitting on the friction traction curve, the friction traction curve can be directly used in the dynamic model.

[0123] S103, establishing a bearing thermal network model that can take into account the lubricating oil temperature, calculating the internal temperature distribution of the bearing and the thermal deformation of each bearing component,

[0124] The calculation process of the internal temperature of the bearing is as follows: set temperature nodes at each key position of the rolling bearing, assuming that the air temperature and the temperature inside the bearing cavity are known, then according to the law of conservation of energy (the heat flowing in and out of each node is equal), 9 heat balance equations can be listed: Node 1: (T1-T 10 ) / R 1cv10 +(T1-T2) / R 2cc1 =0

[0125] Node 2: (T2-T1) / R 1cc2 +(T2-T3) / R 3cc2 =0

[0126] Node 3: (T3-T2) / R 2cc3 +(T3-T4) / R 4cc3 +(T3-T0) / R 0cv3 =Q1 node 4:(T4-T5) / R 5cc4+(T4-T3) / R 4cc3 =0

[0127] Node 5: (T5-T4) / R 5cc4 +(T5-T6) / R 6cc5 +(T5-T0) / R 0cv5 =Q2 node 6:(T6-T5) / R 6cc5 +(T6-T7) / R 7cc6 =0

[0128] Node 7: (T7-T6) / R 7cc6 +(T7-T8) / R 8cc7 +(T7-T0) / R 7cv0 =Q3 node 8:(T8-T7) / R 8cc7 +(T8-T9) / R 9cc8 =0

[0129] Node 9: (T9-T8) / R 9cc8 +(T9-T 10 ) / R 10cv9 =0

[0130] In the above formula, T1 is the outer ring temperature node, T2 is the outer raceway temperature node, T3 is the outer raceway lubricant temperature node, T4 is the rolling element-outer raceway outer node, T5 is the rolling element temperature node, T6 is the rolling element-inner raceway inner node, T7 is the inner raceway lubricant temperature node, T8 is the inner raceway temperature node, T9 is the inner ring temperature node, T 10 is the air temperature node, T0 is the lubricant temperature node; Q1, Q2, Q3 are the friction heat generation rate at the contact between the rolling element and the outer raceway, the friction heat generation rate between the rolling element and the lubricating oil in the bearing cavity, and the friction heat generation rate at the contact between the rolling element and the inner raceway, respectively; R is the thermal resistance value (in the subscript, cv represents heat convection, cc represents heat conduction, and Arabic numerals represent temperature nodes).

[0131] After calculating the temperature T of each bearing component, the calculation formula for the thermal deformation of the bearing component is:

[0132] u t =ΓL(TT α )

[0133] In the above formula, Γ is the thermal expansion coefficient, T α is the ambient temperature, and L is the characteristic length of each bearing element.

[0134] S104, considering the interaction relationship between the various bearing components, a dynamic model of the cylindrical roller bearing is established.

[0135] The dynamic model of cylindrical roller bearing is established based on the GUPTA bearing model.

[0136] S105, coupling the models in the above four steps S101-S104 is achieved by parameter coupling, a thermal-fluid-solid coupling dynamic model of a cylindrical roller bearing is established, and the contact loads between the rotating ring and the stationary ring and the rolling element are calculated respectively. The specific steps are as follows:

[0137] S1051, initializing the temperature, position and speed of each bearing component;

[0138] S1052, calculating the current friction traction coefficient based on the bearing lubrication model in step S102;

[0139] S1053, calculating the geometric approach and contact load between the bearing components and the friction heat generation rate based on the GUPTA dynamic model in step S104;

[0140] S1054, substituting the friction heat generation rate obtained in step S1053 into the bearing thermal network model in step S103 to update the temperature of each bearing component, and obtain the thermal deformation of each bearing component;

[0141] S1055, updating the geometric dimensions of each bearing component based on the elastic deformation of each bearing component calculated by the ring flexible body model in S101 and the thermal deformation of each bearing component calculated in step S1054;

[0142] S1056, substituting the updated temperature obtained in step S1054 into the bearing lubrication model in step S102 to update the viscosity and traction coefficient of the lubricating oil;

[0143] S1057, substituting the geometric dimensions and traction coefficients updated in the above steps S1055 and S1056 into the GUPTA model to update the position and speed of the bearing element;

[0144] S1058, repeatedly iterate the above steps S1051-S1057 until the termination condition is met (the calculation time is reached).

[0145] The above steps complete the establishment of the bearing thermal-fluid-solid coupling dynamic model.

[0146] S2, establish a contact subsurface stress calculation model, and calculate the contact subsurface stress according to the contact load calculated in S1. The specific steps are as follows:

[0147] S201, calculate the contact half-width b, and discretize the entire contact area into N along the rolling element length direction (x direction), contact half-width direction (y direction) and depth direction (z direction) x 、N y and N z points;

[0148] S202, calculate the influence coefficient sequence {D} N (symbol{} N represents a sequence of length N), and through zero padding and wrapping, the sequence {D} N Expands to { D} 2N (symbol{} 2N represents a sequence of length 2N);

[0149] S203, for the sequence { D} 2N Perform fast Fourier transform and get

[0150] S204, calculate the pressure sequence {P} N ; and by padding with zeros, {P} N Expands to { P} 2N ;

[0151] S205, for the sequence { P} 2N Using Fourier transform, we get the sequence

[0152] S206, the sequence in step S203 and step S205 Multiply each unit to get the frequency

[0153] S207, the frequency sequence in step S206 Using inverse Fourier transform, we get the sequence {V} 2N ;

[0154] S208, the sequence {V} obtained from step S206 2N Extract the nth (n∈[0,N-1]) element number and obtain the stress distribution.

[0155] The calculation formula of the influence coefficient in step S202 is as follows:

[0156]

[0157] In the above formula, x + =x m -x ξ +Δ1 / 2,x - =x m -x ξ -Δ1 / 2,y + =y n -y η+Δ2 / 2,y - =y n -y η -Δ2 / 2, ξ and η are the local coordinates of the discrete element, m, n, k are the element numbers of the element in the x, y, z directions, respectively, q, r = x, y or z (i.e., qr is a combination of x, y, z). For the contact in the elastic half space, the stress field uses the function T Nqr The calculation formula is as follows:

[0158]

[0159] T Nyy (x,y,z)=T Nxx}y,x,z)

[0160]

[0161] In the above formula, When z = 0, (sign() is the sign function).

[0162] The above steps complete the contact subsurface stress calculation model.

[0163] S3, calculate the fatigue life of the bearing based on the IH theory, the specific steps are as follows:

[0164] S301, calculating the LP life of the rolling element-rotating ring raceway contact and the rolling element-stationary ring raceway contact respectively based on the LP;

[0165] For rotating rings, the LP life at the rolling element-raceway contact is calculated as follows:

[0166]

[0167] For the stationary ring, the LP life at the rolling element-raceway contact is calculated as follows:

[0168]

[0169] In the above formula, q μj and q vj is the contact load between the rotating ring and the stationary ring (the subscripts μ and v represent the rotating ring and the stationary ring respectively, and j represents the jth rolling element); q cμj and q cvj is the rated dynamic load of the rotating ring and the stationary ring; where q cμj and q cvj The calculation formula is as follows:

[0170]

[0171] In the above formula, D is the rolling element diameter, d m is the bearing pitch diameter, α is the contact angle, l is the length of the rolling element, and Z is the number of rolling elements.

[0172] S302, according to the data obtained in step S301, the IH life L of the contact points between the j-th rolling element and the rotating ring and the stationary ring is obtained. IHμj and L IHvj , and its calculation formula is as follows:

[0173]

[0174] In the above formula, σ VMμj,max and σ VMvj,max are the maximum Von-Mises stress at the contact point between the jth rolling element and the rotating ring and the stationary ring, respectively, VM,lim is the limiting value of Von-Mises stress. and are fixed values, which only change according to the bearing type. In this embodiment, the coefficients of the cylindrical roller bearing are as follows:

[0175]

[0176] S303, according to the data obtained in step S302, the fatigue life of the bearing is obtained, and the calculation is shown as follows:

[0177]

[0178] The above formula is used to superimpose the IH life of each rolling element obtained in step S302 at the contact point with the rotating ring and the stationary ring to obtain the fatigue life of the entire bearing.

[0179] The above steps complete the calculation of bearing fatigue life.

[0180] The above embodiment is only one of the preferred specific embodiments of the present invention. Common changes and substitutions made by those skilled in the art within the scope of the technical solution of the present invention are all included in the protection scope of the present invention.

Claims

1. A method for analyzing the fatigue life of cylindrical roller bearings based on multi-factor coupling, characterized in that: S1, establish the bearing thermal-fluid-solid coupling dynamics model. The specific steps are as follows: S101, establishing a ferrule flexible body model; S102, establishing a bearing lubrication model considering nonlinear viscoelastic characteristics; S103, establishing a bearing thermal network model taking into account lubricating oil temperature; S104, based on the interaction relationship between the various components of the bearing, a dynamic model of the bearing is established using a GUPTA bearing model; S105, coupling the four models in S101-S104 is realized by parameter coupling, and a thermal-fluid-solid coupling dynamic model of the bearing is established; The step S105 establishes a bearing thermal-fluid-solid coupling dynamics model according to the following steps: S1051, initializing the temperature, position and speed of each bearing component; S1052, calculating the current friction traction coefficient based on the bearing lubrication model in step S102; S1053, calculating the geometric approach and contact load between the bearing components and the friction heat generation rate based on the GUPTA dynamic model in step S104; S1054, substituting the friction heat generation rate obtained in step S1053 into the bearing thermal network model in step S103 to update the temperature of each bearing component, and obtain the thermal deformation of each bearing component; S1055, updating the geometric dimensions of each bearing component based on the elastic deformation of each bearing component calculated by the ring flexible body model in S101 and the thermal deformation of each bearing component calculated in step S1054; S1056, substituting the updated temperature obtained in step S1054 into the bearing lubrication model in step S102 to update the viscosity and traction coefficient of the lubricating oil; S1057, substituting the geometric dimensions and traction coefficients updated in the above steps S1055 and S1056 into the GUPTA model to update the position and speed of the bearing element; S1058, repeatedly iterate the above steps S1051-S1057 until the termination condition is met and then exit; S2, establish a contact subsurface stress calculation model; S3, calculates the fatigue life of the bearing based on the IH theory.

2. The method for analyzing the fatigue life of cylindrical roller bearings based on multi-factor coupling according to claim 1 is characterized in that: In the step S101, the flexible body model of the ferrule is established by using curved beam elements based on the finite element method.

3. The method for analyzing the fatigue life of cylindrical roller bearings based on multi-factor coupling according to claim 1 is characterized in that: In step S102, the constitutive model of nonlinear viscoelastic fluid characteristics is expressed as follows: In the above formula, τ is the shear stress, Middle U s is the sliding velocity between the rolling element and the ring, h is the central oil film thickness, G is the shear modulus, τ0 is the yield shear stress, η is the viscosity of the lubricating oil, The calculation formula of the center oil film thickness h is as follows: h=R y 3.05U 0.69 g 0.56 W -0.1 Where U is the dimensionless velocity parameter, U = η0u / (E eq R y ); W is the dimensionless load, W = Q / (E eq R y l); g is a dimensionless material parameter, g = λE eq In the above expression, Q is the contact load, η0 is the reference viscosity, u is the rolling speed, E eq is the equivalent elastic modulus, l is the rolling element length, λ is the viscosity-pressure coefficient of the lubricating oil, R y is the equivalent radius in the rolling direction.

4. The method for analyzing the fatigue life of cylindrical roller bearings based on multi-factor coupling according to claim 3 is characterized in that: In step S102, considering the influence of pressure p and temperature T on shear modulus G, yield shear stress τ0 and lubricating oil viscosity η, the modified calculation formulas of shear modulus G, yield shear stress τ0 and lubricating oil viscosity η are as follows: In the above formula, η0, are the reference viscosity, reference shear modulus and reference yield shear stress respectively; α i , β j , γ k , α G , β G , α τ and β τ is a specific constant; T and T0 are the current temperature and reference temperature respectively; p is the current pressure.

5. The method for analyzing the fatigue life of cylindrical roller bearings based on multi-factor coupling according to claim 4 is characterized in that: The friction traction curve is obtained according to step S102, and the specific steps are as follows: S1021, the constitutive model of nonlinear viscoelastic fluid characteristics uses the Runge-Kutta method to obtain the shear stress τ at the contact surface between the rolling element and the ring; S1022, integrating the shear stress τ over the contact area to obtain the friction traction force; S1023, comparing the friction traction force and the contact load Q to obtain a friction traction coefficient; S1024, by changing the temperature T and sliding speed U s Different values ​​of T and U s The friction traction curve below.

6. The method for analyzing the fatigue life of cylindrical roller bearings based on multi-factor coupling according to claim 5 is characterized in that: In the step S1024, a BP neural network is used to perform function fitting on the friction traction curve.

7. The method for analyzing the fatigue life of cylindrical roller bearings based on multi-factor coupling according to claim 1 is characterized in that: In step S103, the calculation formula for the thermal deformation of each bearing component is obtained based on the temperature T of each bearing component as follows: u t =ΓL(TT α ) In the above formula, Γ is the thermal expansion coefficient, T α is the ambient temperature, and L is the characteristic length of each bearing element.

8. The method for analyzing the fatigue life of cylindrical roller bearings based on multi-factor coupling according to claim 7 is characterized in that: In step S103, temperature nodes are set at various components of the bearing, and the following calculation formula is listed based on the law of conservation of energy: Node 1: (T1-T 10 ) / R 1cv10 +(T1-T2) / R 2cc1 =0 Node 2: (T2-T1) / R 1cc2 +(T2-T3) / R 3cc2 =0 Node 3: (T3-T2) / R 2cc3 +(T3-T4) / R 4cc3 +(T3-T0) / R 0cv3 =Q1 Node 4: (T4-T5) / R 5cc4 +(T4-T3) / R 4cc3 =0 Node 5: (T5-T4) / R 5cc4 +(T5-T6) / R 6cc5 +(T5-T0) / R 0cv5 =Q2 Node 6: (T6-T5) / R 6cc5 +(T6-T7) / R 7cc6 =0 Node 7: (T7-T6) / R 7cc6 +(T7-T8) / R 8cc7 +(T7-T0) / R 7cv0 =Q3 Node 8: (T8-T7) / R 8cc7 +(T8-T9) / R 9cc8 =0 Node 9: (T9-T8) / R 9cc8 +(T9-T 10 ) / R 10cv9 =0 In the above formula, T1 is the outer ring temperature node, T2 is the outer raceway temperature node, T3 is the outer raceway lubricant temperature node, T4 is the rolling element-outer raceway outer node, T5 is the rolling element temperature node, T6 is the rolling element-inner raceway inner node, T7 is the inner raceway lubricant temperature node, T8 is the inner raceway temperature node, T9 is the inner ring temperature node, T 10 is the air temperature node, T0 is the lubricant temperature node; Q1, Q2, Q3 are the friction heat generation rate at the contact between the rolling element and the outer raceway, the friction heat generation rate between the rolling element and the lubricating oil in the bearing cavity, and the friction heat generation rate at the contact between the rolling element and the inner raceway, respectively; R is the thermal resistance value. In the subscripts of the above thermal resistance values, cv represents heat convection, cc represents heat conduction, and Arabic numerals represent temperature nodes.

9. The method for analyzing the fatigue life of cylindrical roller bearings based on multi-factor coupling according to claim 1, characterized in that: In the step S2, the subsurface stress of the contact area is calculated using discrete convolution and fast Fourier transform methods.

10. The method for analyzing the fatigue life of cylindrical roller bearings based on multi-factor coupling according to claim 9, characterized in that: The S2 step obtains the subsurface stress in the contact area according to the following steps: S201, calculate the contact half width b, and discretize the entire contact area into N along the rolling element length direction, i.e., the x direction, the contact half width direction, i.e., the y direction, and the depth direction, i.e., the z direction. x 、N y and N z points; S202, calculate the influence coefficient sequence {D} N , symbol {} N represents a sequence of length N, and through zero padding and wrapping, the sequence {D} N Expands to { D } 2N , symbol {} 2N Represents a sequence of length 2N; S203, for the sequence { D } 2N Perform fast Fourier transform and get S204, calculate the pressure sequence {P} N ; and by padding with zeros, {P} N Expand to S205, for the sequence { P } 2N Using Fourier transform, we get the sequence S206, the sequence in step S203 and step S205 Multiply each unit to get the frequency S207, the frequency sequence in step S206 Using inverse Fourier transform, we get the sequence {V} 2N ; S208, the sequence {V} obtained from step S206 2N Extract the nth (n∈[0,N-1]) element number and obtain the stress distribution.

11. The method for analyzing the fatigue life of cylindrical roller bearings based on multi-factor coupling according to claim 10, characterized in that: The calculation formula of the influence coefficient in step S202 is as follows: In the above formula, x + =x m -x ξ +Δ1 / 2,x - =x m -x ξ -Δ1 / 2,y + =y n -y η +Δ2 / 2,y - =y n -y η -Δ2 / 2, ξ and η are the local coordinates of the discrete unit, m, n, k are the element numbers of the unit in the x, y, z directions respectively, q, r = x, y or z, that is, qr is a combination of x, y, z.

12. The method for analyzing the fatigue life of cylindrical roller bearings based on multi-factor coupling according to claim 11, characterized in that: For contact in an elastic half-space, the stress field is given by the function T Nqr The calculation formula is as follows: T Nyy (x,y,z)=T Nxx (y,x,z) In the above formula, When z = 0, , sign() is the sign function.

13. The method for analyzing the fatigue life of cylindrical roller bearings based on multi-factor coupling according to claim 1, characterized in that: The S3 step calculates the fatigue life of the bearing according to the following steps: S301, calculating the LP life of the rolling element-rotating ring raceway contact and the rolling element-stationary ring raceway contact respectively based on the LP; S302, according to the data obtained in step S301, the IH life of each rolling element at the contact point with the rotating ring and the stationary ring is obtained; S303, obtaining the fatigue life of the bearing according to the data obtained in step S302.

14. The method for analyzing the fatigue life of cylindrical roller bearings based on multi-factor coupling according to claim 13, characterized in that: In the step S301, for the rotating ring, the calculation formula for the LP life at the rolling element-raceway contact is as follows: For the stationary ring, the LP life at the rolling element-raceway contact is calculated as follows: In the above formula, q μj and q vj is the contact load between the rotating ring and the stationary ring, the subscripts μ and v represent the rotating ring and the stationary ring respectively, and j represents the jth rolling element; q cμj and q cvj is the rated dynamic load of the rotating ring and the stationary ring; where q cμj and q cvj The calculation formula is as follows: In the above formula, D is the rolling element diameter, d m is the bearing pitch diameter, α is the contact angle, l is the rolling element length, and Z is the number of rolling elements.

15. The method for analyzing the fatigue life of cylindrical roller bearings based on multi-factor coupling according to claim 14, characterized in that: According to the formula in step S301, the IH life L of the contact between the j-th rolling element and the rotating ring and the stationary ring in step S302 is obtained. IHμj and L IHvj , and its calculation formula is as follows: In the above formula, σ VMμj,max and σ VMvj,max are the maximum Von-Mises stress at the contact point between the jth rolling element and the rotating ring and the stationary ring, respectively, VM,lim is the limiting value of the Von-Mises stress of the material.

16. The method for analyzing the fatigue life of cylindrical roller bearings based on multi-factor coupling according to claim 15, characterized in that: The calculation formula in step S303 is as follows: The above formula is used to superimpose the IH life of each rolling element obtained in step S302 at the contact point with the rotating ring and the stationary ring to obtain the fatigue life of the entire bearing.

Images