Alternative calculation on transient elasto-hydrodynamic lubrication

Shaocheng Zhu (State Key Laboratory of Traction Power, Southwest Jiaotong University, Chengdu, China and University of Lyon, INSA-Lyon, Villeurbanne, France)
Weihua Zhang (State Key Laboratory of Traction Power, Southwest Jiaotong University, Chengdu, China)
Daniel Nelias (University of Lyon, INSA-Lyon, Villeurbanne, France)

Industrial Lubrication and Tribology

ISSN: 0036-8792

Publication date: 12 March 2018

Abstract

Purpose

The purpose of this study is to propose a new method to solve transient elasto-hydrodynamic lubrication (EHL) problem.

Design/methodology/approach

First, the steady-state EHL solution is modified so that the elastic deformation theory is combined with oil film stiffness distribution instead of steady-state Reynolds equation. Second, subsequent dynamic EHL procedure develops, recursively using transient distributed oil film stiffness and damping, where each time-marching solution is iteratively searched by ensuring both oil film force growth and elastic deformation update for each load increment.

Findings

This method increases calculation speed and provides both distributed EHL stiffness and damping for transient regimes.

Originality/value

This method is of interest for fast applications such as rolling bearings or gears.

Keywords

Citation

Zhu, S., Zhang, W. and Nelias, D. (2018), "Alternative calculation on transient elasto-hydrodynamic lubrication", Industrial Lubrication and Tribology, Vol. 70 No. 2, pp. 423-431. https://doi.org/10.1108/ILT-10-2017-0280

Download as .RIS

Publisher

:

Emerald Publishing Limited

Copyright © 2018, Emerald Publishing Limited


1. Introduction

Elasto-hydrodynamic lubrication (EHL) contact qualities in transmitting components (bearing, gear, etc.) and in special dynamic systems with open connections (wheel-rail contact, etc.), can be measured by stiffness and damping parameters. In a dynamic design, integrated (or lumped) parameters are often preferred instead of distributed values owing to their feasibility in multi-body dynamic system, and it must be also recognized that steady-state parameters are used exclusively, as transient values are technically hard to seize.

In the past decades, efforts have gone into improving the accuracy and efficiency of steady-state EHL calculation: direct iterative method for light load in point contact was introduced by Hamrock and Dowson (1976) and an inverse method was developed for heavy loads, first used in line contact by Dowson and Higginson (1966) and later in point contact by Evans and Snidle (1981). Since the Newton–Raphson method was used, relevant research studies in EHL problem have become more popular (Rohde and Oh, 1975; Oh and Rohde, 1977; Houpert and Hamrock, 1986). When a multigrid method was introduced by Lubrecht et al. (1986), incredibly high efficiency in calculation was obtained. Interest in a numerical technique also led to relevant developments of derived methods, such as multilevel, multi-integration method (Venner, 1991), fast Fourier transform (FFT) method (Wang et al., 2009), etc.

For what is mentioned above, It is difficult to apply any type of method to the dynamic systems, from transient EHL problem to subsequent application in vibration system.

Apart from full numerical methods, a simplified method neglecting the Poiseuille term used by Chang (2000), Messe and Lubrecht (2002) and others, brought some efficiency in analysis.

As for integrated EHL behaviors, other efforts were made. A rough and direct integration of Reynolds equation was proposed and was often used in parallel to Hertz relation in bearing dynamic calculation (Walford and Stone, 1983), (Hagiu and Gafitanu, 1997), (Dietl, 1997), (Wensing, 1998), (Nonato and Cavalca, 2010). Recently, another rough and direct integration method was also used; oil film parameters are integrated within different spaces by Wiegert et al. (2013a), Wiegert et al. (2013b), but in series with Hertz stiffness. Both the numerical method and the curve-fitting method were used to extract lumped parameters from distributed EHL parameters by Sarangi et al. (2004) and Chippa and Sarangi (2013), and until now, this could be the most complex way to get EHL damping and stiffness coefficients. Before getting transient lumped parameters, transient distributed parameters should be studied.

