A method for simulating early spalling evolution and life prediction of rolling bearings
By combining Hertz contact theory and elastohydrodynamic lubrication theory with finite element model and damage constitutive equation, the simulation problem of early damage initiation and spalling evolution of rolling bearings was solved, realizing the simulation of early damage mechanism and life prediction of rolling bearings, and improving prediction accuracy and efficiency.
Patent Information
- Application Number
- CN202311240930.X
- Authority / Receiving Office
- CN · China
- Patent Type
- Patents(China)
- Current Assignee / Owner
- Filing Date
- 2023-09-25
- Publication Date
- 2025-10-28
- Estimated Expiration
- 2043-09-25
AI Technical Summary
Existing methods for simulating the damage evolution of rolling bearings lack physical and mechanical mechanisms. Traditional models are complex and lack generalizability, making it difficult to effectively simulate the early damage initiation and spalling process, resulting in inaccurate life predictions.
The Hertz contact theory and elastohydrodynamic lubrication theory were used to calculate the contact stress distribution, establish a two-dimensional planar finite element model of the bearing raceway, generate a Voronoi model, and combine it with Abaqus software to simulate damage evolution. The damage constitutive equation was written in Fortran language to calculate the damage degree and predict the cycle life.
It enables rapid and effective simulation of the early crack and spalling formation process of rolling bearings, improving the accuracy and efficiency of life prediction, and is suitable for accelerated fatigue test verification of rolling bearings throughout their entire life.
Smart Images

