Predicting Friction of Tapered Roller Bearings with Detailed Multi-Body Simulation Models

: In the presented work, a parametric multibody simulation model is presented that is capable of predicting the friction torque and kinematics of tapered roller bearings. For a highly accurate prediction of bearing friction, consideration of solid and lubricant friction is mandatory. For tapered roller bearings in particular, the friction in the contact between the rolling element and raceway is of importance. Friction forces in the contact between the rolling element end face and inner ring rib as well as roller cage pocket contacts are also considered in the model. A large number of tests were carried out to validate the model in terms of the simulated frictional torque. Inﬂuencing variables such as speed, axial load, radial load, and temperature were investigated. The simulation results show good agreement with the measured friction torque, which conﬁrms that the model is well suited to predict frictional torques and therefore the kinematics of tapered roller bearings.


Introduction
As part of the ongoing efforts to increase drive train efficiency, there is a general trend to minimize frictional losses in gearboxes.Besides gear teeth, one of the main aspects is the rolling contact of bearings [1,2].By reducing preload forces in adjusted bearing arrangements, friction losses are often reduced.A disadvantage of this approach is that when it comes to axial operating clearance, the bearing arrangement is subjected to onesided loads, or in the case of an O arrangement, the overall system experiences heating.The clearance has a significant impact on the operational behavior of one or both bearings.It is particularly noticeable in applications where the bearing is subjected to low loads and insufficient lubrication throughout its life cycle and where the rollers are relatively large and heavy.Under these conditions, there is an increased sliding component in the motion of the rolling elements.Combined with rapid speed and load changes, this increases the risk of bearing damage caused by slippage [3][4][5].
For the calculation of bearing losses, there are many approaches available.Over the past years, numerous studies have been carried out and precise contact models developed in order to make the friction losses in highly loaded contacts calculable [6][7][8][9][10][11][12][13][14][15][16].The first empirical equations for bearing friction calculation were formulated by Stribeck [17], Sjovall [18], Lundberg and Palmgren [19][20][21][22].Soon, more complex computer-based models had been developed.Jones [23] and Harris [24] were two of the drivers of the growth in bearing modeling theory.With the continuous increase in demand for higher efficiency and lower power loss, the optimization of rolling bearings has become more and more important.Analytical studies on bearings were often accompanied by experimental investigations to validate the analysis, as is also performed in this work.
The most accurate friction calculation can be carried out with dynamic simulation models of roller bearings.A brief overview of the development of dynamics models can be found in [25][26][27].Nowadays, the most accurate models are the dedicated calculation tools developed by the bearing manufacturers.These detailed dynamic simulation models are also capable of predicting and analyzing conditions in which damages caused by the dynamics of the bearing components, as mentioned above, can occur.Among the best-known tools are BEAST (BEAring Simulation Tool) of the SKF company (Goteborg, Sweden) [28][29][30][31], BRAIN (BeaRing Analysis In NSK) of the NSK company (Maidenhead, England) [32,33], Caba3D (Computer Aided Bearing Analyzer 3D) of the Schaeffler company (Herzogenaurach, Germany) [34], CAGEDYN of the Timken company (North Canton, OH, USA) [35][36][37] and IBDAS used by NTN (Osaka, Japan) [38].Since these models are part of the company's know-how, the information regarding them is scarce or lacking in detail.Moreover, these models are not accessible for general research.In addition to the above-mentioned programs, there are a large number of models for the dynamic simulation of rolling bearings capable of taking into account both dynamic and contact in a single program.Most of these models have been developed for specific problems [39][40][41][42][43][44][45].Other studies using dynamic simulation models use simplified friction calculations and focus on structural deformation or vibration analysis [45][46][47][48][49][50][51][52][53][54][55][56][57].Therefore, their use in general problems could lead to unreliable results.

Materials and Methods
As described in Section 1, the reduction of preload forces in adjusted bearing arrangements can be used to minimize frictional losses in gearboxes.In order to predict the frictional losses of Tapered Roller Bearings (TRBs) and prevent critical operational conditions, a Multi-Body Simulation (MBS) model is needed.Due to the fact that there are no publicly available models that are generally applicable, sufficiently accurate, and validated, a highly detailed MBS model for TRB has been developed in the present paper.The MBS model runs under the program name LaMBDA (Lager MehrkörperBerechnung und DynamikAnalyse).It is described in Section 2.1.
To compare the simulated results with experimentally measured ones and validate them under different operating conditions, friction torque measurements have been performed.The experimental setup used is described in Section 2.2.Simulation results and measurements are presented in Section 3 and discussed in Section 4.

