Bearing Life Calculations in Rotating and Oscillating Applications

In 1999, Houpert [1] provided a model for calculating the roller bearing life of an oscillatory inner ring (IR) as a function of the oscillation angle θ and load zone parameter ɛ. In 2020, Breslau and Schlecht [2] provided a numerical model for calculating the life of an oscillating roller bearing accounting for both the IR and outer ring (OR) and identified some inconsistent results observed in two figures in Ref. [1]. By coincidence, the same inconsistencies were also noted by Menck in the same year in the scope of a project on the same topic on which Houpert and Menck were collaborating. The 1999 published oscillation factor of the oscillating inner ring was found to increase as the oscillation angle θ decreases and to “usually” decrease as the load zone parameter ɛ decreases, except in some case (large θ, low ɛ) where the factor was increasing as ɛ decreases, the latter being inconsistent. It will be demonstrated that the survival probability S of an oscillating ring is defined using $ln(1/S)∼(Qeqp⋅N)e$ and hence the product $Qeqp⋅N$⁠, where p is the life—load exponent, Q_eq is the equivalent load, and N is the number of stress cycles. Both parameters (Q_eq and N) can be defined using either a constant oscillating arc (2θ) or a variable loaded arc H(ψ,θ), see the upper images in Fig. 2. Since these arcs are used at the denominator when defining $Qeqp$ and numerator when defining N, so that the same results are obtained when using either 2θ or H(ψ,θ). The mindless error conducted in 1999 was to use 2θ for defining $Qeqp$ and H(ψ,θ) for defining N. This single error in Ref. [1] has been rectified, and as such the first objective of this article is to explain the 1999 error and 2020 correction as well as to present two updated figures. In addition, the original corrected model is expanded to include the outer stationary ring’s lifetime to obtain the lifetime of the entire bearing. Results are given for both ball and roller bearings of arbitrary dimensions by means of an oscillation factor a_osc.

Before describing bearing life in continuous rotation and oscillatory conditions, the authors found useful to share, in a specific chapter, some general comments about life models and concepts and then continue with standard exponents and life model. Then, analytical relationships for calculating the bearing life in continuous rotations are developed as these are required for explaining the correction factor to apply in oscillatory applications. This will be done with appropriate line contact (LC) sets of exponents (LC Dominik [3] or LC ISO) for roller bearings and using the appropriate set of point contact (PC) ISO exponents for ball bearings. The aforementioned analytical relationships can also be used for explaining and deriving the dynamic capacity (or rating) and the dynamic equivalent load using standard set of exponents (c, h, and e) or any new set.

The oscillation factor will then be described, accounting for both rings, as a function of the oscillation angle θ, the load zone parameter ɛ, geometrical parameter γ, and the osculation factors f_i and f_o used in ball bearings.

Finally, some useful curve-fitted relationships will be given for easy calculation of the final oscillation factor as a function of these input parameters. Moreover, approaches to calculate the required parameter ɛ from the results of the finite element analysis (FEA) are given.

Therefore, the following sections about life in continuous rotation are not new but have the merit of explaining in an analytical manner how to derive miscellaneous bearing life relationships. Readers who are solely interested in the oscillation factor may proceed directly to Sec. 3.3.

Cste is a constant, and c, h, and e are exponents to be defined experimentally. The exponent e is also called the Weibull slope.

The model has since been fine-tuned using elementary small volumes dV to integrate through the volume using a triple integration (accounting for possible pressure spikes at the edge of the roller-race contact and pressure increases on the skin of rough surfaces) and/or using an endurance limit when calculating the final life, leading to possible infinite life at low load when the maximum stress is below the endurance limit, see, for example, the work of Ioannides and Harris [6], Houpert and et al. [7], or Gnagy et al. [8]. Ioannides’ approach served as the basis in current ISO standards, see Refs. [9,10], for example, where the factor a_ISO has been outlined.

In all the previously cited references, the material was considered homogeneous with an implicit uniform distribution of inclusions in the volume, justifying the exponent 1 on the volume V.

In this article, the authors offer some general comments about life models and steel cleanliness effects, even if such comments are not the primary objectives of this article.

Ai developed an advanced life model in Ref. [11] accounting for steel cleanliness when simulating the number of inclusions (in a given stressed volume), inclusions acting as stress raisers. A steel cleanliness life reduction factor was then found to be inversely proportional to the cumulative length of inclusions raised to the power 0.865.

