# Numerical Investigation into the Effect of Geometric Gap Idealisation on Wheel-Rail Rolling Contact in Presence of Yaw Angle

Authors: Boyang An; Jing Wen; Panjie Wang; Ping Wang; Rong Chen; Jingmang Xu
Journal: Mathematical Problems in Engineering
 (2019)

---

## Abstract

For a fast calculation of vehicle-track dynamics and wheel-rail contact mechanics, wheel-rail contact geometric gap is usually idealised in elliptic or nonelliptic form. These two idealisations deviate from the actual one if the lateral combined curvature within the contact patch is not constant or the yaw angle of wheelset exists. The influence of these idealisations on contact solution has not yet been deeply understood, and thus the accuracy of simplified contact modelling applied to vehicle-track dynamics and wheel-rail contact mechanics remains uncertain. This paper presents a numerical methodology to treat 3D wheel-rail rolling contact, in which the asymmetric geometric gap due to yaw angle is fully taken into account. The attention of this work is placed on investigating the effect of geometric gap idealisation on wheel-rail contact force, rolling contact solution, and wear distribution. It can help with the effective wheel-rail contact modelling on the computation of both vehicle-track dynamics and wheel-rail contact mechanics.


---

## Body

## 1. Introduction

As an economic and environmental protection transportation tool, railway keeps attracting attention from governments and industries. The main parts of railway are vehicle and track, which are coupled by the frictional rolling contact of wheel-rail system. That means the load of dozens of tons is bore by contact patch in the level of mm2, resulting in a very high contact stress. For this reason, wheel-rail materials are easily subjected to rolling contact fatigue or severe wear despite improper wheel-rail profile design or being not timely maintained. The set-up of profile design or maintenance strategy mainly depends on accurate wheel-rail contact modelling and vehicle-track dynamics simulation. In addition, the accurate contact model should be simplified appropriately to be integrated with vehicle-track dynamics to meet its low computational cost demand.In railway, wheel-rail profiles are generally defined by a series of straight lines and arc curves in order to facilitate curving passing and reduce wheel-rail wear. Such a design inevitably results in a periodic wheelset motion due to variable rolling radius difference. During this process, the wheel has two main degrees of freedom: lateral shift and yaw rotation. It may result in three kinds of contact shapes. The first is elliptic shape, symmetric with respect to lateral and longitudinal axes. It occurs only when both the lateral and longitudinal combined curvature of wheel-rail profile keep constant within the contact zone. If the lateral combined curvature is variable (a common phenomenon for wheel-rail contact), the second case of nonelliptic shape will appear, which is unsymmetrical with respect to longitudinal axis. The last one is asymmetric with both lateral and longitudinal axes due to the existence of wheelset’s yaw angle. The last two conditions are also termed non-Hertz geometry. For computational efficiency, the first geometric gap is adopted by most vehicle-track dynamic models to get a global force [1–4], and the second one is usually used for a detailed wheel-rail rolling contact solution even if the actual contact shape is asymmetric [5–7]. To the best of authors’ knowledge, a thorough study of wheel-rail rolling contact considering asymmetric geometry is rare. Therefore, the effect of simplification of contact models by ignoring yaw angle on rolling contact is unclear, which motivates us to finish this paper.Combination of Hertz theory [8] and ellipse shape based tangential method such as Kalker’s FASTSIM algorithm [9] or Shen-Hedrick-Elkins model [10] are widely used in vehicle-track dynamics for fast calculation [1–4]. In this case, one basic assumption is the two contacting bodies are represented by second degree polynomials, meaning contact shape is elliptic. This assumption is invalid if abrupt variation of the radii of curvature of the wheel and rail profile appears near the contact point. So the solution of nonelliptic contact attracts many scholars’ interest, yielding several famous contact models such as those developed by Knothe [11], Kik-Piotrowski [12], Ayasse-Chollet [13], and Sichani [14]. In their models, FASTSIM is adapted for nonelliptic contact in terms of tangential treatment. Liu [15] extended Kik-Piotrowski’s approximate method to consider the effect of yaw angle, realizing asymmetric wheel-rail normal contact in fast modelling. However, the tangential solution may be a challenge for their models from a fast modelling viewpoint because there is no an available approach for FASTSIM applied to asymmetric contact cases.The exact solution for asymmetric contact was finished by Kalker [16] based on half-space assumption and implemented in his well-known Contact code. In his model, the potential contact area is discretized into rectangular elements, and it achieves final solution iteratively till all contact boundary conditions satisfied. Kalker’s model was further extended by Li [17] and Vollebregt [18] for conformal contact problem. Another solution strategy is FE method as reported in [19], which can consider arbitrary contact geometry and elastoplastic material. More advanced FE model is developed by Zhao and Li [20] and further applied for computing nonsteady rolling contact due to yaw angle [21]. Similar FE models are also developed by Pletz [22], An [23], Zong and Dhanasekar [24]. However, high time consuming restricts this kind of FE model to be applied for engineering application, which usually needs more than ten hours to finish a case.Before performing the wheel-rail contact computation, the first step is to determine the contact point, on which the contact patch is acting. In the literature, several methods [25–28] are available for the determination of contact point based on the minimum distance principle. Then a key step is to search the spatial geometric gap between the wheel and rail surface near the contact point. In [11–14], lateral geometric gap is fully considered, while the gap in longitudinal direction is assumed to follow semielliptic distribution. That is to say their geometric gap searching method is two-dimensional. It seems few detailed mathematical descriptions for spatial geometric gap are reported in the literature. To fill this gap, this paper presents a spatial geometric gap searching method for 3D wheel-rail rolling contact solution with yaw angle. Furthermore, it is aimed at investigating the role of spatial geometry on vehicle-track dynamics and wheel-rail rolling contact and to which level the 3D geometric gap should be simplified.


## 2. Methodology

### 2.1. Contact Point Detection Method

When the wheel interacts with the rail, contact geometric parameters mainly include roll angle, rolling contact radius, and wheel-rail gap, which are the function of lateral shift and yaw angle. In order to find the location of wheel-rail contact point, it establishes three Cartesian coordinate systems in Figure1. The first one,o1x1y1z1, is the track coordinate system, witho1x1 axis along the rolling direction,o1y1 axis in the direction of the track transverse, ando1z1 axis perpendicular to the track plane. The second one iso2x2y2z2, parallel too1x1y1z1 ando2 coinciding with the gravity of wheelset. This coordinate system is shifted by lateral shift Δy and lifted by double nominal rolling radius 2R0. The last one,o3x3y3z3, is the wheelset coordinate system formed by rotatingo2x2y2z2 through the shared origin pointo3. In Figure 1,ψandφ are yaw angle and roll angle, respectively. They are defined positive when the wheelset rotates theo2z2 axis ando2x2 axis in anticlockwise. Therefore, these two angles in Figure 1 are negative and positive, respectively.C0 is the lowest point of local circle plane marked in red andOC0 is parallel too3z3 axis. The pointC is where the wheel and rail come into contact, and thus the angleθ is defined as the inclined angle ofOC andOC0.Figure 1
Schematic diagram of contact trace method for the left wheel-rail pair.From Figure1, it can be found searching wheel-rail contact point with yaw angle is a three-dimensional problem, which can be time consuming since the wheel surface near the potential contact point should be meshed. In order to accelerate the searching process, Wang [25] developed a semianalytical method. The idea of this method is to transfer the potential 3D wheel surface to a 2D contact trace by analytical derivation and finally determine the contact point by iteration. An example is shown in Figure 2, which shows the contact trace in wheel under the condition of no yaw angle and -50 mrad. The coordinate systems of wheel and rail areo3x3y3z3 ando1x1y1z1. The blue and red lines mean the potential trace where the contact point will appear.Figure 2
Contact trace in wheel tread (red line: -50 mrad; blue line: no yaw angle).For solving the distance between wheel and rail surface, the local contact trace of wheel ino3x3y3z3 should be converted too1x1y1z1 using (1)x1=y3plx-lxRy3tan⁡δy3y1=y3ply+Ry31-lx2lx2lytan⁡δy3+lz1-lx21+tan2δy3+Δyz1=y3plz+Ry31-lx2lx2lztan⁡δy3-ly1-lx21+tan2δy3+2R0where y3p is the lateral coordinate of wheel profile ino3x3y3z3 and δ(y3) is the corresponding contact angle.lx,ly, andlz are direction vectors,(2)lx=-cos⁡φ·sin⁡ψly=cos⁡φ·cos⁡ψlz=sin⁡φOnce having determined the wheel’s contact trace, the corresponding vertical coordinatez1r in rail surface can be easily obtained by interpolating the rail profile withy1 since the rail section is constant in the longitudinal direction.(3)z1r=Iry1whereIr means the interpolation function of rail profile. The contact point can be found through adjusting roll angleφ till both sides’ minimum distance between the wheel and rail in vertical direction getting equality. According to the experience of authors, it needs 5~10 iterations to finish a case.