Multi Body Simulation Model
The MBS model was developed based on an approach that has been established for model development at MEGT (Chair of Machine Elements, Gears and Tribology).It is parameterized and modular and uses self-developed routines for calculation.Routines for such calculations to be used in MBS models have been developed and improved over the past 20 years at MEGT [58][59][60][61].
Within the model, a combination of commercial software and self-programmed calculation routines is used.Depending on the level of detail, all bodies of a rolling bearing are modeled with their real geometries.In this work, the most detailed model structure available is described.In order to define the bodies, markers, forces, and boundary conditions, a graphical user interface is used.It was also developed in this work and can access a database of selected bearings and lubricants.This allows for user-friendly model generation.The force elements of the model use the self-written routines to calculate the individual force components in the contact (normal force, damping, and friction).The flow chart in Figure 1 shows a simplified diagram of the procedure followed within the MBS model to calculate contact forces between the single elements by means of roller and ring raceway contact.In the MBS model geometry, material data and lubricant properties are defined.Together with the state parameters, this data is retrieved from the calculation core.In this case, state parameters mean, for example, the distance  ⃗  between the coordinate systems of the individual bodies and their velocities  ⃗  .Based on these, the contact between the rolling element and raceway is determined (see Section 2.1.1).The result values of the contact routine are contact point  ⃗, contact normal vector  ⃗ and penetration  .With a load-deformation relationship, the effective contact normal force  ⃗ . is determined from the penetration.In the next step, the time-dependent contact state parameters are calculated.These are the relative velocities  ⃗ and sum velocities  ⃗ in the contact point.The relative velocity is further used to calculate the damping forces  ⃗ as described in Section 2.1.2.The normal forces, velocities, and lubricant data are used to determine the lubricant film height ℎ in contact, the specific lubricant film height  and the solid load-bearing ratio  .On their basis frictional forces  ⃗ and resulting torques  ⃗ are calculated.The calculation of frictional forces for each contact in the model is presented in detail in Section 2.1.3.In the last step, the three components of the contact force  ⃗ are summed up and given back to the MBS solver for iteration of the force equilibrium.

Contact Calculation
In order to correctly represent the contact between two bodies, the area of contact must be discretized.Depending on the type of contact, the discretization used in this model can be one-or two-dimensional.In the case of a roller ring raceway contact, where there is line contact, disc models are the state of the art and are also used here.Depending on the geometry pairing, the contact between the rolling element face and the ring rib deviates greatly from the idealized point or line contact.Therefore, a cell model is used to calculate this contact.

Slice Model
A conventional slice model is described in [62].In this work, the Alternative Slicing Technique (AST) is used.It was first presented by Teutsch in [63] and subsequently implemented by Kiekbusch for different types of bearings in [26].The AST model allows radial deformation of the slices toward each other.However, the discs cannot twist against each other.Thus, the deformation of each slice can be determined by its penetration into the contour of the opposite slice of the raceway.Conventional slice models neglect the influence that neighboring slices have on each other.As a result, they do not adequately represent the pressure distribution under bearing loads, especially with tilting.The disk model (AST) used here allows the calculation of excess pressure in the case of edge bearing or strong angular misalignments.In the MBS model geometry, material data and lubricant properties are defined.Together with the state parameters, this data is retrieved from the calculation core.In this case, state parameters mean, for example, the distance → s (t) between the coordinate systems of the individual bodies and their velocities → .
s (t).Based on these, the contact between the rolling element and raceway is determined (see Section 2.1.1).The result values of the contact routine are contact point

Contact Calculation
In order to correctly represent the contact between two bodies, the area of contact must be discretized.Depending on the type of contact, the discretization used in this model can be one-or two-dimensional.In the case of a roller ring raceway contact, where there is line contact, disc models are the state of the art and are also used here.Depending on the geometry pairing, the contact between the rolling element face and the ring rib deviates greatly from the idealized point or line contact.Therefore, a cell model is used to calculate this contact.

Slice Model
A conventional slice model is described in [62].In this work, the Alternative Slicing Technique (AST) is used.It was first presented by Teutsch in [63] and subsequently implemented by Kiekbusch for different types of bearings in [26].The AST model allows radial deformation of the slices toward each other.However, the discs cannot twist against each other.Thus, the deformation of each slice can be determined by its penetration into the contour of the opposite slice of the raceway.Conventional slice models neglect the influence that neighboring slices have on each other.As a result, they do not adequately represent the pressure distribution under bearing loads, especially with tilting.The disk model (AST) used here allows the calculation of excess pressure in the case of edge bearing or strong angular misalignments.

Cell Model
Contact geometries that cannot be reduced to a point or a line are calculated with cell models.One example is the contact between the rolling element's end face and the inner ring rib.Depending on the profile at the rolling element end face and the ring rib, the contact area can take a crescent shape.For contact calculation, a limited area of the ring and the rolling element end face is divided into squares.Based on the body geometry, the penetration between the contact partners is calculated.For this purpose, the contact problem as described in POLONSKY and KEER is solved [64].The surface deformation is described by an influence matrix.It can be determined according to the Boussinesq equation [65].Kiekbusch shows in [60] that the problem can be solved most effectively in dynamics simulation with a combination of FFT (Fast Fourier Transformation) and CG (Conjugate Gradient) solver methods.Accordingly, the method described there is used in this work.

Damping
In dynamics simulation, parametric damping models are often used since a variety of influencing factors make it difficult to determine the exact damping values of the contact points.Material damping occurs as a result of deformations and the non-linear elastic properties of contact partners.An additional damping effect is lubricant film damping, which occurs in the run-in zone of the Elasto-Hydrodynamic Lubrication (EHL) contact between the rolling element and raceway [66].In this study, a parametric model is also used to model the damping force between the rolling element and the other bearing componentsthe ring raceway, ring ribs, and cage pocket.It is determined by two parameters, the maximum damping coefficient d max and the penetration δ max , above which maximum damping is reached.Thus, it represents an easily adjustable modeling variant of the damping.
The direction of the damping force is opposite to the impact velocity, which can be understood as the relative velocity of the two touching bodies in the normal direction The damping value d is calculated as a function of the penetration depth δ.
It is described in the range 0 ≤ δ < δ max with a continuous cubic function.
The damping model is applied to each disc or cell of the contact, depending on the contact point in the bearing.

