© S.O. Farah et al., published by EDP Sciences 2020

1 Introduction

Industrial machines require a good lubrication regime with an appropriate lubricant to avoid wear mechanisms that could occur at the contact interface between two surfaces. In rotating machines, several components such as ball bearings, cams and gears operate under EHD lubrication regimes with non-conformal contacts. The nonconformity of the surfaces in contact means that the contact areas are subjected to a high contact pressure which can be in the range of several gigapascals. With a view to being able to assess whether the machine elements are working with safety conditions or are close to reaching boundary lubrication and ultimate failure, the fluid film thickness is consequently of a prime importance. Such mechanical components are often designed for the lubricant film to be thick enough, more than 1 μm, to avoid contact between the mating surfaces ensuring the expected bearing lifetime. For that purpose several experimental methods are dedicated to measuring the lubricant film thickness. The most widely used one is based on optical interferometry [1–3]. Nevertheless, the electrical technique has previously been used with measurements of electrical resistance [4–6] and capacitance [7–9]. The electrical resistance can indeed be useful when the contact between two surfaces is metallic which limits its use as a conducting contact interface. The electrical capacitance, on the other hand, is specifically linked to the contact interface thickness that accommodates the lubricant film. In the field of ball bearing diagnosis, the electrical techniques have been used to measure the load distribution by capacitance sensors [10] as well as to highlight the angular speed effect on the lubrication regime by means of the electrical resistance [11]. More recently, Jablonka et al. [12] proposed a quantitative evaluation of lubricant film thickness in a radially loaded deep-groove ball bearing using electrical capacitance measurements. The authors used in-situ monitoring tool to identify the lubrication conditions as function to the operating parameters and lubricant properties.

Regarding numerical modeling in the field of ball bearing there are basically two approaches to perform numerical computations, namely the finite element method (FEM) and the multibody method. Numerical tools allow analysis on specific ball bearing parameters, such as kinematics, dynamics, stresses, deformations, lubrication regime and lifetime. Because of the cyclic loadings, bearing rings are the most critical components since the flaking phenomena often appearing on rolling raceways is predictive of the onset of the damage mechanisms. In order to better understand the fatigue behavior of bearing rings, FEM was used to quantify the stress level and its distribution within the material according to the applied mechanical loading [13]. From a computational time point of view, however, FEM modeling remains costly due to the numerical management of the contact. Furthermore, FEM is often performed to describe only a portion of ball bearing. Unlike FEM, the multibody method offers the advantage of allowing a description taking into account the entire bearing, often with the exception of the cage, for modeling the loading distribution and lubrication effects [14,15]. Recently DEM has been used successfully for simulating the dynamic behavior of bearings with tribology enrichments [16] as well as for modeling lubrication regimes by means of an EHD model [11,17]. To achieve this, the digital twin of ball bearings implements a major part of the functionalities according to the state of the art, such as the slice model, the mutual influences of neighboring elements with clearance, elliptical Hertzian contact model, lubricant film estimation, damping forces, centrifugal forces, etc. In the following paper, the ball bearing is simulated with the assumption of EHD lubrication conditions. The coupling between the mechanical pressure and the dielectric behavior of the lubricant film thickness is carried out by means of the electrical capacitance. According to the model of a two parallel-plate capacitor separated by a dielectric material, the capacitance of the lubricant film thickness is extracted numerically. The influence of lubricants properties (viscosity and dielectric constant) on the capacitance is well highlighted. The numerical predictions of fluid film capacitance for the two simulated lubricants provide a promising result in comparison to the experimental data available in the bibliography [12].

The paper is divided into four sections including two sections dedicated to the introduction and conclusion, respectively. The second section deals with a three-dimensional discrete modeling performed on a ball bearing, with a digital twin called MULTICOR3D. The simulated radial ball bearing of 6208 series is designed fully and accurately in accordance with its geometrical conformity. The third section focuses on the mechanical modeling of contact forces based on a smoothed formulation of DEM in the case of an elliptical Hertzian contact with the assumption of piezo-viscous-elastic behavior of the lubricant film. The electrical capacitance of the lubricant is then introduced by means of the parallel-plate capacitor formula depending on the contact pressure and fluid film thickness. The next section deals with the numerical simulations where the effects of operating conditions on the fluid film capacitance are highlighted and discussed. Some comparisons with experimental measurements extracted from the bibliography are also discussed in this section. In the conclusion, we summarize the contributions of this work and we end with a few thoughts on future work.

2 Digitalt win of radial ball bearing

A three-dimensional digital twin of ball bearing (Fig. 1a) is developed by means of MULTICOR3D software using DEM [17]. The radial ball bearing of 6208 series (SKF 6208 TN9/DBGA) has been considered in the present study (Fig. 1b).