Furthermore, it can be demonstrated that for a given steel cleanliness, the number of inclusions, hence also the cumulative length, is proportional to the volume, meaning that one can question the exponent 1 on the volume V.

More importantly, it can also be demonstrated that, when simulating a given high quality steel cleanliness (or inclusion density), the probability of finding one inclusion in a batch of 24 small ball bearings (ball bearings of size 6204 tested in Ref. [6]) is almost nil, thus explaining the perceived infinite life and introduction of an endurance limit. However, when simulating bearings 100 times larger (as found in wind turbine applications, for example), the stressed volume is approximately 1 million times greater and many inclusions will be numerically simulated, causing a true risk of a finite life, so that the concept of endurance limit and infinite life can be questioned too. Bearing endurance tests are conducted on small- to medium-sized bearings for testing current ratings, but very large size bearings are not tested so that current dynamic ratings of large size bearing can be questioned. Note that a higher steel quality is used in large bearings for compensating the risk of inclusions effects on life or rating.

Even when dealing with small-size to medium-size bearing, ratings should be redefined using recently defined exponents (c, h, and e) and life-load exponent p. Large variations of the ratings and calculated bearing life (with the new set of exponents) can be expected when extrapolating these calculations to large size bearings.

Ai did not use any endurance limit and suggested a set of exponents (c = 11.385, h = 0.319, and e = 1.278), leading to a life-load exponent p equal to 4.72 when using Eq. (5) developed later (instead of 10/3 in Dominik’s model or 4 in ISO standards) for matching endurance test results. This exponent is in line with the proposal of Zaretsky et al. of using p = 5 in Ref. [12]. For ball bearings, Londhe suggested p = 4.1 (instead of 3) in Ref. [19].

In the scope of this study, roller and ball bearing rating C and dynamic equivalent load P_eq relationships will be derived using standard exponents, keeping analytically track of the (c, h, and e) exponents and without using an endurance limit (as initially done by Lundberg and Palmgren). Harris and Kotzalas conducted a similar exercise in Ref. [13], splitting this exercise into steps to define first a contact rating Q_c and then an equivalent load Q_e obtained using appropriate summations as also done in ISO 16281.

However, we can employ more concise analytical calculations here (thereby avoiding additional steps) to derive the roller and ball bearing rating and dynamic equivalent load, drawing inspiration from Dominik’s work [14] (for roller bearings only) and Houpert’s approach (also for roller bearings) described in Ref. [3].

Since the primary focus of this article is to define an oscillation factor, this article describes the main concepts used for calculating the bearing life of the inner ring and outer ring and, in particular, identifies the differences to be taken into consideration when studying oscillatory applications. The details presented in this section could be used in the future should the need arise to revisit rating and dynamic equivalent load relationships with any new set of exponents (c, h, e, and p).

The Hertzian exponent n is equal to 10/9 for roller bearings and 1.5 for ball bearings and 2 · ψ_l is the loaded arc (to be compared later with the loaded arc H in oscillatory application from Sec. 3.1.1). J_r and J_a can also be calculated via the integrals L_r and L_a used by Houpert in Ref. [15].

Hertz theory gives for roller-race line contact: τ₀ = 0.25 · P_max (where P_max is the maximum contact pressure) and z = 0.5 · b, where b is the half contact width. Slightly smaller ratios are used when calculating ball-race point contact. Note that τ₀ is described next as τ for allowing the possibility of introducing the index i and o, differentiating the inner from the outer ring, respectively. P_max (therefore also τ) and b (therefore also z) can be calculated as a function of the contact load using Hertz analytical relationships for LC or Hertz numerical results for PC (requiring complex numerical calculations of elliptical integrals), curve fitted by Houpert in Refs. [16,17].

In the following, outer ring results will be expressed as a function of the inner ring ones. For the sake of simplicity, it will first be assumed that the inner ring (index i) is rotating and the outer ring (index o) is stationary, although more general cases (for example, both rings rotating) could also be considered.

Here, E_eq = E/(1 – v²) = 2.26 × 10⁶ N/mm² is the equivalent Young's modulus derived from the Young's modulus E and Poisson’s ratio v of steel, L_i,o is the equivalent roller-race length, and R_{x_eq} is the equivalent contact radius in the rolling direction x. Variables D and d_m are the rolling element diameter and pitch diameter, respectively, and α is the nominal contact angle.