Friction
The mathematical and physical principles for describing friction at the contact points in rolling bearings are complex since different friction phenomena occur depending on the contact geometry and relative motion.Rolling and sliding friction are essential for the contact between the rolling element and raceway.Pure rolling occurs when the surface velocities of the contact partners are equal in terms of magnitude and direction.If there is a tangential relative movement between the two contact partners in the contact area, sliding occurs.Both forms of friction are present in the developed model.Losses resulting from non-linear elastic material behavior are here shown as hysteresis moments.Figure 2 gives an overview of the friction components considered.These are calculated for each slice in the discretized contact.In the following, the implemented approaches to describing the friction components are explained.

Lubricant Friction (Sliding)
In EHL contact, a distinction can be made between a rolling and a sliding component of the frictional force [67].The sliding component results from the shear of the lubricant.The relative movement (sliding) between the rolling element and the ring raceway shears the lubricant.The resulting shear stresses  work against the movement of the rolling element.The force acting on the roler can be written as the integral of the shear stresses over A, which is the contact area between the rolling element and ring raceway.
If a constant equivalent shear stress  is taken as a basis, the relationship can be formulated as follows [67].
For the contact area A, the Hertzian contact area  is used as a basis.This approach can also be formulated in discrete time if  and  are known at each time step.The two models described below are available in LaMBDA for estimating the equivalent shear stress.
Bair and Winer assume that the shear stresses that a lubricant can carry are limited to a certain value, a limiting shear stress  .This value is a characteristic of the lubricant.The viscous component of the shear gradient  can be calculated from the quotient of the lubricants limiting shear stress  .and  as a logarithmic function [68].
The dynamic viscosity  is described by the modulus equation according to Dicke [69].
The parameters  ,  ,  and  of the equation must be determined for each oil from viscosity measurements.The lubricant viscosity  at ambient pressure is calculated according to Vogel with the lubricant-dependent parameters K, B, and C [69].The parameters describe the temperature dependence of the viscosity.The contact pressure is assumed to be p.In the following, the implemented approaches to describing the friction components are explained.

Lubricant Friction (Sliding)
In EHL contact, a distinction can be made between a rolling and a sliding component of the frictional force [67].The sliding component results from the shear of the lubricant.The relative movement (sliding) between the rolling element and the ring raceway shears the lubricant.The resulting shear stresses τ EHL work against the movement of the rolling element.The force acting on the roler can be written as the integral of the shear stresses over A, which is the contact area between the rolling element and ring raceway.
If a constant equivalent shear stress τ EHL is taken as a basis, the relationship can be formulated as follows [67].
For the contact area A, the Hertzian contact area A Hertz is used as a basis.This approach can also be formulated in discrete time if τ EHL and A Hertz are known at each time step.The two models described below are available in LaMBDA for estimating the equivalent shear stress.
Bair and Winer assume that the shear stresses that a lubricant can carry are limited to a certain value, a limiting shear stress τ L .This value is a characteristic of the lubricant.The viscous component of the shear gradient .γ can be calculated from the quotient of the lubricants limiting shear stress τ L .and η as a logarithmic function [68]. .
The dynamic viscosity η is described by the modulus equation according to Dicke [69].
The parameters a 1 , a 2 , b 1 and b 2 of the equation must be determined for each oil from viscosity measurements.The lubricant viscosity η 0 at ambient pressure is calculated according to Vogel with the lubricant-dependent parameters K, B, and C [69].The parameters describe the temperature dependence of the viscosity.The contact pressure is assumed to be p.
The temperature ϑ is given in • C in this equation.For the calculation of the heat of compression in the contact area, Gold et al. [70] propose the following relation: The constants C 1 and C 2 are lubricant dependent.The relative velocity in contact u rel and the lubricant film height h 0 are used to calculate the shear rate .γ [71]. .
The relative velocity is calculated from the surface velocities of the two bodies, the roller and the ring.The equations presented by Moes are used to calculate the lubricant film height.The numerical implementation of these equations is presented in [14] and is also used here.In addition, lubricant viscosity is considerably influenced by thermal effects [44].To take these into account, correction factors φ ϑ are used.
The relationship according to Zhu and Cheng takes into account not only the relative motion in the contact but also the influence of the maximum Hertzian pressure [72].
E is the reduced Young's modulus of both contacting bodies and p 0 is relative pressure.Deviations from perfect rolling, slippage is included as follows.
For the thermal load parameter Γ ZC , Zhu and Cheng [73] take the temperature gradient of the viscosity as a basis in addition to the average conveying velocity u av and the thermal conductivity λ ϑ .
This can be calculated by differentiating the equation according to Vogel [74].
The parameters A V , B V and C V are pressure-and temperature-dependent quantities that must be determined for each lubricant from viscosity measurements [74].

Lubricant Friction (Rolling)
Rolling resistance results from compressing and overrolling the lubricant in EHL contact.It is calculated according to Biboulet and Houpert [75] E and R are the reduced Young's modulus and the reduced radii of the contacting bodies and L the effective contacting length.The dimensionless parameters of these approaches are defined as follows: Load parameter The load parameter W is calculated on the basis of the load Q imposed on one slice of the rolling element.

Solid Rolling Friction
Rolling friction losses due to contact point displacement in the rolling direction caused by elastic deformation of the contact partners are taken into account as a function of normal force, as in Scheuermann [76].
The rolling resistance exponent e r and the rolling friction coefficient c r are materialdependent and can be taken from the literature on the subject or experiment.