The dimensions and geometrical conformity of the 6208 radial ball bearing are provided in Table 1. The number of rolling elements which are likely to be in contact with the inner or outer raceways is Z = 9. Furthermore, nine additional smaller ballsof radius R_z in contact with two toruses are considered to reproduce the cage effects (Fig. 2). Even if the cage is different in shape in comparison to the standard form, the design does not at all affect the contact forces between rolling elements and cage [17].

We note that this description of the cage mainly has two roles. The first is that it ensures a conforming internal radial load distribution,while the second is that it allows intermittent rolling elements/cage contacts to be taken into account when the ball bearing operates. However, details should be supplied about the cage in the sense that the only other possible interactions involving the cage balls may be occur with rolling elements. By considering all bearing components, the proposed modeling gives the opportunity to develop a digital twin without any restriction.

We note that the digital twin design performed by DEM allowed us to numerically achieve standard bearing frequencies [18]. For validation purpose, normal load distribution for several diametral clearances and free contact angle variation as a function of the raceway conformities are shown in Appendix A.

Fig. 1 (a) Digital twin of radial ball bearing 6208, (b) SKF 6208 radial ball bearing.

Fig. 2 Exploded view of digital twin (6208 radial ball bearing).

3 Mechanical modeling of contact forces

The contact stiffness model that we used is based on a smoothed formulation of DEM initially proposed by Cundall and Strack [19]. The contact interactions are obtained by means of Kelvin-Voigt Spring-Dashpot Model and a Coulomb’s law if sliding occurs at the contact (Fig. 3). The main assumption made with DEM is that the particles in contact are assumed perfectly rigid and therefore, there is no elastic deformation of all bearing components. The elasticity is only considered at the contact between, rolling element/raceways (inner or outer), cage ball/toruses and rolling element/cage ball (Fig. 2).

Fig. 3 Contact stiffness model.

3.1 Hertzian contact stiffness model