This paper is oriented toward the study of transient distributed EHL behaviors. Besides the mentioned manners, we have focused on improving the Reynolds equation. With distributed damping and stiffness, the distributed transient EHL behaviors will be also be described.

2. Decomposition of transient EHL force

This paragraph discusses in depth the interacting mechanism on the EHL surface and introduces new principles of simplification.

Generally, oil film stiffness and damping co-exist on EHL surfaces. Distributed behaviors are reduced to three concentrated interactions. The zero effect at exit zone on total load is due to cavitation caused by negative pressure. Oil film performs indifferently at the Hertz zone, as the central approach is dominated by the Hertz deformation. Neglecting two preceding weak interactions, the oil film effect on the contact load is described at the entry zone by a rough and direct integration of the Reynolds equation, in a parallel form of normal supporting load (Walford and Stone, 1983):

(1) Q=4η0RyuyLhmin,iso+6.66η0uZRy1.5Lhmin,iso1.5

The inlet-equivalent relation and the central hertz relation turn one EHL pair to damping and stiffness coefficients, having been widely used to perform bearing dynamic calculation or to formulate bearing stiffness or damping matrices in the past.

However, this simplification fails to consider the deformation outside the Hertz area. Materials compress at Hertz area (primary stiffness), but stretch at entry zone (extra damping produced by “inciting wings”). Grubin (1949) predicted the inlet oil film thickness according to this assumption.

Meanwhile, preceding analytical expression well matches highly loaded contacts, but it fails to serve weakly loaded contacts, as thought by Wiegert et al. (2013a). He proposed a new force model to apply in the oscillatory system with multiple EHL contacts by interpreting each contact as a hydrodynamic force in series with the Hertz force. Hydrodynamic force related to the cavitation position and its height was reported (Wiegert et al., 2013b). However, both film and Hertz coefficients were separately derived according to the Reynolds equation and the Hertz theory, this daring assumption is introducing a linear damping into elastic contact was also specially used in his paper.

Although two preceding methods are performed well, the surface loading is distributed around the entire deformed area, and the interaction between the elastic structure and the lubricant should be considered; this is why Grubin’s prediction was incomplete. Until a full numerical method is introduced into the coupling calculation, obtaining accurate distributed solution can be time-consuming. Chippa and Sarangi worked on this method to extract lumped parameters, but their perturbation method was still based on steady-state results. Wherein, a same oil film perturbation was used along the entire contact surface, still leading to rough estimations of integrated contact coefficients by directly accumulating each parameter.

In this paper, elastic deformation keeps the original calculating strategy, but the lubricant stiffness and the damping forces decompose from the Reynolds equation. The mutual effect between oil film and elastic deformation under transient EHL can be maintained by transient damping and stiffness distributions. Theoretically, this simplification is more accurate than two preceding methods.

Reynolds equation:

(2) x(h3ηdpdx)=12uhx+12ht

Using boundary condition Px|xe=0, transient stiffness and damping distributions are derived:

(3) {dfkfilmdx=12ηuxh(x)(he+h(x))3dfcfilmdx=12ηd(h(x)+he)dtdx(he+h(x))3
(4) {dkfilmdx=12ηux(he+h(x))33(he+h(x))212ηuxh(x)(he+h(x))6=36ηuxh(x)(he+h(x))4+12ηux(he+h(x))3dcfilmdx=12ηx(he+h(x))3

Under steady states, equation (4) can be transformed into:

(5) {h03ηdPhdx+3h02ηdP0dx=12uxh03ηdPhdx=12x
where, Ph = k film, Ph = c film

It agrees perfectly with Chippa’s results based on the perturbation method (Nonato and Cavalca, 2010) (see Appendix).

The perturbation method performs well for one parameter in dynamic equations, but not for multiple (or distributed) parameters simultaneously perturbed, using a same condition. No comparability exists among different variables.

More hidden information of the perturbed Reynolds equation is discovered based on the following two equations, where oh2) and oh2) are neglected:

(6) x((h0+Δh3)η(dP0dx+dPhdxΔh+dPhdxΔḣ+PhdΔhdx+PḣdΔḣdx)=12u(h0+Δh)x+12(h0+Δh)t)
(7) x(h03η(dP0dx+dPhdxΔh+dPhdxΔḣ)+3h02Δhη(dP0dx))+x(h03η(PhdΔhdx+PḣdΔḣdx))=12u(h0+Δh)x+12(h0+Δh)t

Chippa assumes that dΔhdx=dΔḣdx=0 exists under transient condition, so that total stiffness and damping coefficients can be derived by {K=PhdxC=Pḣdx:

(8) x ( h 0 3 η ( d P 0 d x + d P h d x Δ h + d P h ˙ d x h ˙ ) + 3 h 0 2 Δ h η ( d P 0 d x ) ) = 12 u ( h 0 x + Δ h x ) + 12 h ˙

To boil down to equation (8) and to respect differences between original and simplified equations, comparisons among internal terms in the original equation (6) tell the approximating conditions.

dPhdxΔhPhdΔhdx=dlnPhdxdlnΔhdx>1dlnPhdxdlnΔhdx=lnPhi+1PhilnΔhi+1ΔhiΔx=lnPhi+1ΔhiPhiΔhi+1Δx>0Phi+1ΔhiPhiΔhi+1>1

It is found that only Phi+1>Phi from the inlet to exit of the EHL contact can get dΔhdx=dΔḣdx=0hi = Δhi + 1).

Further to this, dynamic distribution of oil film calculated by steady-state distributed stiffness and damping can clarify its non-robustness.

Steady-state equation and fluctuation equation of transient Δh can be collected based on equation (7):

(9) {x(h03ηdP0dx)=12uh0xx((h03ηdPhdx+3h02η(dP0dx))Δh)+x(h03ηdPḣdxḣ)+x(h03η(PhdΔhdx+Pḣdḣdx))=12uΔhx+12ḣ

Phk and Pḣk come from the steady-state EHL solution, and serve the fluctuation equation. To test the usability of Phk and Pḣk, the fluctuation between the continuous steady-state EHL solutions should satisfy the load conservation equation (11).

Transient EHL curves will be generated by calculating Δh. If transient EHL curves show reasonable transitions, Phk and Pḣk from each steady-state EHL solution may be useful.

The oil film force evolves:

(10) Pk=P0k+PhkΔh+Pḣkḣ
where, ḣ=h0tk+Δhtkht1kΔt

By keeping the surface loading unchanged, we obtain:

(11) (PhkΔh+Pḣkḣ)=0

With equations (5), the fluctuation equation (9) is reduced and architectured by perturbed Δh, variable ḣ and some items incorporating steady-state solution (marked with A, B, C, D):

(12) 12xdḣdx+x(h03η(PhdΔhdx+Pḣdḣdx))=12xdḣdx+dΔhdxx(h03ηPh)+h03ηPhd2Δhdx2+dḣdxx(h03ηPḣ)+h03ηPḣd2ḣdx2=(x(h03ηPḣ)+12xA)dḣdx+h03ηPḣBd2ḣdx2+x(h03ηPh)CdΔhdx+h03ηPhDd2Δhdx2=0

Distributed perturbed Δh and variable are substituted by discrete forms:

(13) {dΔhdx=(ΔhtkΔhtk1)Δxkd2Δhdx2=(ΔhtkΔhtk1)(Δhtk1Δhtk2)(Δxk)2dḣdx=dh0tk+Δhtkht1kΔtdx=(h0tkh0tk1)+(ΔhtkΔhtk1)(ht1kht1k1)ΔxkΔtd2ḣdx2=(h0tk2h0tk1+h0tk2)+(Δhtk2Δhtk−1+Δhtk2)(ht1k2ht1k1+ht1k2)(Δxk)2Δt
where, superscript k denotes the film thickness at location xk and subscript t film thickness at time t.

Substituting equation (13) into equation (12):

(14) A(h0tkh0tk1)Δhtk1(ht1kht1k1)Δt+B(h0tk2h0tk1+h0tk2)+(2Δhtk1+Δhtk2)(ht1k2ht1k1+ht1k2)ΔxkΔtCΔhtk1+D2Δhtk1+Δhtk2Δxk=AΔhtkΔtBΔhtkΔxkΔtCΔhtkDΔhtkΔxk(A(h0tkh0tk1)Δt+B(h0tk2h0tk1+h0tk2)ΔxkΔt)E(A(ht1kht1k1)Δt+B(ht1k2ht1k1+ht1k2)ΔxkΔt)F+(AΔt+2BΔxkΔtC+2DΔxk)GΔhtk1+(BΔxkΔt+DΔxk)HΔhtk2=Δhtk(AΔt+BΔxkΔt+C+DΔxk)IΔhtk=EF+GΔhtk1+HΔhtk2I

Assume that the first two increments of the film thickness are zero at each state. The Dimensional results under both steady and transient states are shown in Figure 1.

The figures briefly give the perturbed values, to find that steady-state Phk and Pḣk distributions cannot be used under transient EHL conditions, so that real transient values cannot be perturbed with respect to steady-state solutions.

Under transient states, oil film changes Δh and Δḣ lead to new Ph and P and in turn affect integrated stiffness and damping coefficients. To upgrade the robustness, an extended usage of distributed parameters is performed below.

3. Solution of transient EHL contact

Transient EHL procedure starts from steady-state solution. The Newton–Raphson method has been explored by Houpert to solve the steady-state EHL. To take advantage of the proposed contact stiffness distribution, a new steady-state EHL algorithm is performed instead.

The film thickness equation is:

{h(x)=h+x22R+δkδk=4πEbbp(s)ln|sxsx0|ds

Based on the deformation-matrix 6method:

{δk=2πEx1x2p(s)ln(xks)2ds=2πEi=1npix1x2φi(s)ln(xks)2ds=i=1npiDikDik=12π(XkXk1)ln(XkXi)2(ik)

The film thickness equation turns into:

(15) h(x)=h+x22R+i=1npiDik

The Newton iterative method and the following frame constitute the solver (Figure 2):

Ah0=h+x22R+i=1npiDikBh03ηdPhdx+3h02ηdP0dx=12uxCPk=P0k+PhkD2bHertzN(PhkΔh)=0

The same advantage of the proposed contact stiffness and the damping distributions in solving transient EHL problem is also discovered.

The contact stiffness and the damping distributions [equation (4)] will fluctuate under the transient states, changing the film thicknesses along the contact line. The successive states will be obtained by the continuation method.

Assume that pressure increment along the contact line relative to the previous moment is Δp distribution. Thus, update film thickness and its increment:

(16) {ht(x)=h+x22R+i=1n(pi+Δpi)DikΔh(x)=ht(x)ht1(x)

Real-time stiffness and damping distributions k(x), c(x) [equations(4)] establish:

(17) Δh(x)k(x)+Δḣ(x)c(x)=Δp(x)

Combining equation (16) with equation (17):

(18) [h+x22R+i=1n(pi+Δpi)Dikht1(x)](kt1(x)+1Δtct1(x))=Δpi

Besides equation (18), conservation of the load increment should be also satisfied:

(19) 2bHertzNΔpi=Q(t+Δt)Q(t)

Combining equation (18) with equation (19), the increments of the film thickness and contact pressure under transient conditions can be obtained within a few steps by using the Newton iterative method.

Please note, in the dynamic process:

  • real-time update of contact grids; and

  • zero distributed stiffness under pure squeezing condition, but non-zero under both rolling and squeezing conditions

The Newton iterative method and the following frame constitute the solver:

E h t ( x ) = h * + x 2 2 R + i = 1 n ( p i + Δ p i ) D i k F Δ h ( x ) = h t ( x ) h t 1 ( x ) G Δ p i = P h k Δ h + P h ˙ Δ h ˙ H 2 b H e r t z N Δ p i = Q ( t + Δ t ) Q ( t )

Using equation (18), equations E, H, G can turn to I, where D = [Dik] (Figure 3):

h+x22R+D{pi}{ht1}{hc*}+D{Δpi}={Δpikfilm+cfilmΔtCK(i,i)}
I{Δpi}=inv[I(N)2δHertzPHertz×CIJ=CK×D]×(CK=diag(CK(i,i)×{hc*})CIJ(i,:)=D(i,:)×CK(i,i)

Using the proposed semi-analytical method, the pure squeezing effect is calculated below, and its behaviors are studied after suddenly halting the rotating speed of U.

Case 1: U = 14.205 m/s, R = 0.5 m, E0 = 2.273e + 11, W = (1.0e-5 + 1.0e-5 · t) E0 · R N/m, η0 = 0.08

Dimensional results (Figure 4).

During rapid acceleration and deceleration, film behaviors are governed by two different mechanisms – fluid entrainment and film squeeze.

When load, speed and others fluctuate, the transient squeeze effect occurs inevitably. Figure 5 captures what happens when the squeeze effect takes place during loading, and the immediate disappearance of rolling speeds enables the pure squeeze effect to be specifically considered. Zero time is assumed to halt rolling motions.

Case 2: R = 0.5 m, E0 = 2.273e + 11, W = (1.0e-5 + 1.0e-5 · t) ·E0 · R N/m, η0 = 0.08.

Dimensional results.

Figure 5 advocates how the curves change under the pure squeeze action when the initial oil film state is differently selected. These initial states are generated under different entrainment speeds.

The findings are:

  • Two pressure peaks occur near the Hertzian contact center and outlet. With load increasing, the second pressure peak tends to approach the Hertz pressure curve and the mutual influence (lubricant and elastic bodies) is fading away at a longer loading time.

  • The transient oil film distributes more thickly compared to the steady-state solutions. With load increasing, both film peak and the second pressure peak slowly move toward the contact center.

  • The overall effect of the entrainment speed (exits only at zero moment) on the transient EHL resembles that on steady-state EHL. Distributed damping and stiffness of oil film are excited differently under distinct rolling speeds.

Considering engineering facts, only the squeezing action remains when the contact media stop rolling. The method herein cannot work for a long simulation time since the squeezing effect is fading away with time. Fortunately, some traces are found at the early stage of this simulation that Hertz pressure distribution will replace it eventually.

As for the oil film growth at the contact center, lubricants are explained to fail to escape immediately owing to the sudden stop of rolling motion. Practically, the motion change needs time. The pure squeezing action rarely exists if it evolves from normal lubricating operation. However, according to the transient phenomenon of the oil film tested under the ball longitudinal oscillation (Figure 6), Glovnea and Spikes (2005) mentioned that during rapid deceleration, the formation of fluid entrapment dominated the film thickness over the contact; meanwhile, the depth and shape of this entrapment are dependent on the rate of speed variation. The zero halting time assumed herein is an extreme situation, but also tells a same significant truth of existence of fluid entrapment.

In most cases, both entrainment and squeezing effects participate, even in the machine halting. The corresponding case is calculated below:

R = 0.5m, E0 = 2.273e + 11, W = (1.0e-5 + 1.0e-5·t) E0·R N/m, η0 = 0.08

Dimensional results (Figure 7).

The figure shows:

  • Two pressure peaks occur near the inlet and outlet. The pressure near the inlet point distributes much steeper than steady-state solution. Besides, the more load increases, the steeper is the increase in inlet pressure.

  • The transient oil film distributes less thickly with an increase in load. The oil films also shrink at the inlet point compared to the steady-state solution.

  • With an increase in the entrainment velocity, the pressure tends to distribute close to the contact center. The oil film thickness increases and takes more effect.

The oil film transition between successive states is simulated with oil film stiffness and damping fluctuating at the same contact area. But pressure and thickness distributions involved should be updated in time to reproduce each new transient state.

Please remember, it is still a dummy transition. Different from the simulation on the pure squeeze effect; in this case, the growth in the oil inlet can be well predicted, but the recovery in the oil exit cannot. Owing to rolling, the oil film stiffness and damping distributed near the oil outlet should be partly abandoned. If cavitation can be well located during rolling, the force jump at the oil outlet will naturally disappear in simulation.

4. Conclusion

The use of the perturbation method is extended to prove that steady-state EHL contact damping and stiffness distributions should not be used to analyses the transient conditions, let alone acquire equivalent total contact damping and stiffness coefficients for dynamic analysis.

By regarding the steady-state solution of EHL as an initial state, and continuously acquiring EHL stiffness and damping distributions, the dynamic mechanical behaviors of the oil film are predicted. The method provides a possible insight to understand the distributions of the EHL oil film damping and stiffness, and it can be extended for the transient point EHL problem. The idea proposed can be easily delivered to any engineer.

Figures

Transient oil film thickness distributions generated by steady-state solution

Figure 1

Transient oil film thickness distributions generated by steady-state solution

Reduced algorithm for steady-state EHL solution

Figure 2

Reduced algorithm for steady-state EHL solution

Reduced algorithm for transient EHL solution

Figure 3

Reduced algorithm for transient EHL solution

Comparison of steady-state and transient EHL contact pressures and oil film thicknesses

Figure 4

Comparison of steady-state and transient EHL contact pressures and oil film thicknesses

Transient EHL solutions under different entrainment velocities

Figure 5

Transient EHL solutions under different entrainment velocities

Film thickness profile during rapid deceleration in longitudinal oscillatory motion

Figure 6

Film thickness profile during rapid deceleration in longitudinal oscillatory motion

Transient EHL solutions under different entrainment velocities

Figure 7

Transient EHL solutions under different entrainment velocities

Appendix. Perturbation formulation for stiffness and damping distribution

In Reynolds equation, the main variables are perturbed with respect to its steady-state values, so the perturbed film thickness terms and corresponding pressure changes are given by:

(A1) { h = h 0 + Δ h P = P 0 + P h Δ h + P h ˙ Δ h ˙
where:
(A2) ph=Ph,Pḣ=Pḣ

Substituting perturbed terms (A.1) into Reynolds equation:

(A3) x((h0+Δh)3ηd(P0+PhΔh+PḣΔḣ)dx)=12u(h0+Δh)x+12(h0+Δh)t
where, terms containing Δḣn Δhm,(n + m > 1) are neglected.

According to what is referred in equation (19), three equations are acquired:

(A4) {x(h03ηdP0dx)=12uh0xfor(Δh,Δḣ)0x(h03ηdPhdx+3h02ηdP0dx)=0for(Δh)1x(h03ηdPhdx)12=0for(Δḣ)1

References

Chang, L. (2000), “A simple accurate method to calculate transient EHL film thickness in machine components undergoing operation cycles”, Tribology Transactions, Vol. 43 No. 1, pp. 116-222.

Chippa, S.P. and Sarangi, M. (2013), “On the dynamics of lubricated cylindrical roller bearings, part I: evaluation of stiffness and damping characteristics; part II: linear and nonlinear vibration analysis”, Tribology Transactions, Vol. 56 No. 6, pp. 1087-1108.

Dietl, P. (1997), “Damping and stiffness characteristics of rolling element bearings – theory and experiment”, PhD thesis, Technical University of Vienna, Austria.

Dowson, D. and Higginson, G.R. (1966), Elasto-Hydrodynamic Lubrication, Pergamon Press Ltd, SI, Oxford, p. 236.

Evans, H.P. and Snidle, R.W. (1981), “Inverse solution of Reynolds equation of lubrication under point-contact elasto-hydrodynamic conditions”, Journal of Lubrication Technology, Vol. 103 No. 4, pp. 539-546.

Glovnea, R.P. and Spikes, H.A. (2005), “EHD contacts in low-amplitude oscillatory motion”, Proceedings of World Tribology Congress, Washington, DC, 12-16 September .

Grubin, A.N. (1949), Fundamentals of the Hydrodynamic Theory of Lubrication of Heavily Loaded Cylindrical Surfaces, Ketova, K.F. (Ed.), Central Scientific Research Institute for Technology and Mechanical Engineering, 30. Moscow, D.S.I.R. London Translations No. 337, pp. 115-166.

Hagiu, G.D. and Gafitanu, M.D. (1997), “Dynamic characteristics of high speed angular contact ball bearings”, Wear, Vol. 211 No. 1, pp. 22-29.

Hamrock, B.J. and Dowson, D. (1976), “Isothermal EHD lubrication of points contacts: Part I—theoretical formulation”, Journal of Lubrication Technology, Vol. 98 No. 2, pp. 223-229.

Houpert, L.G. and Hamrock, B.J. (1986), “Fast approach for calculating film thicknesses and pressures in elasto-hydrodynamically lubricated contacts at high loads”, Journal of Tribology, Vol. 108 No. 3, pp. 411-420.

Lubrecht, A.A., Ten Napel, W.E. and Bosma, R. (1986), “Multigrid: an alternative method for calculating film thickness and pressure profiles in elasto-hydrodynamically lubricated line contacts”, ASME Journal of Tribology, Vol. 108 No. 4, pp. 551-556.

Messe, S. and Lubrecht, A.A. (2002), “Approximating EHL film thickness profiles under transient conditions”, Journal of Tribology, Vol. 124 No. 3, pp. 443-447.

Nonato, F. and Cavalca, K. (2010), “On the non-linear dynamic behaviour of elasto-hydrodynamic lubricated point contact”, Journal of Sound and Vibration, Vol. 329 No. 22, pp. 4656-4671.

Oh, K.P. and Rohde, S.M. (1977), “Numerical solution of the point contact problem using the finite element method”, International Journal for Numerical Methods in Engineering, Vol. 11 No. 10, pp. 1507-1518.

Rohde, S.M. and Oh, K.P. (1975), “A unified treatment of thick and thin film elasto-hydrodynamic problems by using higher order element methods”, Proceedings of the Royal Society A: Mathematical, Physical and Engineering Sciences, Vol. 343 No. 1634, pp. 315-331.

Sarangi, M., Majumdar, B.C. and Sekhar, A.S. (2004), “Stiffness and damping characteristics of lubricated ball bearings. Part 1: theoretical formulation”, Proceedings of the Institution of Mechanical Engineers, Part J: Journal of Engineering Tribology, Vol. 218 No. 6, pp. 529-538.

Venner, C.H. (1991), “Multilevel solution of the EHL line and point contact problems”, PhD thesis, University of Twente, Endschende.

Walford, T.L.H. and Stone, B.J. (1983), “The sources of damping in rolling element bearings under oscillating conditions”, Proceedings of the Institution of Mechanical Engineers, Part C: Journal of Mechanical Engineering Science, Vol. 197 No. 4, pp. 225-232.

Wang, Z.J., Wang, W.Z., Wang, H. and Hu, Y.Z. (2009), “Stress analysis on layered materials in point elastohydrodynamic-lubricated contacts”, Tribology Letters, Vol. 35 No. 3, pp. 229-244.

Wensing, J. (1998), On the Dynamics of Ball Bearings, PhD thesis, University of Twente, Enschede.

Wiegert, B., Hetzler, H. and Seemann, W. (2013a), “A simplified elasto-hydrodynamic contact model capturing the nonlinear vibration behaviour”, Tribology International, Vol. 59, pp. 79-89.

Wiegert, B., Hetzler, H. and Seemann, W. (2013b), “An analytical expression of the normal force of hydrodynamic line contacts under transient conditions”, Tribology International, Vol. 61, pp. 32-39.

Acknowledgements

National Key R&D Plan(Security techniques on High speed train-analysis on train security and accessment on comprehensive safety, 2016YFB 1200401-102A), China Railway(PHM analysis on D-Series High-Speed Train, 2016J007-B), NSFC (Vibration behaviour, vibration failure and evaluation research on high speed train, U1234208), CRRC Changchun Railway Vehicles (Key parts vibration intensity of CRH380 (BL)), TPL-SWJTU (Axle-box bearing’s matching research with high speed train and its exploitation,TPL1506), and China Scholarship Council.

Corresponding author

Shaocheng Zhu can be contacted at: 506379989@qq.com