Material Hysteresis
Hysteresis refers to a contact that does not take place ideally elastically.The deformation energy is not completely returned.Part of the deformation energy dissipates into heat [61,77].For the rolling contact between the rolling element and raceway, this means that the contact pressure is asymmetrical.When rolling under load, the contact elements deflect on the run-in side.However, due to the processes described, the deflection does not occur to the same extent.The force Q resulting from the compression is offset relative to the axis of rotation.It therefore causes a moment M T,Hys , which is directed in the opposite direction to the rotary motion of the two bodies.The resulting moment is independent of the velocity.According to Johnson, it can be described as a function of a hysteresis loss factor a v and is only proportional to the contact load and half of the Hertzian contact width b, which corresponds to the disk width in the discretized contact [78].

Solid Sliding Friction
In order to represent the proportion of solid friction in the contact between rolling elements and the raceway, a section-wise defined function is used.It is divided into a cubic part and a constant part.The distinction is made on the basis of the relative velocities in the contact.Using the input parameters v s and v d the friction coefficient-relative velocity curve can be adapted to experimentally determined values.Above the limit of v d , a constant value µ d is assumed for the solid friction.
The parameters a cubic , h cubic , and δ cubic are chosen for the ranges such that the value µ i decreases to the constant value µ d after a degressive slope to v s .
The resulting frictional force is calculated by multiplying the contact normal force with the described friction coefficient.

Mixed Friction
Mixed friction is the transitional area between solid friction and friction due to EHL.Both types of friction are significantly influenced by the surface roughness of the contact bodies.Between the roughness peaks, small micro-EHL structures form, which can only be captured iteratively with great effort.In order to take these effects into account in the dynamics simulation, a simple model according to Zhou and Hoepprich is used [61,73].It is used to approximate the proportion of solid friction in the mixed friction regime φ.This portion can be determined as a first approximation from the asperity load ratio of the surfaces.
The proportion of the normal force Q s transmitted at solid contacts is set in relation to the total load Q here.The parameters B ZH and C ZH describing the roughness of the surfaces and are determined according to [73].The lubricant film thickness parameter Λ is described as the quotient of the film thickness h 0 and the combined standard derivation of surface roughness σ.
Friction in Roller Rib Contact Friction occurring between the rolling element end face and ring rib is called drilling friction due to the rotational movement the roller is performing.This drilling movement can also be understood as sliding when looking at the relative velocity vector of each discretized point in the contact.The closer the discretized point is to the axis of rotation of the roller, the smaller the magnitude of the velocity.The further away it is, the greater its magnitude becomes.In order to describe the friction that occurs during a sliding movement, suitable approaches exist.Here, a distinction is made between sliding friction in solid contact (rolling element and rib) and the resistance of the lubricant due to its shear.

Solid Sliding Friction in Roller Rib Contact
The coefficient of friction for sliding in solid contact between the rolling element end face and the ring ribs is approximated via a cubic function, as it is used in the roller raceway contact (Equation ( 21)).In order to define the continuous function necessary for the MBS simulation, the pole point at v sl = 0 is approximated with a cubic section.The parameters v s , v d µ s and µ d can be used to set the range until a constant friction value is reached.

Lubricant Friction in Roller Rib Contact
The friction component due to sliding in the lubricated contact results from the shear of the lubricant.As in the contact between the rolling element and the ring raceway, Equation ( 4) is therefore used.The shear stresses in the lubricant are represented by the Bair-Winer model.In this case, in contrast to the line contact, the following equation according to [79] is used to calculate the lubricant film height.
The equation was originally derived for an elliptical contact surface and takes into account the velocity, material, and load using the parameters U i,j , (Equation ( 33)) G i,j (Equation ( 34)) and W i,j (Equation ( 35)).It can be used in good approximation also for this kind of contact, even if the contact area deviates strongly from an ellipse in some cases.The elliptic parameter k e is the ratio of the two contact axes a and b.It can be estimated from the reduced radii in x direction R x and y direction R y .
The dimensionless parameters for film height calculation result from the following variables: Velocity parameter For calculation of the material parameter G, the pressure-viscosity coefficient α p is used.The correction of the lubricant film height for a thermal influence and the description of the dynamic viscosity are carried out according to Equations ( 6)- (15).

Mixed Friction
The resulting friction force is calculated by summing the two friction force components (solid sliding friction and EHL sliding friction, or lubricant shear).Following Steinert, they are weighted with a dimensionless key figure for the solid body friction component [80].

Friction in Roller Cage Contact
The contact between the rolling element and cage pocket contributes only slightly to the total frictional torque of the bearing.There are high relative speeds or predominantly sliding between the two surfaces.Coulomb's friction law is therefore used to calculate the frictional force.
The friction value µ used as a basis is approximated as in Equation (21).The polarity of the coefficient of friction at zero velocity is thus avoided.The cage is guided by the rolling elements.Therefore, no further frictional contacts are created.The cage cannot touch the bearing rings.