The contact forces F_n,t in normal and tangential directions according to the local frame at the contact plane, are described with an explicit mechanical model depending on elastic force displacement law, Coulomb friction and viscous damping coefficients (1). $$\text{{}{F}_{n,t}=K_{n,t}{\delta }_{n,t}+C_{n,t}{v}_{n,t}\text{\hspace{0.17em}with\hspace{0.17em}}{F}_t=-\text{min}\left(F_t,\mu F_n\right)\times sgn\left(v_t\right)$$(1)

with K_n,t the normal and tangential stiffnesses, C_n,t the normal and tangential viscous damping coefficients, δ_n,t the normal and tangential relative displacements, v_n,t the normal and tangential relative speeds according to the local frame and μ the Coulomb friction coefficient in the case of a sliding contact. The normal viscous damping coefficient C_n to ensure a mechanical steady state is a function of two damping types, namely hysteretic and EHD, related to the contact shape, roughness, fluid viscosity, temperature, entrainment speed, etc. The tangential viscous damping coefficient C_t is more difficult to evaluate. The ratio $\frac{C_t}{C_n}$ is usually considered to range within the interval [0, 1]. The normal contact stiffness K_n is related to the mechanical characteristics and curvatures of the surfaces in contact according to non-linear model; $F_n\propto {\delta }_n^{\frac{3}{2}}$. For a given rolling element in contact with the inner-raceway (or outer-raceway), with the assumption of an elliptical Hertzian contact, the normal stiffness is written as function of the curvature sum radius and the approximate elliptic integrals as introduced by Hamrock [20]: $$\{\begin{array}{l}{K}_n=\frac{2\pi \kappa G}{1-\nu }\sqrt{\frac{2\overline{E}{R}_{\text{curve}}}{9\overline{F}{ }^3}} {\delta }_n^{\frac{1}{2}}\hfill \\ R_{\text{curve}}=\frac{R_x{R}_y}{R_x+R_y}\hfill \\ R_x=\frac{R_{1x}{R}_{2x}}{R_{2x}\pm R_{1x}}, R_y=\frac{R_{1y}{R}_{2y}}{R_{2y}-R_{1y}}\hfill \end{array}$$(2)

where R_curve is the curvature sum radius, R_x the effective radius in x direction, R_y the effective radius in y direction, κ the ellipticity parameter $\left(\kappa ={\alpha }^{\frac{2}{\pi }}, \alpha =\frac{R_y}{R_x}\right)$, $\overline{F}$ the approximate elliptic integral of first kind $\left(\overline{F}=1+q \text{ln}\alpha , q=\frac{\pi }{2}-1\right)$ and $\overline{E}$ the approximate elliptic integral of second kind $\left(\overline{E}=1+\frac{q}{\alpha }\right)$. The equivalent radius R_x is calculatedby taking into account the cases of convex contact (rolling element/inner-raceway) and concave contact (rolling element/outer-raceway), as given in equation (2). For the sake of clarity, we note that R_b =R_1x = R_1y, R_2x = R_i and R_2y = R_r (Fig. 4). $G=\frac{E}{2\left(1+\nu \right)}$ is the shear modulus of steel in the case of identical materials in contact, E the Young’s modulus (E = 210 GPa) and ν the Poisson’s ratio (ν = 0.3).

When interactions between rolling element and cage ball occur, the contact is assumed to be without friction and the damping coefficient is considered critical. The normal stiffness (3) only is taken into account, and therefore the following equation is considered [21]: $$K_n=\frac{4E_{\text{eq}}\sqrt{R_{\text{eq}}}}{3\left(1-\nu \right)} {\delta }_n^{\frac{1}{2}}$$(3)

with E_eq the equivalent Young’s modulus $\left(\frac{1}{E_{\text{eq}}}=\frac{1-{\nu }^2}{E}+\frac{1-{\nu }_c^2}{E_c}\right)$, E_c the Young’s modulus (E_c = 3.3 GPa) and ν_c the Poisson ratio (ν_c = 0.41) of cage balls made with polyamide 66 (PA66), $R_{\text{eq}}=\frac{r_i{r}_j}{r_i+r_j}$ the equivalent radius of rolling element and cage ball in contact of radius r_i   and r_j, respectively.

To achieve a steady regime of the ball bearing under operating conditions, viscous damping and local friction are added as dissipative forces to the mechanical model (Figure 1). When a contact occurs between rolling element and raceways (inner or outer) two kinds of normal viscous damping could be involved: hysteritic and(or) fluid. The viscous effect due to hysteretic damping, which is related to the plastic strain of surface roughness (4), is expressed by the coefficient C_hyst in the normal direction $\vec{n}$. $$C_{\text{hyst}}=\frac{3{\alpha }_e{K}_n{\delta }_n}{2}$$(4)

where α_e denotes the restitution coefficient of steel [22]. A value of 0.08 s.m⁻¹ is considered to compute the hysteritic damping C_hyst.

Regarding the viscous damping due to the fluid film, Walford and Stone [23] established a damping coefficient C_fluid related to the lubrication film in the normal direction of contact (5). Since the bearing contains non-conformal contacts endowed with EHD lubrication regime, the lubricant film must be thick enough to avoid contact between the mating surfaces. $$C_{\text{fluid}}=\frac{3\pi \eta a}{\sqrt{2}}{\left(\frac{R_b}{h_{\text{min}}}\right)}^{\frac{3}{2}}$$(5)

where η is the dynamic viscosity of the fluid at operating temperature and atmospheric pressure [24], a the semi-major axis of the Hertz elliptical contact [25], h_min the minimum lubricant film thickness of the ball-inner-raceway (or -outer-raceway) contact.

The DEM modeling is by solving Newton’s second law for each rolling element taking into account dynamic effects, such as centrifugal forces and gyroscopic effects, when the bearing is functioning. The mechanical resolution requires an explicit time integration based on the velocity-Verlet scheme [17].

Fig. 4 Curvatures in contact in two orthogonal cross sections of a ball bearing: (a) x−z plane, (b) y−z plane.

3.2 EHD lubrication model

In the field of EHD lubrication, oil (or grease) rheology influences the fluid film thickness under operating conditions [26]. Given that the viscosity of a lubricant increases with respect to pressure, this means that the lubricant exhibits piezo-viscous behavior. Furthermore, temperature also plays a role in lubricant viscosity changes with the consequence of diminishing its viscosity. Starting from these assumptions with a piezo effect of the lubricant, we have implemented an EHD model in the digital twin with piezo-viscous-elastic regime to predict the minimum fluid film thickness h_min. As initially proposed by Hamrock and Anderson [20], the general form of the EHD model is a function of parameters related to the speed, load and material properties. Starting from the previous studies in the field of EHD lubrication, Masjedi and Khonsari [27] improved relatively recently the EHD model by introducing the ellipticity effect. Since changing bearing conformity affects the contact ellipticity, we have implemented the new EHD formulation in the digital twin in which the parameters of the model are a function of the ellipticity. In the case of smooth point-contact EHD, the expression of the minimum fluid film thickness (6) is as follows: $$\begin{array}{lll}\frac{h_{\text{min}}}{R_x}\hfill & =\hfill & 1.637 {\overline{U}}_r^{ 0,711{\kappa }^{-0,023}} {\overline{E}}^{ 0.65{\kappa }^{-0.045}} {\overline{W}}^{ -0.09{\kappa }^{-0.15}}\hfill \\ \hfill & \hfill & \times \left(1-0.974e^{-0.676\kappa }\right)\hfill \end{array}$$(6)

where $\overline{E}, {\overline{U}}_r, \overline{W}$ are three dimensionless parameters: $\overline{E}=\xi E^{\prime }$, with ξ the pressure-viscosity coefficient of the lubricant at operating temperature and atmospheric pressure, $E^{\prime }=\frac{E}{1-{\nu }^2}$ the effective elastic modulus of steel; ${\overline{U}}_r=\frac{\eta U_r}{E^{\prime }{R}_x}$, with U_r the entrainment speed; $\overline{W}=\frac{Q_{\psi }}{E^{\prime }{R}_x^2}$, with Q_ψ the ball-inner-raceway (or -outer-raceway) normal load at an azimuth angle ψ (Fig. A.1).

It is commonly agreed that the identification of the lubrication regime in which a system is running is by the use of the dimensionless fluid parameter Λ_f such defined by (7). Λ_f is the ratio of the minimum fluid film thickness to the surface roughness [28]: $${\Lambda }_f=\frac{h_{\text{min}}(Q_{\psi },\omega )}{\sqrt{{\sigma }_1^2+{\sigma }_2^2}}$$(7)

where σ₁ and σ₂ denote the standard mean roughnesses of the measured surfaces for the inner-raceway (or outer-raceway) and rolling elements, respectively. However, it is admitted that only the raceway roughness is considered, and the rolling element roughness can then be neglected $\left({\sigma }_2\ll {\sigma }_1\right)$. In addition, the schematic representation of the Stribeck curve (Fig. 5) for non-conformal contact shows the dependence of the friction coefficient μ on the dimensionless fluid parameter Λ_f. The friction coefficient μ is therefore dependent on the lubrication regime [29]. In addition, the Stribeck curve mainly defines three lubrication regimes (boundary, mixed-film or fully-fluid film) in relation to fluid parameter Λ_f.

Fig. 5 Schematic representation of Stribeck curve: friction coefficient and fluid film thickness as function of fluid parameter.

3.3 Electrical capacitance model for lubricated Hertzian contact

In the context of the identification of a lubrication film by means of electrical measurements, the fluid film thickness can easily be extracted from the capacitance of the contact [30]. According to the parallel-plate formula the lubricant film thickness is inversely proportionalto the capacitance (8). $$C_{\text{Hertz}}={\epsilon }_0{\epsilon }_r\frac{A_{\text{Hertz}}}{h_c}$$(8)

with C_Hertz the Hertzian contact capacitance between two parallel conducting plates of area A_Hertz separated by a dielectric material where ε_r is the dielectric constant of the fluid and h_c its central thickness. ε₀ is the vacuum permitivity $\left({\epsilon }_0=8.854 \text{pF}.{\text{m}}^{-1}\right)$. The theoretical form of the contact area between the rolling element and raceway for a perfect elliptical Hertzian contact is A_Hertz = πab, with a and b the major and minor semi-axes of the ellipse, respectively. In the case of a smooth surface [27] the central film fluid thickness h_c integrating the ellepticity effect is as follows (9): $$\begin{array}{lll}\frac{h_c}{R_x}\hfill & =\hfill & 3.672 {\overline{U}}_r^{ 0,663{\kappa }^{-0,025}} {\overline{E}}^{ 0.502{\kappa }^{-0.064}} {\overline{W}}^{ -0.045{\kappa }^{-0.18}}\hfill \\ \hfill & \hfill & \times \left(1-0.573e^{-0.74\kappa }\right).\hfill \end{array}$$(9)

From experimental stand point, the qualitative conversion of the measured capacitance to lubricant film thickness involves first the estimation of the capacitance related to the nearby area outside the Hertzian contact and some background capacitance. The latter two are then subtracted from all measurement data. The remaining value corresponds to the contact capacitance due to the fluid film thickness [31]. From a numerical point of view, the contact capacitance is easily computed from expression (8) after predicting the central fluid film thickness by (9).

The dielectric constant ε_r of the lubricants varies with Hertz contact pressure. For nonpolar lubricants ε_r can be calculated from Clausius-Mossotti model [32]: $$\frac{{\epsilon }_r-1}{{\epsilon }_r+2}=\frac{N_A\beta \rho }{3M}$$(10)

where N_A is the Avogadro number, β the molecular polarizability, M the molecular weight and ρ the density. Thefactor $\frac{N_A\beta }{3M}$ can be deduced by considering permittivity and density at atmospheric pressure. According to the formula established by Dowson and Higginson [33], the variation of the lubricant density as a function of the contact pressure is given by following equation (11): $$\rho ={\rho }_0\frac{0.59\cdot {10}^9+1.34\cdot P_{\text{Hertz}}}{0.59\cdot {10}^9+P_{\text{Hertz}}}$$(11)

with ρ₀ the ambient-pressure density and P_Hertz the Hertz contact pressure. In Figure 6a are plotted the dielectric constant evolution of synthetic oil film (Tab. Table 2) with respect to Hertz contact pressure at inner and outer raceways. This can easily be understood by the fact that the higher a pressure, the greater the dielectric constant of the lubricant film. As depicted by the curves in Figure 6b, the maximum Hertz contact pressure at the inner raceway is greater than reached at the outer raceway due to the difference of the ellepticity of both inner and outer raceways. The contact pressures are simulated under dynamic operating conditions, with radial force Fr = 5 kN and shaft angular speed ω = 477.5 rpm.

From electrical point of view, a single rolling element in contact with raceways can be considered as an equivalent circuit composed ofcapacitances between the inner and outer raceway contacts, which are in series. It should be pointed out that the cage components (PA66) are insulating which allows the ball bearing to be considered as a set of capacitances in parallel (Fig. 7). The total capacitance of the ball bearing is calculated by summing the equivalent contact capacitances for all rolling elements located in the loaded zone. $$C_{\text{bearing}}={\sum }_{i=1}^{N_l}\left(\frac{C_{\text{Hertz}}^{\text{inner},i}\times C_{\text{Hertz}}^{\text{outer},i}}{C_{\text{Hertz}}^{\text{inner},i}+C_{\text{Hertz}}^{\text{outer},i}}\right)$$(12)

with N_l the number of rolling elements in the loaded zone of the ball bearing. The contact capacitances $C_{\text{Hertz}}^{\text{inner},i}$ and $C_{\text{Hertz}}^{\text{outer},i}$ between the ith rolling element and inner-outer raceways, respectively, are determined according to the parallel-plate formula (8). It must be emphasized that, unlike the case of experimental measurements, the capacitance model implemented in the digital twin does not take into consideration the contribution of capacitance between inner and outer ring and some background capacitance [31,34].

Fig. 6 (a) Dielectric constant as function of Hertz contact pressure (6208 bearing – PAO VG48 at 40 °C). (b) Hertz contact pressure as function of azimuth angle ψ.

Fig. 7 Flowchart of equivalent circuit with parallel capacitance of ball bearing.

4 Numerical simulations with the digital twin of ball bearing

This section is dedicated to numerical simulations performed on the digital twin of ball bearing. The correlation between contact capacitance and fluid film thickness is investigated when the ball bearing operates under imposed mechanical loading and fixed kinematic conditions. Sensitivity analyses are performed on the main driven parameters, namely the radial load F_r, the diametral clearance P_d and the angular speed ω of the shaft on which the inner ring is mounted. We recall that for benchmarking purpose some numerical simulations considering radial loading, free contact angle and free endplay in ball bearing have been carried out and given in Appendix A.

4.1 Effects of radial load and diametral clearance on electrical capacitance

According to the parallel-plate capacitor formula the capacitance of the film fluid thickness is extracted by considering a single rolling element after one revolution. Simulations are then carried out in order to exhibit the influence of the radial load F_r and diametral clearance P_d on the capacitance of the fluid film thickness. Comparisons with respect to experimental measurements [12] are performed in order to highlight the performances of the digital twin of ball bearing. Two kinds of lubricants are considered based on the experimental study given in [12], namely mineral oil and polyalphaolefin synthetic oil (PAO). In Table 2 are summarized the physical properties of both lubricants.

Figure 8a depicts the sensitivity of the capacitance to the radial load F_r ranging over the interval [1–7 kN] in the case of mineral oil. The capacitance is proportional to F_r, which means that the higher the radial load F_r, the greater the change in fluid film capacitance. From a numerical point of view, when the rolling element loses contact with the rings of the ball bearing, the capacitances in series related to the inner and outer raceways are disabled (Fig. 7). In addition, as can be observed in Figure 8a, the predicted capacitance for F_r = 5 kN is in quite good agreement with the experimental measurements obtained with the same entrainment speed in the case of 6306 ETN9 deep groove bearing with a polymer cage [12].

For comparison purpose we have added to our numerical results the capacitance related to the background and the capacitance due to the space between inner and outer raceways. The equivalent capacitance due to these latter effects is about 30 pF according tothe experimental results given in [12]. Therefore, the additional capacitance is assumed to be of the same order of magnitude as for the simulated ball bearing 6208 series (SKF 6208 TN9/DBGA). Good agreement with experimental measurements can also be observed on the curve fit applied to the numerical predictions for a zero azimuth angle ψ (Fig. 8b) where the radial load reaches its maximum. This confirms that the digital twin of ball bearing is able to predict fluidfilm capacitance fairly close to experimental data.

The next simulation deals with the effect of the load parameter ε on the capacitance distribution in the ball bearing by considering several values of $\frac{P_d}{D}$ ratio. The capacitance distribution in the loaded zone is plotted for three configurations, typically with negative clearance, zero clearance and positive clearance. The values of $\frac{P_d}{D}$ ratio that have been considered are summarized in Table 3.

The effect of the diametral clearance P_d is well highlighted in Figure 9. Because there are more contacts between elements and raceways in the case of P_d < 0, this increases the number of capacitances operating in the electrical model (Fig. 7). As can be observed in Figure 9, for $\frac{P_d}{D}$ ratio of − 0.0035, corresponding to a load parameter ε of 1.32, the capacitances of the electrical model are all activated. Conversely, in the case of positive clearance more than 50% of the capacitances are disabled because of the reduced size of loaded zone (Fig. 9).

Fig. 8 (a) Capacitance distribution in the loaded zone with $\frac{P_d}{D}=0$ and Ur = 0.41 m/s. (b) Capacitance at the maximum radial load position $\left(\psi =0^{°}\right)$ as function of Fr.

Fig. 9 Capacitance distribution in the loaded zone with Fr = 5000 N and Ur = 0.41 m/s.

4.2 Capacitance variation with respect to entrainment speed and fluid film thickness

In this part dedicated to the study of capacitance variation as a function of the entrainment speed U_r and the fluid film thickness h_min, we have considered two radial loads F_r = 1 kN and F_r = 5 kN. The diametral clearance P_d is kept constant and equal to zero $\left(\epsilon =0.5\right)$. First of all, we have plotted the capacitance variation with respect to U_r for the two lubricants listed in Table 2. We note a sharp variation in capacitance for low rolling speed under a given radial load F_r. Clearly, there is a strong correlation between the capacitance and the contact area (F_r) on the one hand and fluid film thickness (U_r) on the other. In addition, an increase in U_r carries a slight decrease in the capacitance (Fig. 10a and b). We note that the capacitance values for the synthetic lubricant are higher then those for the mineral lubricant. The difference between the electrical behaviors of the two lubricants is mainly driven by the chemical composition of each oil. According to the experimental measurements performed on PAO6 (synthetic oil [35]), as expected, is characterized by a very low electrical conductivity, while mineral oil is a rather better electrical conductor.

Once again, we have carried out a comparison of the numerical simulations with the measurements resulting from the experimental study performed with 6306 ETN9 deep groove bearing [12]. According to the previous capacitance distributions ( Fig. 8a) plotted as function of the angular position of one rolling element, the capacitance value at an azimuth angle ψ =0, corresponding to the maximum normal load Q_ψ, provides more valuable indications about the lubrication regime. Figure 10a for mineral oil and Figure 10b for PAO VG48 exhibit quite a good agreement between numerical predictions and the experimental results. In the case of radial load F_r = 1 kN the capacitance predictions as function of entrainment speed are quite close to experimental one for both lubricants. However, for F_r = 5 kN, there are some discrepancies with experimental measurements in particular for PAO VG48 lubricant (Fig. 10b). This is related in part to the fluid film thickness difference at the contact between the rolling elements and raceways. In fact, the asperities effect is more significant at a low entrainment speed U_r combined with a high radial load F_r, which explains why the measured capacitance is higher than the numerical capacitance. In [27] the authors have suggested a parameter (asperity load ratio) in order to quantify the load carried by the asperities. Because of a better electrical metal conductivityin comparison to lubricant, the higher the asperity load ratio, the higher the contact capacitance. Furthermore, the effects of the outside zones close to the Hertzian contact [36] have to be taken into consideration for a better agreement of numerical predictions of contact capacitance with experimental results. As mentioned in [12], the background and ring capacitances are subtracted from all experimental measurements.

To avoid the effect of the lubricant properties, Figure 11 depicts the overall capacitance of the ball bearingas a function of the fluid film thickness at the contact between rolling element and inner raceway. Only the lubricant film thickness at the inner raceway has been considered because of a higher contact pressure. Assuming that the two radial loads provide curves with similar trends, only the case of F_r = 5 kN is simulated. As for the previous capacitance curves (Fig. 10a and b) a sharp decrease in capacitance is observed for low fluid film thickness, which corresponds to low entrainment speed U_r. It stands to reason that from a numerical point of view, it is easy to correlate the capacitance to the fluid film thickness by means of equation (6). In this sense, the digital twin gives benefit over an experimental measurement because it easily allows the capacitance change to be followed according to the lubrication regime when the bearing is operating.

Figure 12a and b depict the dependence of fluid parameter Λ_f on entrainment speed U_r and the correlation between Λ_f and the overall capacitance of the ball bearing. In the interval of entrainment speeds displayed on the x-axis and under a radial load F_r = 5 kN, the lubrication regime in the ball bearing skips from boundary to mixed according to the Stribeck curve (Fig. 5). For a given entrainment speed U_r, we clearly observe a difference in lubrication regime when the ball bearing operates with mineral oil and when it operates with synthetic lubricant PAO (Fig. 12a). Because the two tested lubricants have different viscosities, the ball bearing rapidly reaches the mixed regime with mineral oil in comparison to the PAO. This helps in predicting suitable operating conditions with an appropriate lubricant in order to increase the ball bearing lifetime. In Figure 12b the intervals of boundary and mixed regimes are plotted with respect to the ball bearing capacitance. For the considered ball bearing and radial loading, the change from boundary to mixed regime occurs when capacitance reaches 60 pF.

Fig. 10 Capacitance as function of entrainment speed Ur: (a) mineral oil, (b) synthetic oil PAO VG48.

Fig. 11 Capacitance-dielectric constant ratio as function of fluid film thickness of the element-inner raceway contact.

Fig. 12 (a) Film parameter Λf as function of entrainment speed Ur. (b) Correlation of the fluid parameter Λf and ball bearing capacitance.

5 Conclusion

Throughout this work we have attempted to highlight the relevance of the digital twin when a ball bearing operates under EHD lubrication regimes. A mechanical modeling based on DEM has been carried out with the aim of providing the realistic behavior of a ball bearing in terms of kinematic parameters, contact stiffness and lubricant film effects. By means of a parallel-plate capacitor formula the capacitance of the fluid film thickness is numerically extracted with respect to the EHD lubrication regimes. Furthermore, qualitative and quantitative comparisons of the numerical predictions with experimental data highlight the interest of the digital twin for industrial applications. The digital twin could indeed be a useful predictive tool for recommended lubrication regimes in bearings with a view to enhancing the performance of industrial machines. Nevertheless, the proposed modeling could be improved in order to take into account the effect of asperities in the EHD formulation.

Nomenclature

α_e: Restitution coefficient of steel (s.m⁻¹)

δ_r: Ring radial shift (μm)

δ_n,t: Normal/tangential relative displacement (μm)

η: Dynamic viscosity of the fluid (m.Pa.s)

Λ_f: Fluid number

μ: Friction coefficient

ν: Poisson’s ratio of steel

ω: Shaft/inner-ring angular speed (rad.s⁻¹)

$\overline{E}$: Approximate elliptic integral of second kind

$\overline{F}$: Approximate elliptic integral of first kind

$\overline{E}$: Dimensionneless material parameter

${\overline{U}}_r$: Dimensionneless entrainment speed

${\overline{W}}_r$: Dimensionneless load parameter

ψ: Azimuth angle (rad)

ρ: Density of the lubricant

ρ₀: Ambient-pressure density of the lubricant

σ_1,2: Arithmetic mean roughnesses (μm)

ε: Load parameter

ε₀: Vacuum permitivity (pF.m⁻¹)

ε_r: Dielctric constant

ξ: Viscosity-pressure coefficient of the lubricant (GPa⁻¹)

a: Semi-major axis of the Hertz’s elliptical contact (mm)

b: Semi-minor axis of the Hertz’s elliptical contact (μm)

C_n,t: Normal/tangential viscous damping coefficient (N.s.m⁻¹)

D: Rolling element diameter (mm)

E: Young’s modulus of steel (GPa)

E′: Effective elastic modulus of steel (GPa)

F_r: Radial load (F_r)

G: Shear modulus of steel (GPa)

h_c: Central lubricant film thickness (nm)

h_min: Minimum film fluid thickness (nm)

J_r: Radial integral

k: Ellipticity parameter

K_n,t: Contact normal/tangential stiffness (N.m⁻¹)

P_d: Diametral clearance (μm)

Q_max: Maximum normal load at ψ = 0 (N)

R_b: Rolling element radius (mm)

$R_c^i$: Curvature radius of inner raceway (mm)

$R_c^o$: Curvature radius of outer raceway (mm)

R_i: Inner ring radius (mm)

R_o: Outer ring radius (mm)

R_z: Cage ball radius (mm)

R_curve: Curvature sum radius (mm)

U_r: Entrainment speed (m.s⁻¹)

v_t: Sliding velocity (m.s⁻¹)

Z: Number of rolling elements

Appendix A Normal load and axial free endplay

According to the standard theory, the mechanical state of a ball bearing under a radial load F_r is provided by a dimensionless load parameter, $\epsilon =\frac{1}{2}\left(1-\frac{P_d}{2{\delta }_r}\right)$, related to the diametral clearance P_d and the radial shift δ_r of the inner raceway (Fig. A.1).

The normal load Q_ψ due to radial load F_r at an azimuth angle ψ is function of the load parameter ε such that [37]: $$Q_{\psi }=Q_{\text{max}}{\left[1-\frac{(1-\text{cos}\psi )}{2\epsilon }\right]}^{3/2}$$(A.1)

with Q_max the maximum normal load at ψ = 0. Stribeck [29] proposed that Q_max is function of the radial integral J_r and the number of rolling elements $Z\left(Q_{\text{max}}=\frac{F_r}{Z{J}_r\left(\epsilon \right)}\right)$.

Figure A.2 exhibits the normal load Q_ψ plotted with respect to rolling element position at an azimuth angle ψ. Three values of load parameter ε are considered, 1.7, 0.5 and 0.29, which corresponds to negative clearance $\left(\frac{P_d}{D}=-0.0045\right)$, zero clearance $\left(\frac{P_d}{D}=0\right)$ and positive clearance $\left(\frac{P_d}{D}=0.0045\right)$, respectively. The azimuth angle ψ of the rolling element is varied by driving the inner ring at a slow angular speed $\left(\omega =1 \text{rad}\times {\text{s}}^{-1}\right)$ to avoid any centrifugal effect. One can conclude that the digital twin of ball bearing achieves results in terms of radial load distribution very close to the analytic solution given in [37].

Fig. A.1 Radial shift δr due to the imposed radial load Fr.

Fig. A.2 Variation of $\frac{Q_{\psi }}{F_r}$ ratio as function of azimuth angle ψ under radial load Fr = 5000 N for three values of load parameter ε.

Fig. A.3 Schematic representation of the axial free play.

The next comparison with the analytic solution was made on the free contact angle and free endplay. In fact a radial ball bearing is commonly designed with a diametral clearance which causes both radial and axial degrees of freedom. Under a free axial thrust the rolling element and raceway contact become oblique [37]. Figure A.3 shows the effect of axial play on radial ball bearing when the inner ring moves to the right. The distance between the centers Oⁱ and O^o of the inner and outer raceway curvatures respectively, is: $$C=\underset{B}{\underbrace{\left(\frac{R_c^i+R_c^o}{D}-1\right)}}\times D$$(A.2)

with B the total curvature of the bearing. We note that the numerical simulations have been carried out with identical conformities of inner and outer races which yields total curvature: $B=\left(\frac{2\times R_c^{i,o}}{D}-1\right)$. From Figure A.3, one can express the free contact angle α⁰ made by thevertical line and the line passing through the points of contact between the rolling element and both raceways as follows: $${\alpha }^0={\text{cos}}^{-1}\left(1-\frac{P_d}{2C}\right).$$(A.3)

As mentioned before, because of diametral clearance P_d the bearing is free to move axially without mechanical loading. One can define the free endplay as the maximum relative axial displacement of the outer raceway with respect to the inner raceway: $$P_e=2C\times \text{sin}{\alpha }^0.$$(A.4)

Numerical simulations with digital twin of ball bearing are conducted with several values of parameter B summarized in Table A.1. Raceway curvatures $R_c^{i,o}$ and conformities $\frac{R_c^{i,o}}{D}$ are deduced from the expression of the total curvature B.

As shown in Figure A.4, the simulations achieve results in good agreement with analytical model [37]. In conclusion, the benchmarking performed between the theoretical formulation and numerical simulations shows that DEM is well adapted to ball bearing modeling and provides accurate results.

Fig. A.4 Free contact angle and free endplay for several values of total curvature B.

References

All Tables

All Figures

	Fig. 1 (a) Digital twin of radial ball bearing 6208, (b) SKF 6208 radial ball bearing.
In the text

	Fig. 2 Exploded view of digital twin (6208 radial ball bearing).
In the text

	Fig. 3 Contact stiffness model.
In the text

	Fig. 4 Curvatures in contact in two orthogonal cross sections of a ball bearing: (a) x−z plane, (b) y−z plane.
In the text

	Fig. 5 Schematic representation of Stribeck curve: friction coefficient and fluid film thickness as function of fluid parameter.
In the text

	Fig. 6 (a) Dielectric constant as function of Hertz contact pressure (6208 bearing – PAO VG48 at 40 °C). (b) Hertz contact pressure as function of azimuth angle ψ.
In the text

	Fig. 7 Flowchart of equivalent circuit with parallel capacitance of ball bearing.
In the text

	Fig. 8 (a) Capacitance distribution in the loaded zone with $\frac{P_d}{D}=0$ and Ur = 0.41 m/s. (b) Capacitance at the maximum radial load position $\left(\psi =0^{°}\right)$ as function of Fr.
In the text

	Fig. 9 Capacitance distribution in the loaded zone with Fr = 5000 N and Ur = 0.41 m/s.
In the text

	Fig. 10 Capacitance as function of entrainment speed Ur: (a) mineral oil, (b) synthetic oil PAO VG48.
In the text

	Fig. 11 Capacitance-dielectric constant ratio as function of fluid film thickness of the element-inner raceway contact.
In the text

	Fig. 12 (a) Film parameter Λf as function of entrainment speed Ur. (b) Correlation of the fluid parameter Λf and ball bearing capacitance.
In the text

	Fig. A.1 Radial shift δr due to the imposed radial load Fr.
In the text

	Fig. A.2 Variation of $\frac{Q_{\psi }}{F_r}$ ratio as function of azimuth angle ψ under radial load Fr = 5000 N for three values of load parameter ε.
In the text

	Fig. A.3 Schematic representation of the axial free play.
In the text

	Fig. A.4 Free contact angle and free endplay for several values of total curvature B.
In the text