### 2.2. Determination of Wheel-Rail Contact Gap in 3D

#### 2.2.1. Elliptic Gap

In vehicle system dynamics simulation, wheel-rail contact problem is often treated by Hertz theory in the contact point determined in Section2.1, and the geometric gap is simplified from (4)h=AxL2+ByL2where xL and yL are coordinates of local contact system oLxLyLzL that are located on the tangential plane with its origin in the rigid contact point as shown in Figure 3. The longitudinal axis oLxL is parallel too1x1, vertical axis oLzL points upwards the normal direction, and the lateral axis oLyL is defined to form a right-hand coordinate system.A andB are longitudinal and lateral combined curvature of wheel and rail in contact point,(5)A=12Ry3cB=12Cw+CrHere,Cw andCr are the lateral curvature of the wheel and rail; R(y3c) is the rolling radius of the contact point. It should be mentioned the contact point (x1c,y1c,z1c) in wheel surface searched in Section 2.1 is belonging to the track systemo1x1y1z1, which should be converted to (x3c,y3c,z3c) in the wheelset coordinate systemo3x3y3z3 by (6) to find the corresponding curvature.(6)x3cy3cz3c=cos⁡ψsin⁡ψ0-sin⁡ψcos⁡φcos⁡ψcos⁡φsin⁡φsin⁡ψsin⁡φ-cos⁡ψsin⁡φcos⁡φx1cy1c-Δyz1c-2R0Figure 3
The schematic diagram for wheel-rail gap.


#### 2.2.2. Asymmetric Gap

In order to search the actual wheel-rail contact geometric gap, this part describes an approach by meshing wheel surface in 3D. Since the wheel can be considered as the solid of revolution, its vertical coordinatez3 in any wheel surface point (x3,y3) can be obtained:(7)z3x,y=-Ry32-x32Then, the rotated wheel surface for the contact status can be obtained with the prescribed yaw angleψ and calculated roll angleφ,(8)x1y1z1=cos⁡ψ-sin⁡ψcos⁡φsin⁡ψsin⁡φsin⁡ψcos⁡ψcos⁡φ-cos⁡ψsin⁡φ0sin⁡φcos⁡φx3y3z3+0Δy2R0In this approach, it assumes the contact patch is a flat based on half-space assumption as shown in Figure 3. This assumption is also adopted by many contact models including Kalker’s exact theory [16]. By (9), the coordinates of wheel and rail surface can be converted to local contact system oLxLyLzL.(9)xLiyLizLi=cos⁡θsin⁡θsin⁡δ-sin⁡θcos⁡δ0cos⁡δsin⁡δsin⁡θ-cos⁡θsin⁡δcos⁡θcos⁡δΔxiΔyiΔzi,i=w,rwhere (Δxi,Δyi,Δzi) are relative distance of wheel and rail surface to the contact point in three directions.θ is calculated as follows:(10)θ=arcsin⁡x3cz3cUsually,θ is a tiny value and can be ignored. The wheel-rail gap can finally be obtained from (11)h=zLw-zLr


#### 2.2.3. Nonelliptic Gap

The nonelliptic gapf(yL) along the lateral coordinate is a 2D case of asymmetric gap assuming the yaw angleψ is zero. Its modelling efficiency is sped up since wheel surface meshing is not needed. The geometric gap along thex-direction is assumed to follow the elliptic distribution:(12)h=AxL2+fyL


### 2.3. Solving Asymmetric Wheel-Rail Rolling Contact

As reviewed in Section1, two classes of rolling contact models are available for asymmetric contact problem: half-space based [16–18] and FE based [19–24] models. For consideration of accuracy and efficiency, Kalker’s exact theory [16] is employed in this paper to calculate asymmetric wheel-rail rolling contact. In this work, the wheel and rail are assumed to share the same material properties. Hence, the normal and tangential problem can be treated separately.


#### 2.3.1. Normal Contact Modelling

The normal contact problem can be treated as(13)h=δn-ux,ywhereh is the wheel-rail gap solved in Section 2.2;δn is penetration between the wheel and rail;u(x,y) is the elastic deformation calculated in (14)ux,y=21-v2πE∬AcAIzJzpJzx,ydxdywhereAIzJz is the influence coefficient inz-direction, meaning the displacement in elementI subjected to a unit load in elementJ.pJz is the normal pressure of elementJ.v andE are material’s Poisson’s ratio and elastic modulus.From (14), we can see elastic deformation depends on pressure within contact patch, and on the other hand, the determination of contact boundary also needs to consider elastic deformation. Thus there is no explicit solution for (13). The final solution can only be solved by iteration until pressure of any elementI within the contact patch meet the boundary of (15)min⁡ϕ=12pIzAIzJzpJz+h-δnpIzpIz>0,∀I∈AcIn application to vehicle-track dynamics computation, penetrationδn can be calculated in every time-step; see Shabana’s method [26]. However, penetrationδn is an unknown quantity for the case of prescribed load as presented in this paper. An additional loop of iteration is needed to adjust penetrationδn to bear the prescribed load.(16)∑IMpIzδxδy=N


#### 2.3.2. Tangential Contact Modelling

In the tangential contact modelling, the frictional contact is expressed as a variational inequality; thus the problem is transformed into solving minimum of complementary energyϕ, the product of tractions, and displacements.(17)min⁡ϕ=12pIiAIiJjpJj+WIτ-uIτ′pIτwherei,j stands for any direction (x,y orz), andτ = 1, 2, which meansx- andy-direction.pIi stands for the surface traction in elementI in directioni; uIτ′ is the previous displacement at elementI in directionτ;WIz is the rigid shift at elementI and can be calculated by integration from the previous instance t′ to present onet,(18)WIτ=∫t′tcτ+-1τx3-τ′c3vdtwherev is the running velocity of the vehicle; x3-τ′ is localx- ory- coordinates; cτ and c3 are creepage and spin, respectively. They are calculated as follows:(19)c1=1-Ry3cR0cos⁡ψc2=-sin⁡ψcos⁡φ+δy3cc3=-sin⁡δy3cR0


### 2.4. Wear Depth Prediction

The USFD wear model [29] is introduced to evaluate wear depth. It assumes contact patch keeps constant during the process of wheel passing a section of wheel or rail. Therefore, the wear depth distribution at one section can be calculated by integrating the wear contribution in the longitudinal direction:(20)δtoty=∫-ayayδx,ydxwheredx is the element size inx-direction, and it is 0.2 mm in this paper;δ(x,y) is the wear depth at each element, which can be calculated with wear ratek and material density ρ as(21)δx,y=kρdxThe wear ratek is a function of local frictional powerI(x,y) of each element in the slip zone,(22)k=5.3Ix,yI<10.455.010.4≤I≤77.255.0+61.9×Ix,y-77.2I>77.2(23)Ix,y=qtx,y·sx,ywhere qt(x,y) ands(x,y) are tangential tractions and rigid slip, respectively. More detailed explanation of this method is referred to [30].


## 2.1. Contact Point Detection Method

When the wheel interacts with the rail, contact geometric parameters mainly include roll angle, rolling contact radius, and wheel-rail gap, which are the function of lateral shift and yaw angle. In order to find the location of wheel-rail contact point, it establishes three Cartesian coordinate systems in Figure1. The first one,o1x1y1z1, is the track coordinate system, witho1x1 axis along the rolling direction,o1y1 axis in the direction of the track transverse, ando1z1 axis perpendicular to the track plane. The second one iso2x2y2z2, parallel too1x1y1z1 ando2 coinciding with the gravity of wheelset. This coordinate system is shifted by lateral shift Δy and lifted by double nominal rolling radius 2R0. The last one,o3x3y3z3, is the wheelset coordinate system formed by rotatingo2x2y2z2 through the shared origin pointo3. In Figure 1,ψandφ are yaw angle and roll angle, respectively. They are defined positive when the wheelset rotates theo2z2 axis ando2x2 axis in anticlockwise. Therefore, these two angles in Figure 1 are negative and positive, respectively.C0 is the lowest point of local circle plane marked in red andOC0 is parallel too3z3 axis. The pointC is where the wheel and rail come into contact, and thus the angleθ is defined as the inclined angle ofOC andOC0.Figure 1
Schematic diagram of contact trace method for the left wheel-rail pair.From Figure1, it can be found searching wheel-rail contact point with yaw angle is a three-dimensional problem, which can be time consuming since the wheel surface near the potential contact point should be meshed. In order to accelerate the searching process, Wang [25] developed a semianalytical method. The idea of this method is to transfer the potential 3D wheel surface to a 2D contact trace by analytical derivation and finally determine the contact point by iteration. An example is shown in Figure 2, which shows the contact trace in wheel under the condition of no yaw angle and -50 mrad. The coordinate systems of wheel and rail areo3x3y3z3 ando1x1y1z1. The blue and red lines mean the potential trace where the contact point will appear.Figure 2
Contact trace in wheel tread (red line: -50 mrad; blue line: no yaw angle).For solving the distance between wheel and rail surface, the local contact trace of wheel ino3x3y3z3 should be converted too1x1y1z1 using (1)x1=y3plx-lxRy3tan⁡δy3y1=y3ply+Ry31-lx2lx2lytan⁡δy3+lz1-lx21+tan2δy3+Δyz1=y3plz+Ry31-lx2lx2lztan⁡δy3-ly1-lx21+tan2δy3+2R0where y3p is the lateral coordinate of wheel profile ino3x3y3z3 and δ(y3) is the corresponding contact angle.lx,ly, andlz are direction vectors,(2)lx=-cos⁡φ·sin⁡ψly=cos⁡φ·cos⁡ψlz=sin⁡φOnce having determined the wheel’s contact trace, the corresponding vertical coordinatez1r in rail surface can be easily obtained by interpolating the rail profile withy1 since the rail section is constant in the longitudinal direction.(3)z1r=Iry1whereIr means the interpolation function of rail profile. The contact point can be found through adjusting roll angleφ till both sides’ minimum distance between the wheel and rail in vertical direction getting equality. According to the experience of authors, it needs 5~10 iterations to finish a case.