Friction Torque Measurement
The measurements of the frictional torque were carried out on the friction torque test rig of the MEGT (see Figure 3).The test rig was developed at the MEGT in the framework of Aul's work [58].It has been used in many research projects over the past years and has proven itself reliable [26,[58][59][60]81].The test rig allows the frictional torque of a single test bearing to be measured.The outer ring of the test bearing is mounted on a hydrostatic bearing.A beam load cell is used to hold the rotational degree of freedom of the hydrostatic bearing in position while measuring the force required.With the distance between the axis of rotation of the test bearing and the beam load cell, the total frictional torque can be calculated from the measured force.
The measurements of the frictional torque were carried out on the friction torque test rig of the MEGT (see Figure 3).The test rig was developed at the MEGT in the framework of Aul's work [58].It has been used in many research projects over the past years and has proven itself reliable [26,[58][59][60]81].The test rig allows the frictional torque of a single test bearing to be measured.The outer ring of the test bearing is mounted on a hydrostatic bearing.A beam load cell is used to hold the rotational degree of freedom of the hydrostatic bearing in position while measuring the force required.With the distance between the axis of rotation of the test bearing and the beam load cell, the total frictional torque can be calculated from the measured force.
Furthermore, the test bearing is accessible for additional measurements such as cage or rolling element speed.The bearing can be loaded in the radial and axial directions.The tilting module also allows tilting and/or skewing of the test bearing, which was not used in the investigations presented here.Three sets of measurements were performed.First, the TRB was loaded axially with 6 kN.The bearing was lubricated with an oil bath.The oil level was initially set to half the roller height.The lubricant used was an ISO VG 100 mineral oil without additives (reference oil FVA No. 3-for lubricant data, see Table A1 in Appendix A).The speed was kept constant until the temperature on the outer ring of the bearing reached the desired level.The speed ramp was then run through within a few minutes so that there was no change in temperature during the measurement.The measurement was taken at two temperatures, according to the setup shown in Table 1.Second, a constant combined load was then set.The test sequence is the same as for purely axial loading.Then, the speed was varied at two constant temperatures, according to the setup shown in Table 2.The second series of measurements was carried out under combined radial and axial loads.The radial load was applied in a downward direction, meaning the load zone is located in the lower part of the rolling bearing that is well supplied with lubricant.The other boundary conditions were adopted.Furthermore, the test bearing is accessible for additional measurements such as cage or rolling element speed.The bearing can be loaded in the radial and axial directions.The tilting module also allows tilting and/or skewing of the test bearing, which was not used in the investigations presented here.
Three sets of measurements were performed.First, the TRB was loaded axially with 6 kN.The bearing was lubricated with an oil bath.The oil level was initially set to half the roller height.The lubricant used was an ISO VG 100 mineral oil without additives (reference oil FVA No. 3-for lubricant data, see Table A1 in Appendix A).The speed was kept constant until the temperature on the outer ring of the bearing reached the desired level.The speed ramp was then run through within a few minutes so that there was no change in temperature during the measurement.The measurement was taken at two temperatures, according to the setup shown in Table 1.Second, a constant combined load was then set.The test sequence is the same as for purely axial loading.Then, the speed was varied at two constant temperatures, according to the setup shown in Table 2.The second series of measurements was carried out under combined radial and axial loads.The radial load was applied in a downward direction, meaning the load zone is located in the lower part of the rolling bearing that is well supplied with lubricant.The other boundary conditions were adopted.
For the last series of measurements, a constant preload of 6.5 kN of the bearing was applied.At constant speed, the radial load was increased from 1 to 15 kN.These measurements have been performed under a steady temperature.The setup parameters are listed in Table 3.
A TRB of type 32216 was used for all tests.The geometry of the bearing is shown in Table 4. Figure 4 serves as an illustration of the geometric values.The profile retraction of the rolling elements is described by the parameters a p , c p , d p , k p and r e .They are defined in Teutsch [44] and allow the different standardized profiles of rollers to be described mathematically.
Lubricants 2023, 11, x FOR PEER REVIEW 12 of 22 Teutsch [44] and allow the different standardized profiles of rollers to be described mathematically.

Results
For model validation, the comparison of the frictional torque as an integral parameter of the contact modeling of all internal contacts in the bearing with the measured frictional torque has proved useful.The fictional torque is highly dependent on the load distribution inside the bearing and its dynamics.For this purpose, the total frictional torque is being

Results
For model validation, the comparison of the frictional torque as an integral parameter of the contact modeling of all internal contacts in the bearing with the measured frictional torque has proved useful.The fictional torque is highly dependent on the load distribution inside the bearing and its dynamics.For this purpose, the total frictional torque is being used as the main outcome for this comparison.
As shown in Figure 5, the frictional torque obtained in the simulation (LaMBDA 50 • C) when submitted to an axial load of 6.5 kN and a shaft speed of 500 rpm is about 1.200 Nmm.The value is rising with an increasing shaft speed up to 3800 Nmm.It can be observed that the simulated friction curve follows the measured one (measurement at 50 • C) with a small offset.The simulated frictional torque values are slightly higher.The second measurement shown in Figure 5 was taken at 42 • C. Since the load is the same, the changed frictional torque can be attributed exclusively to the temperature dependence of the lubricant properties.This concerns especially temperature-dependent changes in viscosity.In the simulation, the lubricant properties are included in the calculation of the rolling resistance as well as the losses due to shear of the lubricant.The fact that the simulated frictional torque follows the measured curve for both temperatures very well shows that the modeling of the lubricant properties is quite suitable.