Equation (5) shows that the number of stress cycles N for a given survival probability S_dψ is described by N ∝ Q^−p. So far Q may have been considered constant, but the use of an elementary volume represented by dψ will make it possible to account for a variable load on the rotating ring (via the introduction of an equivalent load Q_eq, in Eq. (6)) or stationary ring (via the product of survival probabilities of all elementary volumes as will be shown in Eq. (11)).

Let us now consider a varying load Q that each elementary volume (defined by dψ) of the rotating ring will face during its excursion through the load zone spanning from –ψ_l to +ψ_l, the loaded arc being 2 · ψ_l. The percentage of occurrence %_j for this elementary volume of facing any load level Q_j is constant. Each load Q_j corresponds to a potential number of stress cycles N_j. The final or weighted number of cycles N that this elementary volume can endure is defined using Miner’s rule and results in the derivation of an equivalent load Q_{eq_rot} for the rotating ring.

The load Q_{eq_rot} can now be considered as constant for all elementary volumes (represented by dψ) of the rotating ring (as opposed to what will be seen in oscillatory applications) and can be used for calculating S_{dψ_i} of the rotating ring according to Eq. (5). $Qeq_rotp$ can also be considered as the mean value of Q^p over the loaded arc of an elementary volume’s movement with 2 · ψ_l appearing at the denominator.

Note how the loaded arc 2 · ψ_l seen by all elementary volumes (appearing at the numerator of N_i and denominator of $Qeq_rotp$⁠) disappears in Eq. (5) when using the product $(Qeq_rotp⋅Ni)e$⁠. This explains why Refs. [4,13] obtain the same result without consideration of the load zone size, replacing 2 · ψ_l with 2 · π, for example.

The final survival probability S_i of the entire rotating ring is obtained using the product of survival probabilities S_{dψ_i} or the integral of ln(1/ S_{dψ_i}) as defined in Eq. (5), explaining the factor 2 · π in front of the mean race radius R_i when defining the volume V_i, see Eq. (9).

B_LC can also be used for defining the life of the outer race as a function of the life of the inner race, see Eq. (31).

The second equation in Eq. (15) is important since it can be used to derive analytically the bearing life N_rev as a function of the survival probability, applied radial load, and load zone, keeping track of any exponents c, h, e, and p selected. It can also be used to derive the dynamic rating and dynamic equivalent load using the selected exponents c, h, e, and p.

When fixing the survival probability to a reference value (S_ref = 0.9), the number of revolutions obtained is often called L₁₀ (in revolutions) and the corresponding radial load F_r is the dynamic rating called C_ref in this article or C_r (or C₁ in this article) at ISO and C₉₀ by Dominik. C₁ and C₉₀ are hence defined using a life of reference (L₁₀ = L_ref = 1 million or 90 million revolutions, respectively), but also with a load zone of reference (2 · ψ_{l_ref} = 180 deg for ISO and 150 deg for Dominik) used for defining all integrals (J_{r_ref}, K_{bi_ref}, and K_{be_ref}) and the ratio $Kbe/Kbie$ used in B_{LC_ref}.

Cste is equal to 1.39 · 10⁻⁴⁹ for matching M = 57.6 (used for defining Dominik C₉₀ corresponding to 90 million revs) or M = 222.18 (used for defining Dominik C₁ corresponding to 1 million revs) with a load zone of reference equal to 150 deg.

The latter ratio can be plotted versus the load zone parameter ɛ or the ratio F_a/F_r/tan α or F_r · tan α/F_a, see Fig. 1 or 3 in Ref. [3], or can be used for defining the load zone factor using both rings (thus including B_LC) or only the rotating inner ring (thus fixing B_LC = 0), see Houpert’s Figs. 6 and 7 in Ref. [3].

Minor differences are observed when using one or two rings for calculating P_eq, so only P_{eq_IR+OR} is shown in Fig. 1.

A degree 6 polynomial relationship is also suggested in the same figure for easily calculating the curve-fitted ratio P_{eq_IR+OR_cf} as a function of the ratio F_r · tan α/F_a.

Note how P_eq increases (and therefore the life decreases) when F_a is either very large or small (the latter causing a narrow load zone) with both cases causing an increase of Q_max.