## 2.2. Determination of Wheel-Rail Contact Gap in 3D

### 2.2.1. Elliptic Gap

In vehicle system dynamics simulation, wheel-rail contact problem is often treated by Hertz theory in the contact point determined in Section2.1, and the geometric gap is simplified from (4)h=AxL2+ByL2where xL and yL are coordinates of local contact system oLxLyLzL that are located on the tangential plane with its origin in the rigid contact point as shown in Figure 3. The longitudinal axis oLxL is parallel too1x1, vertical axis oLzL points upwards the normal direction, and the lateral axis oLyL is defined to form a right-hand coordinate system.A andB are longitudinal and lateral combined curvature of wheel and rail in contact point,(5)A=12Ry3cB=12Cw+CrHere,Cw andCr are the lateral curvature of the wheel and rail; R(y3c) is the rolling radius of the contact point. It should be mentioned the contact point (x1c,y1c,z1c) in wheel surface searched in Section 2.1 is belonging to the track systemo1x1y1z1, which should be converted to (x3c,y3c,z3c) in the wheelset coordinate systemo3x3y3z3 by (6) to find the corresponding curvature.(6)x3cy3cz3c=cos⁡ψsin⁡ψ0-sin⁡ψcos⁡φcos⁡ψcos⁡φsin⁡φsin⁡ψsin⁡φ-cos⁡ψsin⁡φcos⁡φx1cy1c-Δyz1c-2R0Figure 3
The schematic diagram for wheel-rail gap.


### 2.2.2. Asymmetric Gap

In order to search the actual wheel-rail contact geometric gap, this part describes an approach by meshing wheel surface in 3D. Since the wheel can be considered as the solid of revolution, its vertical coordinatez3 in any wheel surface point (x3,y3) can be obtained:(7)z3x,y=-Ry32-x32Then, the rotated wheel surface for the contact status can be obtained with the prescribed yaw angleψ and calculated roll angleφ,(8)x1y1z1=cos⁡ψ-sin⁡ψcos⁡φsin⁡ψsin⁡φsin⁡ψcos⁡ψcos⁡φ-cos⁡ψsin⁡φ0sin⁡φcos⁡φx3y3z3+0Δy2R0In this approach, it assumes the contact patch is a flat based on half-space assumption as shown in Figure 3. This assumption is also adopted by many contact models including Kalker’s exact theory [16]. By (9), the coordinates of wheel and rail surface can be converted to local contact system oLxLyLzL.(9)xLiyLizLi=cos⁡θsin⁡θsin⁡δ-sin⁡θcos⁡δ0cos⁡δsin⁡δsin⁡θ-cos⁡θsin⁡δcos⁡θcos⁡δΔxiΔyiΔzi,i=w,rwhere (Δxi,Δyi,Δzi) are relative distance of wheel and rail surface to the contact point in three directions.θ is calculated as follows:(10)θ=arcsin⁡x3cz3cUsually,θ is a tiny value and can be ignored. The wheel-rail gap can finally be obtained from (11)h=zLw-zLr


### 2.2.3. Nonelliptic Gap

The nonelliptic gapf(yL) along the lateral coordinate is a 2D case of asymmetric gap assuming the yaw angleψ is zero. Its modelling efficiency is sped up since wheel surface meshing is not needed. The geometric gap along thex-direction is assumed to follow the elliptic distribution:(12)h=AxL2+fyL


## 2.2.1. Elliptic Gap

In vehicle system dynamics simulation, wheel-rail contact problem is often treated by Hertz theory in the contact point determined in Section2.1, and the geometric gap is simplified from (4)h=AxL2+ByL2where xL and yL are coordinates of local contact system oLxLyLzL that are located on the tangential plane with its origin in the rigid contact point as shown in Figure 3. The longitudinal axis oLxL is parallel too1x1, vertical axis oLzL points upwards the normal direction, and the lateral axis oLyL is defined to form a right-hand coordinate system.A andB are longitudinal and lateral combined curvature of wheel and rail in contact point,(5)A=12Ry3cB=12Cw+CrHere,Cw andCr are the lateral curvature of the wheel and rail; R(y3c) is the rolling radius of the contact point. It should be mentioned the contact point (x1c,y1c,z1c) in wheel surface searched in Section 2.1 is belonging to the track systemo1x1y1z1, which should be converted to (x3c,y3c,z3c) in the wheelset coordinate systemo3x3y3z3 by (6) to find the corresponding curvature.(6)x3cy3cz3c=cos⁡ψsin⁡ψ0-sin⁡ψcos⁡φcos⁡ψcos⁡φsin⁡φsin⁡ψsin⁡φ-cos⁡ψsin⁡φcos⁡φx1cy1c-Δyz1c-2R0Figure 3
The schematic diagram for wheel-rail gap.


## 2.2.2. Asymmetric Gap

In order to search the actual wheel-rail contact geometric gap, this part describes an approach by meshing wheel surface in 3D. Since the wheel can be considered as the solid of revolution, its vertical coordinatez3 in any wheel surface point (x3,y3) can be obtained:(7)z3x,y=-Ry32-x32Then, the rotated wheel surface for the contact status can be obtained with the prescribed yaw angleψ and calculated roll angleφ,(8)x1y1z1=cos⁡ψ-sin⁡ψcos⁡φsin⁡ψsin⁡φsin⁡ψcos⁡ψcos⁡φ-cos⁡ψsin⁡φ0sin⁡φcos⁡φx3y3z3+0Δy2R0In this approach, it assumes the contact patch is a flat based on half-space assumption as shown in Figure 3. This assumption is also adopted by many contact models including Kalker’s exact theory [16]. By (9), the coordinates of wheel and rail surface can be converted to local contact system oLxLyLzL.(9)xLiyLizLi=cos⁡θsin⁡θsin⁡δ-sin⁡θcos⁡δ0cos⁡δsin⁡δsin⁡θ-cos⁡θsin⁡δcos⁡θcos⁡δΔxiΔyiΔzi,i=w,rwhere (Δxi,Δyi,Δzi) are relative distance of wheel and rail surface to the contact point in three directions.θ is calculated as follows:(10)θ=arcsin⁡x3cz3cUsually,θ is a tiny value and can be ignored. The wheel-rail gap can finally be obtained from (11)h=zLw-zLr


## 2.2.3. Nonelliptic Gap

The nonelliptic gapf(yL) along the lateral coordinate is a 2D case of asymmetric gap assuming the yaw angleψ is zero. Its modelling efficiency is sped up since wheel surface meshing is not needed. The geometric gap along thex-direction is assumed to follow the elliptic distribution:(12)h=AxL2+fyL


## 2.3. Solving Asymmetric Wheel-Rail Rolling Contact

As reviewed in Section1, two classes of rolling contact models are available for asymmetric contact problem: half-space based [16–18] and FE based [19–24] models. For consideration of accuracy and efficiency, Kalker’s exact theory [16] is employed in this paper to calculate asymmetric wheel-rail rolling contact. In this work, the wheel and rail are assumed to share the same material properties. Hence, the normal and tangential problem can be treated separately.


### 2.3.1. Normal Contact Modelling