Figure CN117291075B_ABST
Abstract
Description
Technical Field
[0001] This invention belongs to the field of rotating machinery fault diagnosis, specifically relating to a method for simulating the early spalling evolution and life prediction of rolling bearings. Background Technology
[0002] Main bearings are core supporting components of critical equipment such as aero-engines and ground-based gas turbines, and their failure is one of the leading causes of engine shutdowns and premature engine replacements. The main failure inducing factors of bearings can be divided into two categories: first, pitting corrosion caused by surface stress concentration (e.g., surface contamination, surface roughness, insufficient lubrication); and second, spalling from the initiation of subsurface cracks to the surface. Under conditions of correct alignment, installation, good lubrication, and no excessive overload, the final failure mode is subsurface spalling caused by rolling contact fatigue. However, existing methods, based on probabilistic statistical models, require a large amount of real experimental data for verification and lack a physical and mechanical mechanism for the evolution of rolling bearing damage, thus limiting the generalizability of simulations of the rolling bearing damage evolution process. Regarding mechanical and physical models, traditional crack propagation models based on fracture mechanics cannot simulate the initiation of initial bearing damage and its propagation path beneath the subsurface. Damage accumulation evolution failure models based on continuous damage mechanics are theoretically more advantageous; however, existing models are complex to establish and solve, lacking a unified framework. Therefore, establishing an efficient model for bearing damage initiation and evolution is essential. Summary of the Invention
[0003] Purpose of the invention: To address the problem of early damage initiation and spalling evolution in rolling bearings, this invention proposes a method for simulating the early spalling evolution and life prediction of rolling bearings, providing an important technical approach for simulating the formation of early damage and predicting the life of rolling bearings.
[0004] This invention includes the following steps:
[0005] Step 1: Based on the load conditions of the rolling bearing, calculate the contact stress distribution between the inner and outer raceways and the rolling elements using the bearing size parameters, Hertz contact and elastohydrodynamic lubrication theory, and obtain the maximum contact stress and half-width length.
[0006] Step 2: Based on the maximum contact stress and half-width length obtained in Step 1, calculate the elastohydrodynamic lubrication contact pressure distribution, establish a two-dimensional planar finite element model of the bearing raceway, and generate a Voronoi model to simulate irregularly shaped metal grains. Discretize the pressure distribution and apply it to the model surface to simulate the rolling process of the rolling element.
[0007] Step 3: Calculate the degree of damage to the bearing after it has been subjected to cyclic loads;
[0008] Step 4: Calculate the damage degree of the material elements of the two-dimensional plane finite element model of the bearing raceway under multiple cyclic loads, and store and update the data until crack initiation and evolution to surface spalling stop, thus obtaining the cycle life.
[0009] Step 1 includes: using Hertz contact theory, and based on the load on the bearing, obtaining the stress distribution P in the two-dimensional contact region:
[0010]
[0011] Where b is the half-width of the contact area, P max denoted as the maximum contact stress. x represents the x-coordinate of the distribution area.
[0012] Step 1 also includes: calculating the maximum Hertz contact stress P under line contact. max :
[0013]
[0014] Where l is the line contact length and w is the maximum contact load.
[0015] In step 1, the half-width b of the contact area is calculated using the following formula:
[0016]
[0017] Parameter b * The curvature and function Σρ can be obtained by looking up tables (Reference: T.A. Harris, M.M. Kotzalas, et al. Rolling Bearing Analysis: Basic Concepts of Bearing Technology [M]. Machinery Industry Press, 2010) and calculated based on relevant bearing parameters.
[0018] In step 1, the contact stress distribution considering lubrication is calculated using elastohydrodynamic lubrication theory:
[0019] Define the dimensionless load parameter W, velocity parameter U, and material parameter G as follows:
[0020]
[0021] Where Q represents the load per unit area, and α is the Barus viscous compressive strength, typically ranging from 1 × 10⁻⁶. -8 ~3×10 -8 Pa -1 .
[0022] In step 1, assuming a constant temperature and a Newtonian fluid environment, the Reynolds equation for one-dimensional linear contact dimensionless elastohydrodynamic lubrication is expressed as:
[0023]
[0024] Among them, P, H, Let X be the dimensionless pressure, oil film thickness, density, viscosity, and coordinate, respectively: P = p / P max H = hR / b 2 ; X = x / b; d is the differential symbol;
[0025] λ is the dimensionless velocity coefficient, expressed as λ = 3π 2 U / 4W 2 ;
[0026] b is the Hertz contact radius; P max ρ is the maximum Hertz contact pressure; ρ is the density of the lubricating oil (kg / m³). 3 ), ρ0 is the initial density of the lubricating oil at room temperature and atmospheric pressure, p is the oil film pressure (MPa), η is the fluid viscosity (Pa·s), u s Where is the fluid inlet velocity (m / s), h is the oil film thickness (m), and R represents the equivalent radius of curvature of the contacting object.
[0027] In step 1, the equation for calculating the dimensionless film thickness is:
[0028]
[0029] X in and X out Let dX' represent the points where the lubricating oil enters and leaves the coordinate system, respectively, and dX' represent the derivative with respect to the horizontal axis.
[0030] In step 1, the initial film thickness H0 is calculated using Dawson's formula:
[0031] H0 = 2.65G 0.54 U 0.7 W -0.13 (7)
[0032] The viscosity equation, density equation, and load balance equation are as follows:
[0033]
[0034]
[0035]
[0036] Where p0 is a constant related to pressure viscosity, taken as 1.98 × 10⁸ Pa; η0 is the initial viscosity of the lubricating oil; z represents the viscosity-pressure exponent, with a value of 0.6; X in and X outLet dX represent the coordinates of the entry and exit points of the dimensionless contact region, respectively, and let dX represent the differential. Solve the simultaneous equations (5) to (10) according to the Newton-Raphson method or multigrid method to obtain the final elastohydrodynamic lubrication pressure distribution P.
[0037] Step 2: Based on the maximum contact stress and half-width length obtained in Step 1, calculate the elastohydrodynamic lubrication contact pressure distribution, establish a two-dimensional planar finite element model of the bearing raceway based on Abaqus software, and generate a Voronoi model to simulate irregularly shaped metal grains. Discretize the pressure distribution and apply it to the model surface to simulate the rolling process of the rolling element.
[0038] Step 3: Use Fortran language to write the damage constitutive equation of the bearing material element in the secondary development function of Abaqus software, and calculate the degree of damage to the bearing after being subjected to cyclic load.
[0039] Step 4: Calculate the damage degree of the material elements of the two-dimensional plane finite element model of the bearing raceway under multiple cyclic loads, and store and update the data until crack initiation and evolution to surface spalling stop, thus obtaining the cycle life.
[0040] Step 2 includes: Within the contact area between the rolling element and the raceway, the problem is simplified to a plane strain problem along the axial direction, i.e., the rolling direction, with respect to the Hertz contact. The rolling direction is along the contact half-width b. Within the contact area, the stress P of the contact half-width b is distributed on the contact surface. The rolling process of the bearing rolling element from the raceway is simulated using a discretized form from left to right. A two-dimensional finite element model with a length of 10b and a depth of 6b is established to simulate crack initiation and damage evolution. A larger quadrilateral mesh (typically 20μm) is used below 1.5b on the surface; a smaller quadrilateral mesh (typically 10μm) is used above 1.5b on the subsurface. Different mesh shapes are beneficial for simulating the random propagation effect of subsurface cracks and spalling.
[0041] In step 3, the damage accumulation equation for the material element is:
[0042]
[0043] Where Δτ represents the range of shear stress variation experienced by the material during contact; D is the degree of damage to the material element; N is the number of cycles; and the material parameters m and σ... rThese are parameters related to bearing steel materials; according to the test (reference: Shigeo, Shimizu, Kazuo, et al. Probabilistic Stress-Life (PSN) Study on Bearing Steel Using Alternating Torsion Life Test[J]. Tribology Transactions, 2009, 52(6): 807-816), m = 11.1, σ r =5979 MPa. d is the differential symbol;
[0044] Theoretically, the damage to a rolling bearing accumulates with each rolling element's cyclic rolling. However, the number of cycles for a rolling bearing typically exceeds one million, and calculating each cycle individually requires significant computational power, making it difficult to implement quickly. Therefore, in the simulation, the damage increment is assumed to consist of two or more cyclic blocks on average. Within each cyclic block, the maximum damage increment ΔD of a certain element in the material remains constant, and the shear stress variation range Δτ is fixed within a given cyclic block. Then, the number of cycles ΔN in the k-th cycle segment is determined. k Represented as:
[0045]
[0046] Wherein, the subscript k,crit represents the maximum value of dD / dN in the material element of the k-th cycle block;
[0047] The damage degree of each material element in the loop block is then updated as follows:
[0048]
[0049] Where j represents the number of each unit in the material; This represents the damage increment of the j-th unit;
[0050] At the beginning of the next loop block, i.e., the (k+1)th loop, we have:
[0051]
[0052] in, This represents the damage level of the j-th unit in the k-th cycle. This represents the damage degree of the j-th element in the material at the start of cycle k+1. Let be the elastic modulus of the j-th element in the k-th cycle. This represents the elastic modulus of the j-th element in the material at the start of cycle k+1.
[0053] Step 3 only performs one "cyclic block" calculation. However, the damage evolution of rolling bearings requires multiple "cyclic blocks" to simulate the crack initiation and propagation process. Therefore, in step 4, the above steps need to be cyclically loaded to achieve the effect of multiple rolling simulations. This continues until crack initiation leads to surface spalling, at which point the calculation stops. The damage degree of the material element under multiple cyclic loads is calculated by combining Abaqus script commands to perform multiple cyclic loading simulations in step 3 to achieve the effect of multiple rolling simulations, and the data is stored and updated. The calculation stops when the damage degree D of each material element reaches 0.99, until crack initiation evolves to surface spalling, thus obtaining its cycle life.
[0054] In at least one embodiment of the present invention, the invention is applied to accelerated fatigue tests of rolling bearings throughout their entire lifespan to verify its effectiveness and accuracy. The invention possesses at least the following beneficial technical effects:
[0055] This invention provides a method for simulating the early spalling evolution process of rolling bearings and predicting their cycle life. It can quickly, effectively and conveniently simulate the early crack and spalling formation process of rolling bearings, providing an important technical approach for understanding the damage evolution mechanism of rolling bearings and predicting their cycle life. Attached Figure Description
[0056] The present invention will be further described in detail below with reference to the accompanying drawings and specific embodiments, and the advantages of the present invention in the above and / or other aspects will become clearer.
[0057] Figure 1 This is a flowchart illustrating the early damage evolution process of rolling bearings according to the present invention.
[0058] Figure 2 This is a schematic diagram of a 6206 type deep groove ball bearing.
[0059] Figure 3 This is a schematic diagram of the ABLT-1A bearing fatigue testing machine and its loading.
[0060] Figure 4 This is a schematic diagram of the experimental results of inner ring spalling in three sets of rolling bearings.
[0061] Figure 5 This is a schematic diagram of the elastohydrodynamic lubrication pressure distribution obtained under experimental conditions.
[0062] Figure 6a This is a schematic diagram of the simulation and prediction results of the present invention (crack initiation, N = 2.141 × 10⁻⁶). 8 ).
[0063] Figure 6b This is a schematic diagram of the simulation and prediction results of the present invention (crack propagation, N = 2.167 × 10⁻⁶).8 ).
[0064] Figure 6c This is a schematic diagram of the simulation and prediction results of the present invention (crack propagation, N = 2.169 × 10⁻⁶). 8 ).
[0065] Figure 6d This is a schematic diagram of the simulation and prediction results of the present invention (peeling formation, N = 2.172 × 10⁻⁶). 8 ). Detailed Implementation
[0066] like Figure 1 As shown, this invention provides a method for simulating the early spalling evolution and life prediction of rolling bearings, comprising the following steps:
[0067] Step 1: Based on the load conditions of the rolling bearing, calculate the contact stress distribution between the inner and outer raceways and the rolling elements using the bearing size parameters, Hertz contact and elastohydrodynamic lubrication theory, and obtain the maximum contact stress and half-width length.
[0068] Step 2: Based on the maximum contact stress and half-width length obtained in Step 1, calculate the elastohydrodynamic lubrication contact pressure distribution, establish a two-dimensional planar finite element model of the bearing raceway, and generate a Voronoi model to simulate irregularly shaped metal grains. Discretize the pressure distribution and apply it to the model surface to simulate the rolling process of the rolling element.
[0069] Step 3: Calculate the degree of damage to the bearing after it has been subjected to cyclic loads;
[0070] Step 4: Calculate the damage degree of the material elements of the two-dimensional plane finite element model of the bearing raceway under multiple cyclic loads, and store and update the data until crack initiation and evolution to surface spalling stop, thus obtaining the cycle life.
[0071] Step 1 includes: using Hertz contact theory, and based on the load on the bearing, obtaining the stress distribution P in the two-dimensional contact region:
[0072]
[0073] Where b is the half-width of the contact area, P max denoted as the maximum contact stress. x represents the x-coordinate of the distribution area.
[0074] Step 1 also includes: calculating the maximum Hertz contact stress P under line contact. max :
[0075]
[0076] Where l is the line contact length and w is the maximum contact load.
[0077] In step 1, the half-width b of the contact area is calculated using the following formula:
[0078]
[0079] Parameter b * The curvature and function Σρ can be obtained by looking up tables (Reference: T.A. Harris, M.M. Kotzalas, et al. Rolling Bearing Analysis: Basic Concepts of Bearing Technology [M]. Machinery Industry Press, 2010) and calculated based on relevant bearing parameters.
[0080] In step 1, the contact stress distribution considering lubrication is calculated using elastohydrodynamic lubrication theory:
[0081] Define the dimensionless load parameter W, velocity parameter U, and material parameter G as follows:
[0082]
[0083] Where Q represents the load per unit area, and α is the Barus viscous compressive strength, typically ranging from 1 × 10⁻⁶. -8 ~3×10 -8 Pa -1 .
[0084] In step 1, assuming a constant temperature and a Newtonian fluid environment, the Reynolds equation for one-dimensional linear contact dimensionless elastohydrodynamic lubrication is expressed as:
[0085]
[0086] Where P, H, ρ, η, and X are the dimensionless pressure, oil film thickness, density, viscosity, and coordinates, respectively, and are: P = p / P max H = hR / b 2 ; X = x / b; d is the differential symbol;
[0087] λ is the dimensionless velocity coefficient, expressed as λ = 3π 2 U / 4W 2 ;
[0088] b is the Hertz contact radius; P max ρ is the maximum Hertz contact pressure; ρ is the density of the lubricating oil (kg / m³). 3 ), ρ0 is the initial density of the lubricating oil at room temperature and atmospheric pressure, p is the oil film pressure (MPa), η is the fluid viscosity (Pa·s), u s Where is the fluid inlet velocity (m / s), h is the oil film thickness (m), and R represents the equivalent radius of curvature of the contacting object.
[0089] In step 1, the equation for calculating the dimensionless film thickness is:
[0090]
[0091] X in and X out Let dX' represent the points where the lubricating oil enters and leaves the coordinate system, respectively, and dX' represent the derivative with respect to the horizontal axis.
[0092] In step 1, the initial film thickness H0 is calculated using Dawson's formula:
[0093] H0 = 2.65G 0.54 U 0.7 W -0.13 (7)
[0094] The viscosity equation, density equation, and load balance equation are as follows:
[0095]
[0096]
[0097]
[0098] Where p0 is a constant related to pressure viscosity, taken as 1.98 × 10⁸ Pa; η0 is the initial viscosity of the lubricating oil; z represents the viscosity-pressure exponent, with a value of 0.6; X in and X out Let dX represent the coordinates of the entry and exit points of the dimensionless contact region, respectively, and let dX represent the differential. Solve the simultaneous equations (5) to (10) according to the Newton-Raphson method or multigrid method to obtain the final elastohydrodynamic lubrication pressure distribution P.
[0099] Step 2: Based on the maximum contact stress and half-width length obtained in Step 1, calculate the elastohydrodynamic contact pressure distribution. Establish a two-dimensional planar model of the bearing raceway using Abaqus software, and generate a Voronoi model simulating irregularly shaped metal grains. Discretize the pressure distribution and apply it to the model surface to simulate the rolling process of the rolling elements.
[0100] Within the contact region between the rolling elements and the raceway, the problem is simplified to a plane strain problem along the rolling direction with respect to the Hertz contact. In a very small contact region, a stress P with a contact half-width of b is distributed on the contact surface. The rolling process of the bearing rolling elements from the raceway is simulated using a discretized form from left to right. In Abaqus, a two-dimensional finite element model with a length of 10b and a depth of 6b is established to simulate crack initiation and damage evolution. To ensure computational efficiency, a 20μm quadrilateral mesh is used below 1.5b on the surface. A 10μm quadrilateral mesh is used above 1.5b on the subsurface. The different mesh shapes are beneficial for simulating the random propagation effects of subsurface cracks and spalling.
[0101] Step 3: Using Fortran in the secondary development function of Abaqus software, the damage constitutive equation for the bearing material element is written to calculate the damage degree of the bearing after being subjected to cyclic loading. The damage accumulation equation for the material element is:
[0102]
[0103] In the formula, Δτ represents the range of shear stress variation experienced by the material during the contact process. D represents the damage degree of the material element, and N represents the number of cycles. Material parameters m and σ r These are parameters related to bearing steel materials. Based on testing, m = 11.1, σ... r = 5979 MPa. Assume that within each cycle block, the maximum damage increment ΔD of a certain element in the material remains constant, and within a defined "cycle block," the range of shear stress variation Δτ is constant. Then, the number of cycles ΔN in the k-th cycle segment... k Represented as:
[0104]
[0105] In the formula, the subscript k ,crit This represents the maximum dD / dN value in the material element of the k-th cycle block. The damage degree of each material element in this "cycle block" can then be updated as follows:
[0106]
[0107] In the formula, j represents the number of each unit in the material. ΔD k j This represents the damage increment of the j-th unit. At the start of the next loop block, i.e., the (k+1)-th loop, we have:
[0108]
[0109] In the formula, This represents the damage level of the j-th unit in the k-th cycle. This represents the damage degree of the j-th element in the material at the start of cycle k+1, where N is the number of vertical cycles. Let be the elastic modulus of the j-th element in the k-th cycle. This represents the elastic modulus of the j-th element in the material at the start of cycle k+1.
[0110] Step 4: Calculate the damage degree of the material elements under multiple cyclic loading. By combining Abaqus script commands, perform multiple cyclic loading simulations on the above step 3 to achieve the effect of multiple rolling simulations, and store and update. Stop the calculation and update when the damage degree D of each material element reaches 0.99, until the crack initiates and evolves to the surface and forms spalling, and then stop to obtain its cycle life.
[0111] Example:
[0112] This invention is verified using accelerated fatigue testing of rolling bearings throughout their entire lifespan. A 6206 type deep groove ball bearing is used as the test object. Figure 2 As shown in Table 1, its physical dimensional parameters are as follows. The ABLT-1A bearing fatigue testing machine was used. The loading method was hydraulic loading, with a maximum radial load capacity of 30kN and an axial load capacity of 10kN. Figure 3 As shown. In the experiment, four identical rolling bearings were placed in the test head each time. To accelerate bearing failure and to make the rolling bearing speed close to that of a real aircraft engine, the test speed was 12000 r / min. Kunlun L-HM 46 anti-wear hydraulic oil was used to fully lubricate the bearings. The oil density was 0.85 kg / L, and its kinematic viscosity was 45.61 mm. 2 / s (40℃). To ensure that the contact fatigue failure mechanism of the rolling bearing is not altered, the equivalent dynamic load on the bearing must be less than or equal to half of its rated dynamic load. The rated dynamic load of the bearing used in the test is approximately 19.5kN. In the test, a 10kN radial load is applied to the weights and transmitted to the test head through hydraulic pressure. Figure 3 It can be seen that each bearing is subjected to a radial load of 5kN. Therefore, its rolling contact fatigue failure mechanism will not change.
[0113] Table 1
[0114]
[0115] The operating state of the bearing was estimated by monitoring its real-time vibration signals. A shutdown threshold was set at three times the initial normal operating vibration value. The test was terminated when the effective value of the bearing vibration signal samples exceeded the shutdown threshold by 30 times. Through three different sets of tests, three rolling bearings with natural single-point spalling damage on the inner ring were obtained, and their spalling conditions are as follows: Figure 4As shown, spalling pits are clearly visible to the naked eye in the inner ring of the bearing, and no pitting, dents, or wear were found in other locations on the inner and outer ring raceways, indicating that the bearing was well lubricated during operation. The total number of cycles for the three rolling bearings was approximately 2.1 × 10⁻⁶. 8 2.64×10 8 and 5.57×10 8 .
[0116] It can be calculated that under a radial load of 5 kN, the maximum contact stress between the rolling elements and the inner and outer rings of the 6206 type rolling bearing is 2.89 GPa, and the theoretical contact half-width b (along the rolling direction) is 0.19 mm. At the corresponding rotational speed, its elastohydrodynamic lubrication pressure distribution is as follows: Figure 5 As shown. Applying EHL pressure to the surface of the rolling bearing within a contact radius b = 0.19 mm, the resulting bearing damage evolution is shown in [the table / reference needed]. Figure 6a , Figure 6b , Figure 6c , Figure 6d As can be seen, the secondary surface crack propagates parallel to the secondary surface before extending to the surface and forming spalling. Its final simulated cycle life is 2.172 × 10⁻⁶. 8 Compared with the test cycle life results of the three groups of inner ring spalling, the errors were -3.4%, +17.7%, and +61.0%, respectively. The simulated cycle life prediction results are in good agreement with the test results of bearings No. 1 and No. 2. The error with the cycle life test results of bearing No. 3 is larger, but considering the dispersion of bearing life and the lag in the effective value threshold monitoring during the test, and since the predicted results are within the same order of magnitude as the test results, the overall simulation results are relatively accurate.
[0117] In its specific implementation, this application provides a computer storage medium and a corresponding data processing unit. The computer storage medium is capable of storing a computer program, which, when executed by the data processing unit, can run the invention's content regarding a method for simulating early spalling evolution and life prediction of rolling bearings, as well as some or all of the steps in various embodiments. The storage medium can be a magnetic disk, optical disk, read-only memory (ROM), or random access memory (RAM), etc.
[0118] Those skilled in the art will clearly understand that the technical solutions in the embodiments of the present invention can be implemented using computer programs and their corresponding general-purpose hardware platforms. Based on this understanding, the technical solutions in the embodiments of the present invention, or the parts that contribute to the prior art, can be embodied in the form of computer programs, i.e., software products. These computer program software products can be stored in a storage medium and include several instructions to cause a device containing a data processing unit (which may be a personal computer, server, microcontroller, MUU, or network device, etc.) to execute the methods described in various embodiments or certain parts of the embodiments of the present invention.
[0119] This invention provides a method for simulating the early spalling evolution and life prediction of rolling bearings, taking into account elastohydrodynamic lubrication pressure. Many methods and approaches exist for implementing this technical solution; the above description is merely a preferred embodiment of the invention. It should be noted that those skilled in the art can make various improvements and modifications without departing from the principles of this invention, and these improvements and modifications should also be considered within the scope of protection of this invention. All components not explicitly stated in this embodiment can be implemented using existing technologies.
Claims
1. A method for simulating the early spalling evolution and life prediction of rolling bearings, characterized in that, Includes the following steps: Step 1: Based on the load conditions of the rolling bearing, calculate the contact stress distribution between the inner and outer raceways and the rolling elements to obtain the maximum contact stress and half-width length. Step 2: Based on the maximum contact stress and half-width length obtained in Step 1, calculate the elastohydrodynamic lubrication contact pressure distribution, establish a two-dimensional planar finite element model of the bearing raceway, and generate a Voronoi model to simulate irregularly shaped metal grains. Discretize the pressure distribution and apply it to the model surface to simulate the rolling process of the rolling element. Step 3: Calculate the degree of damage to the bearing after it has been subjected to cyclic loads; Step 4: Calculate the damage degree of the material elements of the two-dimensional plane finite element model of the bearing raceway under multiple cyclic loads, and store and update the data until the cracks initiate and evolve to the point of surface spalling, thus obtaining the cycle life. Step 1 includes: calculating the maximum Hertz contact stress P under line contact conditions. max : Where l is the line contact length; w is the maximum contact load; In step 1, the half-width b of the contact area is calculated using the following formula: Parameter b * The curvature and function Σρ were calculated by looking up tables and based on bearing parameters, respectively. In step 1, the contact stress distribution considering lubrication is calculated using elastohydrodynamic lubrication theory: Define the dimensionless load parameter W, velocity parameter U, and material parameter G as follows: Where Q represents the load per unit area, and α is the Barus viscosity coefficient; In step 1, assuming a constant temperature and a Newtonian fluid environment, the Reynolds equation for one-dimensional linear contact dimensionless elastohydrodynamic lubrication is expressed as: Among them, P, H, Let X be the dimensionless pressure, oil film thickness, density, viscosity, and coordinates, respectively: d is the differential symbol; λ is the dimensionless velocity coefficient, expressed as λ = 3π 2 U / 4W 2 ; b is the Hertz contact radius; P max The maximum Hertz contact pressure; ρ is the density of the lubricating oil, ρ0 is the initial density of the lubricating oil at atmospheric pressure and normal temperature, p is the oil film pressure, η is the fluid viscosity, and u s Where is the fluid inlet velocity, h is the oil film thickness, and R represents the equivalent radius of curvature of the contacting object.
2. The method according to claim 1, characterized in that, Step 1 includes: using Hertz contact theory, and based on the load on the bearing, obtaining the stress distribution P in the two-dimensional contact region: Where b is the half-width of the contact area, P max denoted as the maximum contact stress, and x as the abscissa of the distribution area.
3. The method according to claim 2, characterized in that, In step 1, the equation for calculating the dimensionless film thickness is: X in and X out Let dX' represent the points where the lubricating oil enters and leaves the coordinate system, respectively, and dX' represent the derivative with respect to the horizontal axis.
4. The method according to claim 3, characterized in that, In step 1, the initial film thickness H0 is calculated using Dawson's formula: H0=2.65G 0.54 AT 0.7 IN -0.13 (7) The viscosity equation, density equation, and load balance equation are as follows: Where p0 is a constant related to pressure viscosity; η0 is the initial viscosity of the lubricating oil; z represents the viscosity-pressure exponent; X in and X out Let dX be the coordinates of the entry and exit points of the dimensionless contact region, respectively, and let dX represent the differential. Solve equations (5) to (10) simultaneously to obtain the final elastohydrodynamic lubrication pressure distribution P.
5. The method according to claim 4, characterized in that, Step 2 includes: within the contact area between the rolling element and the raceway, the problem is simplified to a plane strain problem along the axial direction, i.e., the rolling direction, with respect to the Hertz contact. The rolling direction is along the half-width b of the contact area. Within the contact area, the stress distribution P of the half-width b of the contact area is on the contact surface. The process of the bearing rolling element rolling from the raceway is simulated by a discretized form from left to right. A two-dimensional finite element model with a length of 10b and a depth of 6b is established to simulate the crack initiation and damage evolution process. A larger quadrilateral mesh is used below 1.5b on the surface, and a smaller quadrilateral mesh is used above 1.5b on the subsurface.
6. The method according to claim 5, characterized in that, In step 3, the damage accumulation equation for the material element is: Where Δτ represents the range of shear stress variation experienced by the material during contact; D is the degree of damage to the material element; N is the number of cycles; and the material parameters m and σ... r These are parameters related to bearing steel materials; d is the differential symbol; In the simulation, the average damage increment is assumed to consist of two or more cyclic blocks. Within each cyclic block, the maximum damage increment ΔD of a certain element in the material remains constant, and the range of shear stress variation Δτ is constant within a defined cyclic block. Therefore, the number of cycles ΔN in the k-th cycle segment is... k Represented as: Wherein, the subscript k,crit represents the maximum value of dD / dN in the material element of the k-th cycle block; The damage degree of each material element in the loop block is then updated as follows: Where j represents the number of each unit in the material; This represents the damage increment of the j-th unit; At the beginning of the next loop block, i.e., the (k+1)th loop, we have: in, This represents the damage level of the j-th unit in the k-th cycle. This represents the damage degree of the j-th element in the material at the start of cycle k+1. Let be the elastic modulus of the j-th element in the k-th cycle. This represents the elastic modulus of the j-th element in the material at the start of cycle k+1.