The ratio P_{eq_Catalog}/F_r is fixed to 1 in any “Catalog” approach as soon as the estimated ratio F_a/F_r/tan α is smaller than 1.5 and is equal to 0.4 · (1 + F_a/F_r/tan α) otherwise. Such a simplification can result in very erroneous overestimation of the bearing life at low axial force (or narrow load zone). Also, in a Catalog approach, the final axial force F_{a_1,2} in each row is not calculated precisely, but rather estimated using the axial force equilibrium accounting for both induced (by the radial force) axial forces (calculated using approximately ±1.26 · F_{r_1,2} · tan α if ɛ_ref = 0.5) and the external axial force, meaning the final axial force is not very accurate in the unseated row (having a narrow load zone).

As a reminder, the parameter B has been introduced for including the stationary outer ring. A similar exercise could be conducted for introducing another parameter (for example, C) accounting for the roller set effects on bearing life and rating using (1 + B + C) instead of (1 + B), but such an exercise is out of the scope of this study only describing standards.

In addition, similar relationships could be derived for describing the axial dynamic rating and dynamic equivalent load. Slightly different relationships can be derived when accounting for a rotating outer race with a zero or nonzero inner race rotation.

Finally, an inconsistency can be noted when using the ISO roller bearing life practices, where the final L₁₀ life is calculated using (C₁/P_eq)^10/3, hence p = 10/3, while the rating and equivalent load are defined with p = 4.

When analytically developing a relationship for ln(1/S_i) (in a similar manner to Eq. (15), but not explicitly shown in this article), the product $CAi⋅CBi1−h⋅CPic$ will be used when tracing the final exponents applied on (1 ± γ) and osculation ratios.

Only minor differences (due to the curve-fitted relationships used not only by Houpert but also by Palmgren) are observed.

The relationships described earlier are not only important for analytically deriving bearing life models in rotating applications but also for understanding the development of life models in oscillating applications.

The reader is also reminded that the life of the rotating inner ring is finally calculated using the product $(Qeq_rotp⋅Ni)e$ in which the loaded arc 2 · ψ_l cancels out. Moreover, for a given ɛ, the equivalent load $Qeq_rot$ as well as the loaded arc 2 · ψ_l of an elementary volume on the rotating race are constant for any elementary volume dV_i on the rotating inner ring.

In Sec. 4, when investigating oscillatory applications, it will become apparent that the equivalent load of the oscillating ring and the loaded arc H of an elementary volume on the oscillating race will be a function of ɛ, ψ, and the oscillation angle θ.

Based on Sec. 3, we will now derive a factor a_osc to correct the lifetime for oscillatory behavior.

An oscillatory inner ring oscillating between ±θ at an oscillating frequency freq_i will be studied first. The loaded zone is still defined by the angle ψ_l as shown in the upper left corner of Fig. 2. Each elementary volume dV_i is defined by its position ψ, nil in the 12 o’clock position and equal to 90 deg or π/2 rad in the 9 o’clock position in Fig. 2.

As explained in Secs. 4.1.2 and 4.1.3, the calculation of the loaded arc H(ψ) that an elementary volume will face during one oscillation is conceptually needed for correctly calculating the number of stress cycles that each elementary volume dV_i (located at the orbital angle ψ on the oscillating inner ring) will endure as well as its equivalent load Q_eq(ψ) but is not essential because H(ψ) will cancel out in the final derivation of the oscillation factor. However, it is useful to explain how H(ψ) can be derived for obtaining accurate calculation of the lower and upper bounds of the required integrals used in Eq. (35), Sec. 4.1.3.

It is shown in Fig. 2 that the loaded arc H(ψ) is a function of θ, but it is now also a function of the initial orbital angle ψ and load zone angle ψ_l.

An analytical relationship for H(ψ), only applicable in the case θ > ψ_l, is given in Table 1 of Ref. [1] and can now be extended to all cases (see updated Table 1), including the case not considered analytically in 1999 where two partial incursions of dV_i in the load zone (by its two ends) occur as shown in the lower left corner of Fig. 2 (with an example corresponding to ψ_l = 150, θ = 120, and ψ = 160 deg):

The lower right corner of Fig. 2 shows how H(ψ) varies as a function of ψ (with ψ_l = 150 and θ = 120 deg), with the possibility of observing a second plateau corresponding to H = 180 deg in Fig. 2 when a second incursion in the load zone occurs (by the other end of the load zone).

With the lower and upper bounds accurately defined, the required integrals can be calculated using only a few points (101 points for the results in this article, using an Excel table) and an accurate integration method.