The normal contact problem can be treated as(13)h=δn-ux,ywhereh is the wheel-rail gap solved in Section 2.2;δn is penetration between the wheel and rail;u(x,y) is the elastic deformation calculated in (14)ux,y=21-v2πE∬AcAIzJzpJzx,ydxdywhereAIzJz is the influence coefficient inz-direction, meaning the displacement in elementI subjected to a unit load in elementJ.pJz is the normal pressure of elementJ.v andE are material’s Poisson’s ratio and elastic modulus.From (14), we can see elastic deformation depends on pressure within contact patch, and on the other hand, the determination of contact boundary also needs to consider elastic deformation. Thus there is no explicit solution for (13). The final solution can only be solved by iteration until pressure of any elementI within the contact patch meet the boundary of (15)min⁡ϕ=12pIzAIzJzpJz+h-δnpIzpIz>0,∀I∈AcIn application to vehicle-track dynamics computation, penetrationδn can be calculated in every time-step; see Shabana’s method [26]. However, penetrationδn is an unknown quantity for the case of prescribed load as presented in this paper. An additional loop of iteration is needed to adjust penetrationδn to bear the prescribed load.(16)∑IMpIzδxδy=N


### 2.3.2. Tangential Contact Modelling

In the tangential contact modelling, the frictional contact is expressed as a variational inequality; thus the problem is transformed into solving minimum of complementary energyϕ, the product of tractions, and displacements.(17)min⁡ϕ=12pIiAIiJjpJj+WIτ-uIτ′pIτwherei,j stands for any direction (x,y orz), andτ = 1, 2, which meansx- andy-direction.pIi stands for the surface traction in elementI in directioni; uIτ′ is the previous displacement at elementI in directionτ;WIz is the rigid shift at elementI and can be calculated by integration from the previous instance t′ to present onet,(18)WIτ=∫t′tcτ+-1τx3-τ′c3vdtwherev is the running velocity of the vehicle; x3-τ′ is localx- ory- coordinates; cτ and c3 are creepage and spin, respectively. They are calculated as follows:(19)c1=1-Ry3cR0cos⁡ψc2=-sin⁡ψcos⁡φ+δy3cc3=-sin⁡δy3cR0


## 2.3.1. Normal Contact Modelling

The normal contact problem can be treated as(13)h=δn-ux,ywhereh is the wheel-rail gap solved in Section 2.2;δn is penetration between the wheel and rail;u(x,y) is the elastic deformation calculated in (14)ux,y=21-v2πE∬AcAIzJzpJzx,ydxdywhereAIzJz is the influence coefficient inz-direction, meaning the displacement in elementI subjected to a unit load in elementJ.pJz is the normal pressure of elementJ.v andE are material’s Poisson’s ratio and elastic modulus.From (14), we can see elastic deformation depends on pressure within contact patch, and on the other hand, the determination of contact boundary also needs to consider elastic deformation. Thus there is no explicit solution for (13). The final solution can only be solved by iteration until pressure of any elementI within the contact patch meet the boundary of (15)min⁡ϕ=12pIzAIzJzpJz+h-δnpIzpIz>0,∀I∈AcIn application to vehicle-track dynamics computation, penetrationδn can be calculated in every time-step; see Shabana’s method [26]. However, penetrationδn is an unknown quantity for the case of prescribed load as presented in this paper. An additional loop of iteration is needed to adjust penetrationδn to bear the prescribed load.(16)∑IMpIzδxδy=N


## 2.3.2. Tangential Contact Modelling

In the tangential contact modelling, the frictional contact is expressed as a variational inequality; thus the problem is transformed into solving minimum of complementary energyϕ, the product of tractions, and displacements.(17)min⁡ϕ=12pIiAIiJjpJj+WIτ-uIτ′pIτwherei,j stands for any direction (x,y orz), andτ = 1, 2, which meansx- andy-direction.pIi stands for the surface traction in elementI in directioni; uIτ′ is the previous displacement at elementI in directionτ;WIz is the rigid shift at elementI and can be calculated by integration from the previous instance t′ to present onet,(18)WIτ=∫t′tcτ+-1τx3-τ′c3vdtwherev is the running velocity of the vehicle; x3-τ′ is localx- ory- coordinates; cτ and c3 are creepage and spin, respectively. They are calculated as follows:(19)c1=1-Ry3cR0cos⁡ψc2=-sin⁡ψcos⁡φ+δy3cc3=-sin⁡δy3cR0


## 2.4. Wear Depth Prediction

The USFD wear model [29] is introduced to evaluate wear depth. It assumes contact patch keeps constant during the process of wheel passing a section of wheel or rail. Therefore, the wear depth distribution at one section can be calculated by integrating the wear contribution in the longitudinal direction:(20)δtoty=∫-ayayδx,ydxwheredx is the element size inx-direction, and it is 0.2 mm in this paper;δ(x,y) is the wear depth at each element, which can be calculated with wear ratek and material density ρ as(21)δx,y=kρdxThe wear ratek is a function of local frictional powerI(x,y) of each element in the slip zone,(22)k=5.3Ix,yI<10.455.010.4≤I≤77.255.0+61.9×Ix,y-77.2I>77.2(23)Ix,y=qtx,y·sx,ywhere qt(x,y) ands(x,y) are tangential tractions and rigid slip, respectively. More detailed explanation of this method is referred to [30].


## 3. Results and Discussion

In order to model obvious non-Hertz properties of contact patch, the standard profile S1002CN-CN60 of wheel-rail system with 1:40 inclination is considered. This profile combination is widely used in Chinese high-speed railway. The materials of wheel and rail are assumed to be identical, whose Poisson’s ratio and elastic modulus are 0.28 and 206 GPa, respectively. The friction between the wheel-rail interface is 0.5 for a dry condition.


### 3.1. Normal Contact

This section intends to investigate the effect of geometrical gap assumptions on normal contact modelling. A set of lateral shift cases with the prescribed load of 83.3 kN are conducted. With the increase of lateral displacement, the contact patch shift towards the rail gauge as shown in Figure4, of which the nonelliptic gap is employed.Figure 4
Variation of nonelliptic contact shape and pressure with lateral shift varying from -3 to 7 mm (from left to right) in 1 mm increment.


#### 3.1.1. Geometric Gap, Contact Shape, and Pressure

Typical nonelliptic contact conditions where the wheelset has a lateral displacement of -2, 0, and 6 mm, and a similar elliptic case of Δy = 3 mm is chosen for a detailed comparison. Figure 5 illustrates the comparison of contact shape, geometric gap (atx= 0), and pressure. Two yaw angles, namely, 5 and 25 mrad, are considered to express asymmetric gap to different level. It is mentioned that 25 mrad means a relatively large yaw angle for the actual vehicle operation, which can only appear for the small curving (below 100 m) operational condition.Figure 5
Comparison of contact shape (left), gap (middle), and pressure (right) under different lateral shift.
(a)
Δy = -2 mm (b)
Δy = 0 mm (c)
Δy = 3 mm (d)
Δy = 6 mmFor Δy = -2 mm case in Figure 5(a), the lengths of contact patch in these four geometric gaps are almost the same, and their widths only differ a bit. By contrast, the pressure variations within the contact patch for elliptic shape and the other three differ a lot since an abrupt change of curvature appears aty = -8.5 mm. In the second case, i.e., wheelset central position, the contact shape and pressure predicted by elliptic gap and the other three are violently different because the gaps of the latter ones in positive direction differ a lot from the elliptic one. In this case, the pressure for the elliptic gap is higher due to smaller contact area. In Figure 5(c), contact patch and pressure distribution in these four gaps agree very well since they all follow the elliptic form. For the last case, the contact patch is composed of two merging ellipses and pressure appears two peaks. The result based on elliptic gap assumption fails to predict both contact patch and pressure as expected.It is further found the influence of yaw angle on contact shape and pressure in the first three shift cases can be ignored, though a large yaw angle of 25 mrad has been enforced. By contrast, the effect of this yaw angle for case Δy = 6 mm is significant. This phenomena may be explained as contact points for the first three cases located in the rail top and wheel tread that are relatively ‘flat’, while that in the last case locates near the rail gauge and wheel flange that are ‘sharp’. In the condition of ‘sharp’, the curvature varies sharply and is sensitive to yaw angle like in Figure 5(d).


#### 3.1.2. Curve of Penetration versus Contact Force

Then we evaluate the performance of geometric gap simplification on resultant force, which is very important for vehicle dynamics simulation. Figure6 compares the penetration versus contact force curve for the four shift cases. Since the result of 5 mrad yaw angle is almost the same as that using nonelliptic gap, only the result for 25 mrad is presented for the asymmetric case. For both Δy = -2 and 3 mm, the curves predicted by different gap assumptions coincide very well, although in the former case the pressure distribution are different. In the other two shift cases, there is an obvious difference among elliptic shape and the other two. However, it may not significantly affect their performance on predicting normal force as explained in the following.Figure 6
Penetration versus contact force curve under different lateral shift.
(a)
Δy = -2 mm (b)
Δy = 0 mm (c)
Δy = 3 mm (d)
Δy = 6 mmBefore a vehicle-track dynamics computation, it needs a preprocessing to balance the prescribed load and then the penetrationδi is obtained depending on contact model for instance in Figure 6(d). If the vehicle is running on a steady-state, the contact force and penetration will keep constant. Otherwise, the contact force will vibrate due to irregularity between the wheel-rail interfaces. The dynamic forceFd can be obtained from (24)Fd=fδi+Δδ≈Fs+Δδ·tan⁡αiwhere subscripti = 1~3 meaning the elliptic, nonelliptic, and asymmetric gap;Fs is the prescribed load; Δδ is irregularity size; α is the local angle as shown in Figure 6(d).For the case Δy = 6 mm, α1 is smaller than the other two while α2 and α3 are nearly the same. So the dynamic forces predicted by nonelliptic and asymmetric gaps will be similar and that using elliptic gap may be underestimated.