Results
For model validation, the comparison of the frictional torque as an integral parameter of the contact modeling of all internal contacts in the bearing with the measured frictional torque has proved useful.The fictional torque is highly dependent on the load distribution inside the bearing and its dynamics.For this purpose, the total frictional torque is being used as the main outcome for this comparison.
As shown in Figure 5, the frictional torque obtained in the simulation (LaMBDA 50 °C) when submitted to an axial load of 6.5 kN and a shaft speed of 500 rpm is about 1.200 Nmm.The value is rising with an increasing shaft speed up to 3800 Nmm.It can be observed that the simulated friction curve follows the measured one (measurement at 50 °C) with a small offset.The simulated frictional torque values are slightly higher.The second measurement shown in Figure 5 was taken at 42 °C.Since the load is the same, the changed frictional torque can be attributed exclusively to the temperature dependence of the lubricant properties.This concerns especially temperature-dependent changes in viscosity.In the simulation, the lubricant properties are included in the calculation of the rolling resistance as well as the losses due to shear of the lubricant.The fact that the simulated frictional torque follows the measured curve for both temperatures very well shows that the modeling of the lubricant properties is quite suitable.The same experiments were performed with combined axially and radially loaded bearings.This load situation produces significant contact forces in the roller-raceway contact and between roller and rib.Therefore, it is suitable to validate the contact and friction models in the roller-rib contact as well.Both curves, simulation and experiment, show an increase in frictional torque with higher shaft speed and a decrease with lower temperature (Figure 6).The same experiments were performed with combined axially and radially loaded bearings.This load situation produces significant contact forces in the roller-raceway contact and between roller and rib.Therefore, it is suitable to validate the contact and friction models in the roller-rib contact as well.Both curves, simulation and experiment, show an increase in frictional torque with higher shaft speed and a decrease with lower temperature (Figure 6).Contrary to expectation, the measured frictional torque is noticeably lower than it is under pure axial load.This results from the changed load distribution in the bearing.Whereas with purely axial loads, the load zone extends over the entire bearing Contrary to expectation, the measured frictional torque is noticeably lower than it is under pure axial load.This results from the changed load distribution in the bearing.Whereas with purely axial loads, the load zone extends over the entire bearing circumference, with additional radial loads, it is located only in the lower part of the bearing.Consequently, fewer rolling elements and thus a smaller contact length contribute to the friction.To investigate this influence further, the radial bearing load is varied at a constant preload.The results from experiment and simulation are presented in Figure 7. Contrary to expectation, the measured frictional torque is noticeably lower than it is under pure axial load.This results from the changed load distribution in the bearing.Whereas with purely axial loads, the load zone extends over the entire bearing circumference, with additional radial loads, it is located only in the lower part of the bearing.Consequently, fewer rolling elements and thus a smaller contact length contribute to the friction.To investigate this influence further, the radial bearing load is varied at a constant preload.The results from experiment and simulation are presented in Figure 7.As previously assumed, the measurements carried out (in light blue) show that the total frictional torque of the test bearing decreases with increasing radial load.Again, the comparison between experiment and model shows a very good agreement.The simulation model reproduces all the internal contacts of the bearing individually and in detail, which means that the lower frictional torque is displayed with a radial load component.With the high-resolution contacts and the detailed friction calculation described in Section 2.1, the model is very good at predicting the frictional torque.Furthermore, the model provides an explanation for the falling friction torque curve as it allows an insight into the inside of the bearing and shows the load distribution.
As a second proof of the assumption, the load distribution in the rolling bearing is used.The load distribution for the load cases-pure axial load and combined load with a high radial portion-is shown in Figure 8.As previously assumed, the measurements carried out (in light blue) show that the total frictional torque of the test bearing decreases with increasing radial load.Again, the comparison between experiment and model shows a very good agreement.The simulation model reproduces all the internal contacts of the bearing individually and in detail, which means that the lower frictional torque is displayed with a radial load component.With the high-resolution contacts and the detailed friction calculation described in Section 2.1, the model is very good at predicting the frictional torque.Furthermore, the model provides an explanation for the falling friction torque curve as it allows an insight into the inside of the bearing and shows the load distribution.
As a second proof of the assumption, the load distribution in the rolling bearing is used.The load distribution for the load cases-pure axial load and combined load with a high radial portion-is shown in Figure 8.It shows that under pure axial load (left), all rolling elements in the bearing are in the load zone.Thus, all contacts also contribute to the friction.In the load case with an additional radial load of 16 kN, the load zone shifts.It is now located in the lower area of the rolling bearing, which means that only 11 rolling elements are still carrying.This can be seen from the fact that only 12 points are still visible in the diagram.11 represents the load-carrying rolling elements; all other points are at 0 N since the remaining rolling elements do not experience any load.Accordingly, fewer contacts (effective length) contribute to the friction.It shows that under pure axial load (left), all rolling elements in the bearing are in the load zone.Thus, all contacts also contribute to the friction.In the load case with an additional radial load of 16 kN, the load zone shifts.It is now located in the lower area of the rolling bearing, which means that only 11 rolling elements are still carrying.This can be seen from the fact that only 12 points are still visible in the diagram.11 represents the loadcarrying rolling elements; all other points are at 0 N since the remaining rolling elements do not experience any load.Accordingly, fewer contacts (effective length) contribute to the friction.