Therefore, the equivalent load Q_eq and loaded arc H are not constant and vary as a function of ψ, as shown in Fig. 3 using a given example.

The equivalent load distribution in Fig. 3 can be compared to Breslau’s definition in Ref. [2] defined with 2 · θ (instead of H) at the denominator.

The survival probability S_{dψ_i} of the elementary small volume is obtained as in Eq. (5) using the product $[Qeqp(ψ)⋅Ni(ψ)]e$⁠, in which the loaded arc H cancels out since it is used at the denominator when defining Q_eq^p and numerator when defining N_i. The final results are therefore identical irrespective of whether H is defined correctly (conceptually at least with the number of stress cycles being proportional to a variable loaded arc H(ψ)), defined as constant and proportional to 2 · θ (as is the case in Refs. [2,18]) or defined as constant and equal to 1, for example (which makes no sense).

The only error in Ref. [1] is that 2 · θ was used for defining Q_eq^p, while H was correctly used for defining N_i, meaning that cancellation of the denominator and numerator was not possible.

The following steps are described in Ref. [1], where the final survival probability S_i of the oscillating inner ring is obtained by integrating ln(1/S_{dψ_i}) over the entire inner ring volume V_i, leading to the use of a double integral.

These relationships are identical to Eq. (32) in Ref. [1] (where the factor 2 · π was missing—a typographical error—but used in the calculations).

When fixing the acceptable accuracy to 10%, the aforementioned approximation can be used for θ < 36 deg when ɛ > 0.05. The accuracy improves to 2.34% when limiting θ to 10 deg and improves for any θ values as ɛ increases.

On the stationary outer ring, the number of stress cycles that any volume dV_o endures over a time period t_o can be defined using the cage frequency freq_c and is simply equal to 2 · Z · (2 · θ)/(2 · π) · freq_c · t_o, whereas it was Z · freq_c · t_o in continuous rotation.

We will now apply the aforementioned calculations to a variety of cases. Starting with the inner ring factor a_{osc_IR}, we will then turn to the factor a_{osc_IR+OR} for the entire bearing and consider the effect of γ and osculation.

It is now possible to correct Figs. 5 and 6 published in Ref. [1] and replace them with the following Figs. 4 and 5. In Fig. 4, the oscillation factor of the rotating IR can be seen as increasing when the oscillation angle decreases. Furthermore, by using any given oscillation angle, the oscillation factor decreases as ɛ decreases or increases as ɛ increases for reaching asymptotically a plateau when ɛ > 4.

Figure 5 shows the same results as Fig. 4, but divided by a_Harris. The results for θ = 180 are not visible in Fig. 5, because they are constantly equal to 1, therefore being identical to the results for θ = 360.

For the sake of completeness, Figs. 6 and 7 are presented using the ISO line contact (for roller bearings) and point contact (for ball bearings) sets of exponents, respectively, observing a smaller effect of ɛ on the final factor, the lowest ratio being of the order of 0.8 for ISO instead of 0.5 for Dominik.

Appropriate curve-fitted relationships for calculating a_{osc_IR} are presented in Sec. 4.5.1.

Equation (42) can be used for LC or PC using the appropriate relationships and values for B.

Equation (42) is already used in Ref. [13] to produce Figs. 11.28 and 11.29, unfortunately with an incorrect value of a_{osc_IR} in some cases.

Figure 11.28 should, for example, be replaced Fig. 8 (corresponding to L_o = L_i and γ = 0.1), γ = 0.1 being a typical geometrical ratio (γ = D·cos_a/d_m) used when deriving Fig. 11.28. The final oscillation factor continuously decreases as ɛ decreases, while this was not always the case in Fig. 11.28, especially when θ was large. Note that Fig. 8 is also shown in the appendix using the same format as in Ref. [13] (A_osc, here called a_{osc_IR+OR_LC}, instead of a ratio relative to Harris), to compare with Figs. 11.28 and 11.29.

Again, for the sake of completeness, the ISO LC and PC results are shown in Figs. 9 and 10 using γ = 0.1 (and f_i = f_o when defining B_PC).

Appropriate curve-fitted relationships are presented in Sec. 4.5.2.

The previously described models can also be used for scanning on γ and/or the osculation ratio f_o, for example, see Fig. 11 (using θ = 1 deg and f_i = f_o).