### 3.2. Tangential Contact

In this section, an elaborated analysis of shear stress, stick-slip division, microslip, and wear distribution will be presented with full consideration of yaw angle. From (19) and Figure 5, we can see yaw angle devotes to two main changes of contact solution: lateral creepage and asymmetric contact. Therefore, both of them are taken into consideration in this part.Two shift cases, Δy= -2 and 6 mm, are chosen for evaluating the effect of yaw angle on tangential modelling. The reasons for selecting these two cases include the following:(1)
In the former case, both longitudinal creepage and spin are low as the contact point occurring in the rail top, thus facilitating showing the role of lateral creepage on tangential contact.(2)
For the latter case, the contact shape can be changed obviously and then its effect on tangential contact can be considered.


#### 3.2.1. Shear Stress, Slip and Wear Distribution

The stick-slip division, shear stress distribution, and its direction for two shift cases are shown in Figure7. Only results of elliptic and asymmetric shape are included since that by nonelliptic shape is similar to the asymmetric one. For each figure, the slip zone is encircled by the red line and the rest part means stick zone. The vector points the direction of shear stress, and its length is proportional to the magnitude of shear stress.Figure 7
Comparison of shear stress distribution and stick-slip division for two shift cases.
(a)
Δy = -2 mm (b)
Δy = 6 mmIn Figure7(a), the contact shape is nearly not influenced by yaw angle within 5 mrad. It can be found the increase of yaw angle make the slip zone gradually dominates. In the meantime, the directions of shear stress change from longitudinal to transversal. Although the contact shapes are different, the change tendencies for these two cases are the same. The difference of location and magnitude of shear stress is due to different pressure distribution as shown in Figure 5(a).In Figure7(b), relative large longitudinal creepage and spin make the slip zone dominating the contact patch. So it can easily reach a full slip condition once yaw angle appears. However, in this case, it needs a larger yaw angle (or lateral creepage) to change the direction of shear stress. Although the contact shapes under these two gap assumptions are far different, the evolution of stick-slip ratio and shear stress’ direction is similar.For the two shift cases, their relative slip distributions along longitudinal (y=0) direction are illustrated in Figure 8. In the case Δy = -2 mm, the magnitude and distribution of relative slip coincide very well for the two contact shape assumptions. It is found that the increasing yaw angle can enlarge relative slip continuously even though it has achieved a full slip condition (see results of 10 mrad). In Figure 8(b), the relative slip can also achieve a good agreement though their contact shapes vary significantly. This is because the relative slip mainly depends on rigid slip as defined in (19) and the difference is initiated from deformation. Such a difference can be reduced if the full slip condition is reached since the deformation caused by saturated shear stress is similar.Figure 8
Relative slip comparison for (a) Δy = -2 mm and (b) Δy = 6 mm.
(a) (b)The wear distributions predicted by different geometrical gap assumptions are compared in Figure9. In these results, lateral creepage keep the same for three gap models. In addition, the running speed of vehicle is 300 km/h. We can find three main phenomena: (1) there is a significant difference in terms of wear distribution and maximum wear depth among elliptic solution and the others as a result of shear stress and slip in Figures 7 and 8; (2) the yaw angle can continuously worse wear depth due to slip increase; (3) the difference between nonelliptic and asymmetric gap assumption is nonnegligible under case Δy = 6 mm and 25 mrad.Figure 9
Wear distributions under the same combination of creepage and spin.
(a)
Δy = -2 mm (yaw angle = 0, 1, 5 mrad from left to right) (b)
Δy = 6 mm (yaw angle = 0, 5, 25 mrad from left to right)


#### 3.2.2. Creep Curve

Three gap simplifications in Figure6 are employed in this part to evaluate their performance on creep curve. Only pure longitudinal creepage or spin is considered ignoring lateral creepage initiated from yaw angle and thus facilitates understanding the effect of asymmetric contact shape on creep force estimation.Figures10(a) and 10(b) illustrate the longitudinal creep force curve due to pure longitudinal creepage for the three contact cases. Here, creep force is normalized by the normal force. As can been seen, the nonelliptic force estimation agrees with the reference for both shift cases, while the elliptic estimates a lower initial slope of the force build-up for the latter shift case.Figure 10
Creep curve estimation for three contact shape under (a), (b) pure longitudinal and (c), (d) pure spin condition.
(a)
Δy = -2 mm (b)
Δy = 6 mm (c)
Δy = -2 mm (d)
Δy = 6 mmIn addition to the pure creepage case, creep curve for the pure spin condition is also investigated. In Figures10(c) and 10(d), the saturated longitudinal and lateral creep force were estimated by both elliptic and nonelliptic gap deviates from the reference. This is attributed to double or single symmetry of the patch adopted separately in these two gap models. In the case of Δy = 6 mm, the estimation error is large due to obvious difference among the contact shapes predicted by three gap assumptions.


### 3.3. Discussion of This Paper

From aforementioned results, we can find yaw angle plays a double effect on wheel-rail rolling contact solution, i.e., lateral creepage and contact shape. It yields a question about suited contact modelling in vehicle-track dynamics and wheel-rail contact mechanics. Note that discussion below is based on the assumption wheel-rail profiles that are smooth. The irregularity case is beyond the applicability of half-space method employed in this paper.


#### 3.3.1. Effective Geometry Gap and Contact Model for Wheel-Rail Contact

For wheel-rail contact modelling, the effect of yaw angle on normal contact such as shape size and pressure is negligible even for a highly nonelliptic contact condition such as in Figure5(b). By contrast, tangential solutions such as stick-slip division and traction direction can be significantly influenced by yaw angle due to the appearance of lateral creepage; see Figure 7. If the true lateral creepage is taken into account; solution using the nonelliptic geometry gap assumption can coincide well with that in the asymmetric contact case (the true solution). This is verified from the comparison of wear performance (see Figure 9) that consider the combined effect of yaw angle on contact shape and lateral creepage. So in most wheel-rail contact mechanics cases, mainly including wheel tread-rail top contact, nonelliptic geometry gap simplification is sufficient. For the severe contact condition, an exact asymmetric geometry gap is needed. The term ‘severe contact’ here mainly means contact locates on flange-rail gauge contact where the contact shape is not an ellipse or similar like shift case Δy = 6 mm with yaw angle of 25 mrad.For modelling approach, Kalker’s Contact is used in this paper. It is, however, generally considered improper for online vehicle-track dynamics modelling as well as wear prediction for relatively high computational cost. As an alternative solution, simplified nonelliptic contact models [12–14] are regarded as fast and suitable for general nonelliptic cases, while it is considered not sufficiently exact for highly nonelliptic case [31]. The best choice may be a mixed contact models strategy for different contact cases, that is, (1) Hertz theory [8] and FASTSIM [9] used in similar elliptic case for instance Δy = 3 mm; (2) simplified nonelliptic models like KP [12], STRIPES [13] and ANALYN+FaStrip [14] used to treat general nonelliptic cases such as Δy = -2 mm; and (3) employ Contact to solve highly nonelliptic or asymmetric case (see Δy = 0 and 6 mm). To apply this strategy, a calibration work of nonelliptic level should be performed. A simple measure may be the ratio of left and right contact boundary determined by a scaling penetration similar to that in KP model [12], identifying to use which contact model. This work is within our scope of further study.


#### 3.3.2. Effective Geometry Gap for Vehicle-Track Interaction

According to Figures6 and 10, the contact solution difference reflected on contact force is however acceptable for elliptic and nonelliptic gap simplification as explained in Section 3.1.2.For a higher accuracy, nonelliptic gap should be taken into account especially for a highly nonelliptic case such as the examples of Δy = 0 and 6 mm in Figure 5. It is observed in these two cases, the elliptic contacts are narrow along the running direction since at the rigid contact point curvatureA is less thanB, while the true contact shape is flat since the gap along the transversal direction does not follow an elliptic form. The shape sizes are, therefore, very different.One exception is the case of Δy = -2 mm where curvatureA is larger thanB at the rigid contact point. In this case, the contact force predicted by elliptic gap still agrees well with the true result, though the true shape is not a strict ellipse. The reason behind this phenomenon is that the semilength of ellipse in running direction can agree with nearly all contact cases since longitudinal curvature of the wheel keeps constant. Thus the area of ellipse whose length is smaller than its width can achieve an agreement with the true one; see example of Figure 5(a). It should be mentioned that the difference between the elliptic one and the true one in Figures 5(b) and 5(d) is due to different penetration as explained in Section 2.3.1.Therefore, we recommend elliptic gap used for vehicle-track dynamic modelling in the following two conditions: the gap follows the elliptic form or similar (see Figure5(c)), and curvatureA is larger thanB at the rigid contact point (see Figure 5(a)). The level of ‘similar’ and ‘larger’ should be calibrated in the future. Otherwise, nonelliptic gap must be used instead. Only large spin case the asymmetric gap is necessary, which is usually for vehicle derailment analyses.