Discussion
The preliminary work shows that the dynamics simulations with LaMBDA are well suited to simulate the frictional torque of TRBs at low oil levels or minimum quantities of lubrication.Hydraulic losses due to a higher oil level are not calculated by default with LaMBDA.However, equations exist to take the losses into account.The method, according to Liebrecht, among others, can be regarded as the state of the art [82].
First comparative calculations between LaMBDA simulations and measurements that can be found in the literature, i.e., performed by Liebrecht [82] and Gonda [83], are shown in Figure 9.The measurements had been taken at a test rig where the test bearing and the supporting bearings have separate oil reservoirs to allow the measurement of torques due to the latter bearings independently from the losses of the support bearing [82,83].While the friction torque simulated with LaMBDA (gray) with minimum quantity lubrication agrees well with the experimentally determined friction torques (light blue), the simulation with a fully flooded test volume (yellow) deviates strongly from the measured friction torque (orange).This means that the state-of-the-art approach (LaMBDA with Liebrecht) greatly underestimates the hydraulic losses in the bearing.It shows that under pure axial load (left), all rolling elements in the bearing are in the load zone.Thus, all contacts also contribute to the friction.In the load case with an additional radial load of 16 kN, the load zone shifts.It is now located in the lower area of the rolling bearing, which means that only 11 rolling elements are still carrying.This can be seen from the fact that only 12 points are still visible in the diagram.11 represents the load-carrying rolling elements; all other points are at 0 N since the remaining rolling elements do not experience any load.Accordingly, fewer contacts (effective length) contribute to the friction.

Discussion
The preliminary work shows that the dynamics simulations with LaMBDA are well suited to simulate the frictional torque of TRBs at low oil levels or minimum quantities of lubrication.Hydraulic losses due to a higher oil level are not calculated by default with LaMBDA.However, equations exist to take the losses into account.The method, according to Liebrecht, among others, can be regarded as the state of the art [82].
First comparative calculations between LaMBDA simulations and measurements that can be found in the literature, i.e., performed by Liebrecht [82] and Gonda [83], are shown in Figure 9.The measurements had been taken at a test rig where the test bearing and the supporting bearings have separate oil reservoirs to allow the measurement of torques due to the latter bearings independently from the losses of the support bearing [82,83].While the friction torque simulated with LaMBDA (gray) with minimum quantity lubrication agrees well with the experimentally determined friction torques (light blue), the simulation with a fully flooded test volume (yellow) deviates strongly from the measured friction torque (orange).This means that the state-of-the-art approach (LaMBDA with Liebrecht) greatly underestimates the hydraulic losses in the bearing.The equations formulated by Liebrecht are derived from experiments with a tapered roller bearing of type 32208 with vertical shaft alignment.They take into account, in addition to the operating conditions, the flow characteristics and simplify the internal bearing geometry.For the operating points investigated, they are simplified as follows: While length l RB and diameter d RB of the rolling elements and as well as raceway diameter d ORL resp.d IRL , effective surface A ORL resp.A IRL and d m the mean rolling bearing diameter are geometrical values, the oil density ρ, the dynamic viscosity η and the kinematic viscosity v are lubricant properties.The velocities of the single components are taken into account via the cage speed n C .The oil level is taken into account indirectly via the effective surfaces.
Better results for hydraulic losses can be obtained from the studies in [84] using a Computational Fluid Dynamics (CFD) model developed in an OpenFOAM ® environment.The highly specialized CFD simulation obtained agrees with the experiment with only insignificant deviations (see Figure 9 on the right side).In addition, the aforementioned simulations with OpenFOAM ® have been validated with respect to lubricant flows in a 32312-A TRB [85][86][87].Combining the calculated losses from both LaMBDA and CFD simulations, the total bearing losses can be predicted very well (dark blue curve, Figure 9, left).These results show good agreement with the measured total bearing losses (orange).The comparison of the listed calculation methods shows that the state-of-the-art approach predicts the hydraulic losses in a bearing only to a certain extent.With highly specialized CFD simulations, it is possible to determine these losses very precisely [88].
However, there is a need for further research to combine the two methods in order to determine the influence of hydraulic losses on the bearing kinematics.For this purpose, generally applicable equations would have to be derived from the CFD simulation, with which the flow-dependent forces of the lubricant on the bodies in the rolling bearing can be calculated.The forces can then be taken into account in the dynamics simulation and included in the force equilibrium so that they influence the movement-kinematics-of the roller bodies and cage.

Conclusions
The developed dynamic model LaMBDA for tapered roller bearings provides very good results with respect to the simulated frictional torque.This has been shown by comparison with friction torque measurements on a single-bearing test rig.Lubricant properties depending on temperature, load, and stresses are correctly represented in the model.The highly accurate contact models and friction descriptions are the basis for the high expressiveness of the model.
Nevertheless, not all losses that occur in a rolling bearing are taken into account.If rolling bearings in a gearbox are to be simulated, the losses due to the oil outside the immediate contact must also be included.Further research is needed in this regard.As shown in the discussion, CFD simulations can make an important contribution to the description of these losses.CFD simulations can also provide information on the lubricant pumping mechanism of the TRB and the lubrication supply, which are not considered in the model yet.A mathematical-physical description of the losses can be derived from the complex flow simulations.In multi-body simulation, these approaches can contribute to the holistic system optimization of gearboxes.Studies on the implementation of these losses in the MBS require further investigation.

Figure 1 .
Figure 1.Flow chart of contact force calculation between rolling elements and ring raceways [26].

Figure 1 .
Figure 1.Flow chart of contact force calculation between rolling elements and ring raceways [26].