## 3.1. Normal Contact

This section intends to investigate the effect of geometrical gap assumptions on normal contact modelling. A set of lateral shift cases with the prescribed load of 83.3 kN are conducted. With the increase of lateral displacement, the contact patch shift towards the rail gauge as shown in Figure4, of which the nonelliptic gap is employed.Figure 4
Variation of nonelliptic contact shape and pressure with lateral shift varying from -3 to 7 mm (from left to right) in 1 mm increment.


### 3.1.1. Geometric Gap, Contact Shape, and Pressure

Typical nonelliptic contact conditions where the wheelset has a lateral displacement of -2, 0, and 6 mm, and a similar elliptic case of Δy = 3 mm is chosen for a detailed comparison. Figure 5 illustrates the comparison of contact shape, geometric gap (atx= 0), and pressure. Two yaw angles, namely, 5 and 25 mrad, are considered to express asymmetric gap to different level. It is mentioned that 25 mrad means a relatively large yaw angle for the actual vehicle operation, which can only appear for the small curving (below 100 m) operational condition.Figure 5
Comparison of contact shape (left), gap (middle), and pressure (right) under different lateral shift.
(a)
Δy = -2 mm (b)
Δy = 0 mm (c)
Δy = 3 mm (d)
Δy = 6 mmFor Δy = -2 mm case in Figure 5(a), the lengths of contact patch in these four geometric gaps are almost the same, and their widths only differ a bit. By contrast, the pressure variations within the contact patch for elliptic shape and the other three differ a lot since an abrupt change of curvature appears aty = -8.5 mm. In the second case, i.e., wheelset central position, the contact shape and pressure predicted by elliptic gap and the other three are violently different because the gaps of the latter ones in positive direction differ a lot from the elliptic one. In this case, the pressure for the elliptic gap is higher due to smaller contact area. In Figure 5(c), contact patch and pressure distribution in these four gaps agree very well since they all follow the elliptic form. For the last case, the contact patch is composed of two merging ellipses and pressure appears two peaks. The result based on elliptic gap assumption fails to predict both contact patch and pressure as expected.It is further found the influence of yaw angle on contact shape and pressure in the first three shift cases can be ignored, though a large yaw angle of 25 mrad has been enforced. By contrast, the effect of this yaw angle for case Δy = 6 mm is significant. This phenomena may be explained as contact points for the first three cases located in the rail top and wheel tread that are relatively ‘flat’, while that in the last case locates near the rail gauge and wheel flange that are ‘sharp’. In the condition of ‘sharp’, the curvature varies sharply and is sensitive to yaw angle like in Figure 5(d).


### 3.1.2. Curve of Penetration versus Contact Force

Then we evaluate the performance of geometric gap simplification on resultant force, which is very important for vehicle dynamics simulation. Figure6 compares the penetration versus contact force curve for the four shift cases. Since the result of 5 mrad yaw angle is almost the same as that using nonelliptic gap, only the result for 25 mrad is presented for the asymmetric case. For both Δy = -2 and 3 mm, the curves predicted by different gap assumptions coincide very well, although in the former case the pressure distribution are different. In the other two shift cases, there is an obvious difference among elliptic shape and the other two. However, it may not significantly affect their performance on predicting normal force as explained in the following.Figure 6
Penetration versus contact force curve under different lateral shift.
(a)
Δy = -2 mm (b)
Δy = 0 mm (c)
Δy = 3 mm (d)
Δy = 6 mmBefore a vehicle-track dynamics computation, it needs a preprocessing to balance the prescribed load and then the penetrationδi is obtained depending on contact model for instance in Figure 6(d). If the vehicle is running on a steady-state, the contact force and penetration will keep constant. Otherwise, the contact force will vibrate due to irregularity between the wheel-rail interfaces. The dynamic forceFd can be obtained from (24)Fd=fδi+Δδ≈Fs+Δδ·tan⁡αiwhere subscripti = 1~3 meaning the elliptic, nonelliptic, and asymmetric gap;Fs is the prescribed load; Δδ is irregularity size; α is the local angle as shown in Figure 6(d).For the case Δy = 6 mm, α1 is smaller than the other two while α2 and α3 are nearly the same. So the dynamic forces predicted by nonelliptic and asymmetric gaps will be similar and that using elliptic gap may be underestimated.


## 3.1.1. Geometric Gap, Contact Shape, and Pressure

Typical nonelliptic contact conditions where the wheelset has a lateral displacement of -2, 0, and 6 mm, and a similar elliptic case of Δy = 3 mm is chosen for a detailed comparison. Figure 5 illustrates the comparison of contact shape, geometric gap (atx= 0), and pressure. Two yaw angles, namely, 5 and 25 mrad, are considered to express asymmetric gap to different level. It is mentioned that 25 mrad means a relatively large yaw angle for the actual vehicle operation, which can only appear for the small curving (below 100 m) operational condition.Figure 5
Comparison of contact shape (left), gap (middle), and pressure (right) under different lateral shift.
(a)
Δy = -2 mm (b)
Δy = 0 mm (c)
Δy = 3 mm (d)
Δy = 6 mmFor Δy = -2 mm case in Figure 5(a), the lengths of contact patch in these four geometric gaps are almost the same, and their widths only differ a bit. By contrast, the pressure variations within the contact patch for elliptic shape and the other three differ a lot since an abrupt change of curvature appears aty = -8.5 mm. In the second case, i.e., wheelset central position, the contact shape and pressure predicted by elliptic gap and the other three are violently different because the gaps of the latter ones in positive direction differ a lot from the elliptic one. In this case, the pressure for the elliptic gap is higher due to smaller contact area. In Figure 5(c), contact patch and pressure distribution in these four gaps agree very well since they all follow the elliptic form. For the last case, the contact patch is composed of two merging ellipses and pressure appears two peaks. The result based on elliptic gap assumption fails to predict both contact patch and pressure as expected.It is further found the influence of yaw angle on contact shape and pressure in the first three shift cases can be ignored, though a large yaw angle of 25 mrad has been enforced. By contrast, the effect of this yaw angle for case Δy = 6 mm is significant. This phenomena may be explained as contact points for the first three cases located in the rail top and wheel tread that are relatively ‘flat’, while that in the last case locates near the rail gauge and wheel flange that are ‘sharp’. In the condition of ‘sharp’, the curvature varies sharply and is sensitive to yaw angle like in Figure 5(d).


## 3.1.2. Curve of Penetration versus Contact Force

Then we evaluate the performance of geometric gap simplification on resultant force, which is very important for vehicle dynamics simulation. Figure6 compares the penetration versus contact force curve for the four shift cases. Since the result of 5 mrad yaw angle is almost the same as that using nonelliptic gap, only the result for 25 mrad is presented for the asymmetric case. For both Δy = -2 and 3 mm, the curves predicted by different gap assumptions coincide very well, although in the former case the pressure distribution are different. In the other two shift cases, there is an obvious difference among elliptic shape and the other two. However, it may not significantly affect their performance on predicting normal force as explained in the following.Figure 6
Penetration versus contact force curve under different lateral shift.
(a)
Δy = -2 mm (b)
Δy = 0 mm (c)
Δy = 3 mm (d)
Δy = 6 mmBefore a vehicle-track dynamics computation, it needs a preprocessing to balance the prescribed load and then the penetrationδi is obtained depending on contact model for instance in Figure 6(d). If the vehicle is running on a steady-state, the contact force and penetration will keep constant. Otherwise, the contact force will vibrate due to irregularity between the wheel-rail interfaces. The dynamic forceFd can be obtained from (24)Fd=fδi+Δδ≈Fs+Δδ·tan⁡αiwhere subscripti = 1~3 meaning the elliptic, nonelliptic, and asymmetric gap;Fs is the prescribed load; Δδ is irregularity size; α is the local angle as shown in Figure 6(d).For the case Δy = 6 mm, α1 is smaller than the other two while α2 and α3 are nearly the same. So the dynamic forces predicted by nonelliptic and asymmetric gaps will be similar and that using elliptic gap may be underestimated.


## 3.2. Tangential Contact