→p
, contact normal vector → n and penetration δ.With a load-deformation relationship, the effective contact normal force → F N .is determined from the penetration.In the next step, the time-dependent contact state parameters are calculated.These are the relative velocities → u rel and sum velocities → u sum in the contact point.The relative velocity is further used to calculate the damping forces → F D as described in Section 2.1.2.The normal forces, velocities, and lubricant data are used to determine the lubricant film height h in contact, the specific lubricant film height Λ and the solid load-bearing ratio φ.On their basis frictional forces → F T and resulting torques → M T are calculated.The calculation of frictional forces for each contact in the model is presented in detail in Section 2.1.3.In the last step, the three components of the contact force → F Σ are summed up and given back to the MBS solver for iteration of the force equilibrium.

Figure 2 .
Figure 2. Friction calculation in roller raceway contact of the TRB MBS model [26].

Figure 2 .
Figure 2. Friction calculation in roller raceway contact of the TRB MBS model [26].

Figure 3 .
Figure 3. Scematic view and CAD model of the friction torque test bench at MEGT [42].

Figure 3 .
Figure 3. Scematic view and CAD model of the friction torque test bench at MEGT [42].

Figure 4 .
Figure 4. Graphical illustration of the geometry parameters of the TRB used in this study.On the left the sectional view of the TRB and on the right a rolling element is shown.The dotted lines represent the symmetry lines of the bearing and the rolling element.

Figure 4 .
Figure 4. Graphical illustration of the geometry parameters of the TRB used in this study.On the left the sectional view of the TRB and on the right a rolling element is shown.The dotted lines represent the symmetry lines of the bearing and the rolling element.

Figure 4 .
Figure 4. Graphical illustration of the geometry parameters of the TRB used in this study.On the left the sectional view of the TRB and on the right a rolling element is shown.The dotted lines represent the symmetry lines of the bearing and the rolling element.

Figure 5 .
Figure 5.Comparison of the measured and LaMBDA-calculated frictional torque of a TRB of type 32216 under purely axial load of 6 kN and oil bath lubrication with reference oil FVA No. 3.

Figure 5 .
Figure 5.Comparison of the measured and LaMBDA-calculated frictional torque of a TRB of type 32216 under purely axial load of 6 kN and oil bath lubrication with reference oil FVA No. 3.

Lubricants 2023 , 22 Figure 6 .
Figure 6.Comparison of the measured and LaMBDA-calculated frictional torque of a TRB of type 32216 under combined axial load of 6.5 kN and radial load of 6 kN with oil bath lubrication with reference oil FVA No. 3.

Figure 6 .
Figure 6.Comparison of the measured and LaMBDA-calculated frictional torque of a TRB of type 32216 under combined axial load of 6.5 kN and radial load of 6 kN with oil bath lubrication with reference oil FVA No. 3.

Figure 6 .
Figure 6.Comparison of the measured and LaMBDA-calculated frictional torque of a TRB of type 32216 under combined axial load of 6.5 kN and radial load of 6 kN with oil bath lubrication with reference oil FVA No. 3.

Figure 7 .
Figure 7.Comparison of the measured and LaMBDA-calculated frictional torque of a TRB of type 32216 under combined load and oil bath lubrication with reference oil FVA No. 3 under variation of the radial load at 50 °C.

Figure 7 .
Figure 7.Comparison of the measured and LaMBDA-calculated frictional torque of a TRB of type 32216 under combined load and oil bath lubrication with reference oil FVA No. 3 under variation of the radial load at 50 • C.

Figure 8 .
Figure 8. Load distribution of the TRB type 32216 used.Left: with purely axial load of 6 kN, right: with combined load of 6 kN axial and 16 kN radial.

Figure 8 .
Figure 8. Load distribution of the TRB type 32216 used.Left: with purely axial load of 6 kN, right: with combined load of 6 kN axial and 16 kN radial.

Figure 8 .
Figure 8. Load distribution of the TRB type 32216 used.Left: with purely axial load of 6 kN, right: with combined load of 6 kN axial and 16 kN radial.

Figure 9 .
Figure 9.Comparison of the simulated frictional torque of a TRB of type 32208 under oil sump lubrication with reference oil FVA No. 3 at 50 °C and an axial load of 1 kN with the measured frictional torque on a single-bearing test rig.Simulated frictional torque with LaMBDA and experiment (left), experiment and CFD simulation [84] (right).Bearing data is provided inTable A2 in Appendix B) and lubricant data in Table A1 in Appendix A).

Figure 9 .
Figure 9.Comparison of the simulated frictional torque of a TRB of type 32208 under oil sump lubrication with reference oil FVA No. 3 at 50 • C and an axial load of 1 kN with the measured frictional torque on a single-bearing test rig.Simulated frictional torque with LaMBDA and experiment (left), experiment and CFD simulation [84] (right).Bearing data is provided inTable A2 in Appendix B) and lubricant data in Table A1 in Appendix A).

Table 1 .
Test setup for friction torque measurement of a TRB type 32216 at pure constant axial load.

Table 1 .
Test setup for friction torque measurement of a TRB type 32216 at pure constant axial load.

Table 2 .
Test setup for friction torque measurement of a TRB type 32216 at combined, constant axial and radial load.

Table 3 .
Test setup for friction torque measurement of a TRB type 32216 at combined load.

Table 4 .
Geometrical data of a TRB type 32216.

Table A2 in
Appendix B) and lubricant data in Table A1 in Appendix A).

Table A2 .
Geometrical data of a TRB type 32208.
→ sDisplacement between two coordinate systems .→ s Relative velocity between two coordinate systems