In this section, an elaborated analysis of shear stress, stick-slip division, microslip, and wear distribution will be presented with full consideration of yaw angle. From (19) and Figure 5, we can see yaw angle devotes to two main changes of contact solution: lateral creepage and asymmetric contact. Therefore, both of them are taken into consideration in this part.Two shift cases, Δy= -2 and 6 mm, are chosen for evaluating the effect of yaw angle on tangential modelling. The reasons for selecting these two cases include the following:(1)
In the former case, both longitudinal creepage and spin are low as the contact point occurring in the rail top, thus facilitating showing the role of lateral creepage on tangential contact.(2)
For the latter case, the contact shape can be changed obviously and then its effect on tangential contact can be considered.


### 3.2.1. Shear Stress, Slip and Wear Distribution

The stick-slip division, shear stress distribution, and its direction for two shift cases are shown in Figure7. Only results of elliptic and asymmetric shape are included since that by nonelliptic shape is similar to the asymmetric one. For each figure, the slip zone is encircled by the red line and the rest part means stick zone. The vector points the direction of shear stress, and its length is proportional to the magnitude of shear stress.Figure 7
Comparison of shear stress distribution and stick-slip division for two shift cases.
(a)
Δy = -2 mm (b)
Δy = 6 mmIn Figure7(a), the contact shape is nearly not influenced by yaw angle within 5 mrad. It can be found the increase of yaw angle make the slip zone gradually dominates. In the meantime, the directions of shear stress change from longitudinal to transversal. Although the contact shapes are different, the change tendencies for these two cases are the same. The difference of location and magnitude of shear stress is due to different pressure distribution as shown in Figure 5(a).In Figure7(b), relative large longitudinal creepage and spin make the slip zone dominating the contact patch. So it can easily reach a full slip condition once yaw angle appears. However, in this case, it needs a larger yaw angle (or lateral creepage) to change the direction of shear stress. Although the contact shapes under these two gap assumptions are far different, the evolution of stick-slip ratio and shear stress’ direction is similar.For the two shift cases, their relative slip distributions along longitudinal (y=0) direction are illustrated in Figure 8. In the case Δy = -2 mm, the magnitude and distribution of relative slip coincide very well for the two contact shape assumptions. It is found that the increasing yaw angle can enlarge relative slip continuously even though it has achieved a full slip condition (see results of 10 mrad). In Figure 8(b), the relative slip can also achieve a good agreement though their contact shapes vary significantly. This is because the relative slip mainly depends on rigid slip as defined in (19) and the difference is initiated from deformation. Such a difference can be reduced if the full slip condition is reached since the deformation caused by saturated shear stress is similar.Figure 8
Relative slip comparison for (a) Δy = -2 mm and (b) Δy = 6 mm.
(a) (b)The wear distributions predicted by different geometrical gap assumptions are compared in Figure9. In these results, lateral creepage keep the same for three gap models. In addition, the running speed of vehicle is 300 km/h. We can find three main phenomena: (1) there is a significant difference in terms of wear distribution and maximum wear depth among elliptic solution and the others as a result of shear stress and slip in Figures 7 and 8; (2) the yaw angle can continuously worse wear depth due to slip increase; (3) the difference between nonelliptic and asymmetric gap assumption is nonnegligible under case Δy = 6 mm and 25 mrad.Figure 9
Wear distributions under the same combination of creepage and spin.
(a)
Δy = -2 mm (yaw angle = 0, 1, 5 mrad from left to right) (b)
Δy = 6 mm (yaw angle = 0, 5, 25 mrad from left to right)


### 3.2.2. Creep Curve

Three gap simplifications in Figure6 are employed in this part to evaluate their performance on creep curve. Only pure longitudinal creepage or spin is considered ignoring lateral creepage initiated from yaw angle and thus facilitates understanding the effect of asymmetric contact shape on creep force estimation.Figures10(a) and 10(b) illustrate the longitudinal creep force curve due to pure longitudinal creepage for the three contact cases. Here, creep force is normalized by the normal force. As can been seen, the nonelliptic force estimation agrees with the reference for both shift cases, while the elliptic estimates a lower initial slope of the force build-up for the latter shift case.Figure 10
Creep curve estimation for three contact shape under (a), (b) pure longitudinal and (c), (d) pure spin condition.
(a)
Δy = -2 mm (b)
Δy = 6 mm (c)
Δy = -2 mm (d)
Δy = 6 mmIn addition to the pure creepage case, creep curve for the pure spin condition is also investigated. In Figures10(c) and 10(d), the saturated longitudinal and lateral creep force were estimated by both elliptic and nonelliptic gap deviates from the reference. This is attributed to double or single symmetry of the patch adopted separately in these two gap models. In the case of Δy = 6 mm, the estimation error is large due to obvious difference among the contact shapes predicted by three gap assumptions.


## 3.2.1. Shear Stress, Slip and Wear Distribution

The stick-slip division, shear stress distribution, and its direction for two shift cases are shown in Figure7. Only results of elliptic and asymmetric shape are included since that by nonelliptic shape is similar to the asymmetric one. For each figure, the slip zone is encircled by the red line and the rest part means stick zone. The vector points the direction of shear stress, and its length is proportional to the magnitude of shear stress.Figure 7
Comparison of shear stress distribution and stick-slip division for two shift cases.
(a)
Δy = -2 mm (b)
Δy = 6 mmIn Figure7(a), the contact shape is nearly not influenced by yaw angle within 5 mrad. It can be found the increase of yaw angle make the slip zone gradually dominates. In the meantime, the directions of shear stress change from longitudinal to transversal. Although the contact shapes are different, the change tendencies for these two cases are the same. The difference of location and magnitude of shear stress is due to different pressure distribution as shown in Figure 5(a).In Figure7(b), relative large longitudinal creepage and spin make the slip zone dominating the contact patch. So it can easily reach a full slip condition once yaw angle appears. However, in this case, it needs a larger yaw angle (or lateral creepage) to change the direction of shear stress. Although the contact shapes under these two gap assumptions are far different, the evolution of stick-slip ratio and shear stress’ direction is similar.For the two shift cases, their relative slip distributions along longitudinal (y=0) direction are illustrated in Figure 8. In the case Δy = -2 mm, the magnitude and distribution of relative slip coincide very well for the two contact shape assumptions. It is found that the increasing yaw angle can enlarge relative slip continuously even though it has achieved a full slip condition (see results of 10 mrad). In Figure 8(b), the relative slip can also achieve a good agreement though their contact shapes vary significantly. This is because the relative slip mainly depends on rigid slip as defined in (19) and the difference is initiated from deformation. Such a difference can be reduced if the full slip condition is reached since the deformation caused by saturated shear stress is similar.Figure 8
Relative slip comparison for (a) Δy = -2 mm and (b) Δy = 6 mm.
(a) (b)The wear distributions predicted by different geometrical gap assumptions are compared in Figure9. In these results, lateral creepage keep the same for three gap models. In addition, the running speed of vehicle is 300 km/h. We can find three main phenomena: (1) there is a significant difference in terms of wear distribution and maximum wear depth among elliptic solution and the others as a result of shear stress and slip in Figures 7 and 8; (2) the yaw angle can continuously worse wear depth due to slip increase; (3) the difference between nonelliptic and asymmetric gap assumption is nonnegligible under case Δy = 6 mm and 25 mrad.Figure 9
Wear distributions under the same combination of creepage and spin.
(a)
Δy = -2 mm (yaw angle = 0, 1, 5 mrad from left to right) (b)
Δy = 6 mm (yaw angle = 0, 5, 25 mrad from left to right)


## 3.2.2. Creep Curve

Three gap simplifications in Figure6 are employed in this part to evaluate their performance on creep curve. Only pure longitudinal creepage or spin is considered ignoring lateral creepage initiated from yaw angle and thus facilitates understanding the effect of asymmetric contact shape on creep force estimation.Figures10(a) and 10(b) illustrate the longitudinal creep force curve due to pure longitudinal creepage for the three contact cases. Here, creep force is normalized by the normal force. As can been seen, the nonelliptic force estimation agrees with the reference for both shift cases, while the elliptic estimates a lower initial slope of the force build-up for the latter shift case.Figure 10
Creep curve estimation for three contact shape under (a), (b) pure longitudinal and (c), (d) pure spin condition.
(a)
Δy = -2 mm (b)
Δy = 6 mm (c)
Δy = -2 mm (d)
Δy = 6 mmIn addition to the pure creepage case, creep curve for the pure spin condition is also investigated. In Figures10(c) and 10(d), the saturated longitudinal and lateral creep force were estimated by both elliptic and nonelliptic gap deviates from the reference. This is attributed to double or single symmetry of the patch adopted separately in these two gap models. In the case of Δy = 6 mm, the estimation error is large due to obvious difference among the contact shapes predicted by three gap assumptions.


## 3.3. Discussion of This Paper

From aforementioned results, we can find yaw angle plays a double effect on wheel-rail rolling contact solution, i.e., lateral creepage and contact shape. It yields a question about suited contact modelling in vehicle-track dynamics and wheel-rail contact mechanics. Note that discussion below is based on the assumption wheel-rail profiles that are smooth. The irregularity case is beyond the applicability of half-space method employed in this paper.


### 3.3.1. Effective Geometry Gap and Contact Model for Wheel-Rail Contact

For wheel-rail contact modelling, the effect of yaw angle on normal contact such as shape size and pressure is negligible even for a highly nonelliptic contact condition such as in Figure5(b). By contrast, tangential solutions such as stick-slip division and traction direction can be significantly influenced by yaw angle due to the appearance of lateral creepage; see Figure 7. If the true lateral creepage is taken into account; solution using the nonelliptic geometry gap assumption can coincide well with that in the asymmetric contact case (the true solution). This is verified from the comparison of wear performance (see Figure 9) that consider the combined effect of yaw angle on contact shape and lateral creepage. So in most wheel-rail contact mechanics cases, mainly including wheel tread-rail top contact, nonelliptic geometry gap simplification is sufficient. For the severe contact condition, an exact asymmetric geometry gap is needed. The term ‘severe contact’ here mainly means contact locates on flange-rail gauge contact where the contact shape is not an ellipse or similar like shift case Δy = 6 mm with yaw angle of 25 mrad.For modelling approach, Kalker’s Contact is used in this paper. It is, however, generally considered improper for online vehicle-track dynamics modelling as well as wear prediction for relatively high computational cost. As an alternative solution, simplified nonelliptic contact models [12–14] are regarded as fast and suitable for general nonelliptic cases, while it is considered not sufficiently exact for highly nonelliptic case [31]. The best choice may be a mixed contact models strategy for different contact cases, that is, (1) Hertz theory [8] and FASTSIM [9] used in similar elliptic case for instance Δy = 3 mm; (2) simplified nonelliptic models like KP [12], STRIPES [13] and ANALYN+FaStrip [14] used to treat general nonelliptic cases such as Δy = -2 mm; and (3) employ Contact to solve highly nonelliptic or asymmetric case (see Δy = 0 and 6 mm). To apply this strategy, a calibration work of nonelliptic level should be performed. A simple measure may be the ratio of left and right contact boundary determined by a scaling penetration similar to that in KP model [12], identifying to use which contact model. This work is within our scope of further study.


### 3.3.2. Effective Geometry Gap for Vehicle-Track Interaction

According to Figures6 and 10, the contact solution difference reflected on contact force is however acceptable for elliptic and nonelliptic gap simplification as explained in Section 3.1.2.For a higher accuracy, nonelliptic gap should be taken into account especially for a highly nonelliptic case such as the examples of Δy = 0 and 6 mm in Figure 5. It is observed in these two cases, the elliptic contacts are narrow along the running direction since at the rigid contact point curvatureA is less thanB, while the true contact shape is flat since the gap along the transversal direction does not follow an elliptic form. The shape sizes are, therefore, very different.One exception is the case of Δy = -2 mm where curvatureA is larger thanB at the rigid contact point. In this case, the contact force predicted by elliptic gap still agrees well with the true result, though the true shape is not a strict ellipse. The reason behind this phenomenon is that the semilength of ellipse in running direction can agree with nearly all contact cases since longitudinal curvature of the wheel keeps constant. Thus the area of ellipse whose length is smaller than its width can achieve an agreement with the true one; see example of Figure 5(a). It should be mentioned that the difference between the elliptic one and the true one in Figures 5(b) and 5(d) is due to different penetration as explained in Section 2.3.1.Therefore, we recommend elliptic gap used for vehicle-track dynamic modelling in the following two conditions: the gap follows the elliptic form or similar (see Figure5(c)), and curvatureA is larger thanB at the rigid contact point (see Figure 5(a)). The level of ‘similar’ and ‘larger’ should be calibrated in the future. Otherwise, nonelliptic gap must be used instead. Only large spin case the asymmetric gap is necessary, which is usually for vehicle derailment analyses.


## 3.3.1. Effective Geometry Gap and Contact Model for Wheel-Rail Contact

For wheel-rail contact modelling, the effect of yaw angle on normal contact such as shape size and pressure is negligible even for a highly nonelliptic contact condition such as in Figure5(b). By contrast, tangential solutions such as stick-slip division and traction direction can be significantly influenced by yaw angle due to the appearance of lateral creepage; see Figure 7. If the true lateral creepage is taken into account; solution using the nonelliptic geometry gap assumption can coincide well with that in the asymmetric contact case (the true solution). This is verified from the comparison of wear performance (see Figure 9) that consider the combined effect of yaw angle on contact shape and lateral creepage. So in most wheel-rail contact mechanics cases, mainly including wheel tread-rail top contact, nonelliptic geometry gap simplification is sufficient. For the severe contact condition, an exact asymmetric geometry gap is needed. The term ‘severe contact’ here mainly means contact locates on flange-rail gauge contact where the contact shape is not an ellipse or similar like shift case Δy = 6 mm with yaw angle of 25 mrad.For modelling approach, Kalker’s Contact is used in this paper. It is, however, generally considered improper for online vehicle-track dynamics modelling as well as wear prediction for relatively high computational cost. As an alternative solution, simplified nonelliptic contact models [12–14] are regarded as fast and suitable for general nonelliptic cases, while it is considered not sufficiently exact for highly nonelliptic case [31]. The best choice may be a mixed contact models strategy for different contact cases, that is, (1) Hertz theory [8] and FASTSIM [9] used in similar elliptic case for instance Δy = 3 mm; (2) simplified nonelliptic models like KP [12], STRIPES [13] and ANALYN+FaStrip [14] used to treat general nonelliptic cases such as Δy = -2 mm; and (3) employ Contact to solve highly nonelliptic or asymmetric case (see Δy = 0 and 6 mm). To apply this strategy, a calibration work of nonelliptic level should be performed. A simple measure may be the ratio of left and right contact boundary determined by a scaling penetration similar to that in KP model [12], identifying to use which contact model. This work is within our scope of further study.


## 3.3.2. Effective Geometry Gap for Vehicle-Track Interaction

According to Figures6 and 10, the contact solution difference reflected on contact force is however acceptable for elliptic and nonelliptic gap simplification as explained in Section 3.1.2.For a higher accuracy, nonelliptic gap should be taken into account especially for a highly nonelliptic case such as the examples of Δy = 0 and 6 mm in Figure 5. It is observed in these two cases, the elliptic contacts are narrow along the running direction since at the rigid contact point curvatureA is less thanB, while the true contact shape is flat since the gap along the transversal direction does not follow an elliptic form. The shape sizes are, therefore, very different.One exception is the case of Δy = -2 mm where curvatureA is larger thanB at the rigid contact point. In this case, the contact force predicted by elliptic gap still agrees well with the true result, though the true shape is not a strict ellipse. The reason behind this phenomenon is that the semilength of ellipse in running direction can agree with nearly all contact cases since longitudinal curvature of the wheel keeps constant. Thus the area of ellipse whose length is smaller than its width can achieve an agreement with the true one; see example of Figure 5(a). It should be mentioned that the difference between the elliptic one and the true one in Figures 5(b) and 5(d) is due to different penetration as explained in Section 2.3.1.Therefore, we recommend elliptic gap used for vehicle-track dynamic modelling in the following two conditions: the gap follows the elliptic form or similar (see Figure5(c)), and curvatureA is larger thanB at the rigid contact point (see Figure 5(a)). The level of ‘similar’ and ‘larger’ should be calibrated in the future. Otherwise, nonelliptic gap must be used instead. Only large spin case the asymmetric gap is necessary, which is usually for vehicle derailment analyses.


## 4. Conclusions

This paper presents a 3D geometrical gap searching method, which constructs a bridge between contact point detection method and 3D non-Hertz rolling contact theory. Such a methodology can be used to evaluate asymmetric contact due to yaw angle. From presented results, the following conclusions can be drawn.(1)
For the contact between wheel tread and rail top, the effect of yaw angle on contact shape and pressure is negligible, while this effect for a wheel flange/rail gauge contact case is evident especially for a large yaw angle.(2)
Yaw angle plays an important role on tangential contact solution as well as wear distribution even a full slip condition is achieved.(3)
The slip velocity along the running direction can hardly be affected by contact shape, since its value is mainly dominated by creepage and spin.(4)
From an exact modelling viewpoint, the elliptic gap idealisation produces an inaccurate contact solution for a nonelliptic contact case. By contrast, it can always provide an acceptable resultant normal/creep force including the nonelliptic case where curvatureA is larger thanB.(5)
Contact mechanics solution from the nonelliptic gap idealisation can agree well with the true solution in most wheel-rail contact cases with yaw angle. However, a highly nonelliptic contact case with a large yaw angle or spin always needs an exact asymmetric gap considered.The present gap searching method is based on half-space assumption, i.e., the contact patch is located on a tangential plane. Approach solving conformal contact (curved contact patch) should be considered in the further study.



---

*Source: 9895267-2019-04-02.xml*