Influence of white etching layer on rolling contact behavior at wheel-rail interface

The existence of narrow and brittle white etching layers (WELs) on the rail surface is often linked with the formation of rail defects such as squats and studs, which play the key roles in rail surface degradation and tribological performance. In the present study, a systematic investigation on stress/strain distribution and fatigue life of the WEL during wheel-rail rolling contact was conducted based on a numerical model considering the realistic wheel geometry. This is the first study considering the influence of rail materials, loading pressure, frictional condition, WEL geometry (a/b), and slip ratio (Sr) in the practical service conditions at the same time. The results revealed much higher residual stress in WEL than in rail matrix. Stress changes along the rail depth matched with the previously reported microstructure evolutions. The current work revealed that the maximum difference in contact stress between the wheel passages of rail matrix and the WEL region (noted as stress variation) rises with the increase of loading pressure, the value of a/b, and Sr; but drops with the friction coefficient (μ). In addition, a critical length-depth ratio of 5 for a/b has been found. The fatigue parameter, FP, of the WEL decreased quickly with the length-depth ratio when it was less than 5 and then increased slightly when it was larger than 5. This study also revealed that the fatigue life of the WEL was reduced for high strength head hardened (HH) rail compared with standard carbon (SC) rail.


Introduction
Squat is one of the typical rolling contact fatigue (RCF) defects of rails [1][2][3][4][5][6][7][8][9], which plays a key role in rail surface degradation and tribological performance. The first report of squat can be dated back to as early as 1950s in Japan [8]. Then Johnson [7], Bower and Johnson [2], and other researchers [1,3,6] from European countries have further contributed to establish the fundamental understandings on the mechanics of initiation and development of surface-initiated squat cracks in rails. In 1980s, Clayton et al. [4,5] have studied the surface damage of rails (particularly RCF and squats) from a metallurgical view. Since 1990s, preventative maintenance against rail squats have been gradually proposed and developed by Smallwood et al. [9], Kalker et al. [10], and Grassie [11,12]. So far, a number of investigations on squats have already been reported, mainly based on metallurgical evaluations or numerical analysis [13][14][15][16][17][18]. It has been found that squat formation on rail surface can be affected by many factors and a more comprehensive understanding is still necessary.
Based on the recent researches [13,14,[19][20][21][22][23], the occurrence of narrow and brittle white etching layer (WEL), which is named due to the white reflection after being etched in 2%-5% HNO 3 in ethanol, is thought to promote the formation of rail surface squats. Steenbergen and Dollevoet [24] proposed a theory for the origination and physical nature of squat defects on train rails and they pointed out two metallurgical principles of crack initiation: edge delamination of the WELs embedded in the rail surface and transverse fracture of the WELs. Grassie et al. [25] proposed the concept of "stud", which develops without severe plastic flow and associates with WELs at the top layer of rail at a lower borne tonnage. Compared with a typical RCF squat, a stud may show less tendency to cause rail fracture. Metallurgical examination and synchrotron study of ex-service damaged rail by Al-Juboori et al. [19] revealed a close relationship between WELs and squat initiation. Up to now, two theories are widely acknowledged to explain the formation mechanism of the WELs: one is the super saturation of carbon during repeated severe plastic deformation on rail surfaces [26,27]; the other is thermal induced martensite phase transformation due to heating above the austenitization temperature region and then followed by rapid cooling [21,23,28]. Very recently, stratified layers were reported by Li et al. [13] and Messaadi and Steenbergen [29], where a "brown etching layer (BEL)" formed immediately beneath the WEL with comparable properties. The formation mechanism of BEL was similar to the WEL, which was caused by either phase transformation or severe plastic deformation [13,29].
In order to obtain a detailed understanding on the nature and formation mechanism of WELs, a number of laboratory experiments and metallurgical examination of ex-serviced rails have been carried out. For instance, Makino et al. [30] investigated the formation of WELs using a set of two disk type test pieces (TPs) which was applicable to approximate Hertzian rolling contact. Carroll and Beynon [31] fabricated the WELs using a twin disc test under elastic-plastic deformation. These two studies revealed that WEL thickness was influenced by the hardenability of original rail material and the crack morphology in the vicinity of WEL was dependent on the plastic deformation. In Ref. [20], Baumann and co-authors found that the WEL formation on rail could be assisted or even produced by high thermal stress during mechanical ball milling. Besides, Vargolici et al. [32] found that the WELs were very brittle and about three times harder than the traditional rail steel. With the help of synchrotron X-ray diffraction, it has been found that the dislocation density in the WELs was about 10 12 cm −2 and compressive residual stress was over 700 MPa [33,34]. From metallurgical observations of ex-service rail samples, Al-Juboori et al. [19] found two distinguishable types of WELs based on different operational conditions, namely WEL at a heavy braking track region containing martensite and retained austenite, and WEL at a track region under steady traffic speed consisting of nanocrystalline martensite. Recently, Li et al. [13] has found a brown etching transition layer between the WEL and matrix, which may play an important role in crack initiation.
In addition to the experimental observations, a number of numerical simulations have also been conducted to understand the formation mechanism and fatigue properties of WEL during wheel rail contact. Bernsteiner et al. [21] simulated the temperature distribution within and below the rail surface in wheel/rail contact zone. Their results suggested that the austenitization temperature of rail steel could be reached under certain conditions and the thermal induced WEL was possible. Fatigue life of WEL was evaluated by Seo et al. [35] using a finite element simulation. It has been found that the shortest fatigue life located in the leading point of WEL. Kato et al. [36] also studied the WEL fatigue life by an elastic-plastic finite element model. They reported that more cracks initiated from the WEL than rail matrix and the maximum stresses decreased with the WEL size. However, it is worth noting that both these studies [35,36] have a significant drawback and they have introduced a down-scale setup, unpractical small train wheel, in order to save the simulation time and fine meshing of WEL. The downscaled wheel geometry can cause great errors in contact area and stress distribution in both wheel and rail according to the Hertz contact theory [7], which will then significantly affect stress variations around the WEL and the related RCF properties. According to Refs. [37][38][39][40], both the lubricated condition and slip ratio (Sr) played vital roles in wheel and rail wear. Even though the influences of friction and Sr on wheel/rail interface have been extensively studied [21,33,35,[37][38][39][40][41], the effects of friction and Sr on the RCF properties of the WEL were not mentioned in these reports and are still missing.
This study aims to provide a systematic numerical investigation on the wheel-rail rolling contact with consideration of a WEL on rail surface. Influences of loading pressure, friction coefficient (μ), WEL geometry (a/b), and Sr on stress distribution and fatigue behavior of both rail and WEL have been studied considering the realistic wheel geometry. Multi-wheel passages were introduced to investigate the RCF life of rail and WELs. According to Refs. [42][43][44], rail material property also played a critical role in the wear and corrugation in practical train service conditions and therefore harder rail grades were widely introduced to replace softer rail grades. Unfortunately, influence of different rail materials with consideration of the WEL has never been mentioned up to now. Therefore, two typical rail materials (head hardened (HH) and standard carbon (SC)) with large difference in strength were specifically examined in this study in order to understand the influence of rail material. The present modelling was verified by comparing the simulated residual stresses in both WEL and rail matrix with those experimental measurements by synchrotron X-ray diffraction method [33,34]. The microstructure evolution along rail depth direction observed in Ref. [19] was successfully explained with the simulated stress gradients in this study. In addition, rolling contact fatigue lives of both the WEL and rail matrix have been investigated in detail under different conditions.

Material and properties
The studied wheel material was 0.7% carbon steel with pearlite microstructure, having a yield strength of about 1.1 GPa and the Vickers hardness of 330 HV at room temperature [36]. Two typical rail materials were studied: one was a SC rail steel and another one was a HH rail steel. The corresponding yield strength was about 507 and 800 MPa, respectively. Compared to the wheel and rail materials, the WELs were much harder and more brittle. Vickers hardness of the WELs was reported to vary between 550 and 1,200 HV [19-23, 32-36, 45]. It is worth noting that the hardness of WELs can be significantly affected by its origin (plastic deformation or thermal induced), microstructure, and phase constitutions. In the present study, a thermal induced WEL observed on an ex-service rail surface was considered and its hardness was measured to about 840 HV [19]. The yield strength and reduction of area (RA) of the WEL was about 1.39 GPa and 1.3%, respectively, obtained from a micro-tensile testing [36]. Table 1 lists detailed material properties of the wheel, rail, and WEL used in the finite element simulations in this paper. These parameters were summarized from Refs [19,35,36], in which rail steels (SC and HH) had the same pearlite microstructure, similar carbon range (0.65%-0.82%), and similar hardness (300-400 HV).  [22] and Australia, respectively. The first type of WEL ( Fig. 1(a)) was obtained from U71Mn rail which was located in a curved track region from Shenyang to Jilin city of China. The formation mechanism of this WEL was suggested as thermally induced martensite [22]. The second type of WEL ( Fig. 1(b)) was obtained from a damaged rail provided by Sydney Trains from re-railing sites in New South Wales, Australia. According to our previous work [19], the second type of WEL was induced by either severe plastic deformation or thermo-mechanically phase transformation. However, the formation mechanism of WEL was not primary interest in this study, the RCF behaviors of rail surface with pre-existence of WEL was mainly investigated. Metallurgical observations showed that the WELs formed on rail surface in different shapes and sizes [33]. In the area where the wheel contacted steady with the rail, the WELs appear coherently as bands or strips along the rail surface, such kind of WEL feature can be found in Refs. [19,20,29,33,46]. However, the WELs also can be formed in isolated patches [19,20,22,33,35] or discrete islands [24]. Crack failures were more likely to occur at the front or tail of the WELs rather than the middle of the WELs. Hence, such discrete patch or island like WELs may have even worse effects on the rail service life than those continuous bands or strips like WELs. In the longitudinal crosssection of rail surface, the microstructure features of the patch WEL can be seen in Figs. 1(a) and 1(b), a clear boundary between the WEL and pearlite matrix can be found. Moreover, those discrete patch or island like WELs have a common geometry feature, namely a circular arc as highlighted by the red lines. Such feature has also been observed in other reported microstructural studies of the WELs [19,20,22,24,35]. The present study is targeted at those discrete patch or island like WELs and its influence on rolling contact behavior at wheel-rail interface.

Geometry of the WEL
An assumption has therefore been made to simplify the discrete patch or island like WELs geometries into circular arcs on the longitudinal cross-section as illustrated by Fig. 1(c), where X axis indicates the rolling/longitudinal direction, and Y axis indicates the normal/vertical direction. Assuming the maximum length of 2a and the maximum thickness of b, WEL can be mathematically described by the following equation: where the radius R can be calculated by As shown in Fig. 1(c), BC is the tangent line at the intersecting point B between the WEL head and rail surface. The angle ABC can be calculated by In addition, the angle ABD can be calculated by From Eqs. (2)-(4), the ratio of a/b is a critical parameter in describing the WEL geometry. This approximation of arc-shaped WEL in plane-strain condition has already been demonstrated to be satisfactory in Refs. [35,36].
In order to understand influence of the a/b, six typical geometries with a/b varying between 2 and 12.5 based on reported WEL microstructures have been simulated. The details can be found in Table 2.

Boundary and loading conditions
It has been widely accepted that, the wheel-rail contact can be simplified to a two-dimensional (2D) plane strain condition on a longitudinal cross-section in order to save computational cost because the contact region is very small compared to the geometry dimensions of both wheel and rail [35,36]. Figure 2 shows the schematic illustration of the model, where axes X, Y, and Z indicate the rolling/longitudinal, normal/vertical, and lateral directions, respectively. In the present study, we have developed the model based on ANZR1 wheel with a practical wheel radius of 460 mm on the X-Y plane and a representative width of 1 mm on the Y-Z plane; AS60 rail of 1,500 mm on the X-Y plane and a representative width of 1 mm on the Y-Z plane. As a large curvature of 190 mm exists on AS60 rail top along the Y-Z plane, a 2D assumption could be accepted to analyze the wheel/rail contact behavior on the X-Y plane in order to save computational cost as reported in previous studies [35,36]. However, it should be noted that the 2D simplification has a drawback because of the presence of WEL/matrix interface along rail transverse direction. Even though the stress and strain distributions along the rail depth direction may not be significantly affected, but a three-dimensional (3D) comprehensive model is still essential for a more accurate understanding in the future. The simulations were conducted using commercial finite element software ANSYS/LS-DYNA.
The mechanical solver applied in this study is based on Lagrangian formulation for wheel/rail contact problem, this solver is an explicit time integration scheme due to its stability in achieving accurate solutions with efficient computation. For explicit solution, initial contact conditions such as dynamic contact forces will inevitably occur when the wheel runs over the rail [15]. To avoid the unsteady fluctuation of contact force at initial contact point, the wheel rotated over a distance from initial contact to the targeted WEL zone; hence, a steady contact condition would be achieved for the further analysis.
Different from the down-scale setup in the previous simulations [35,36], a full-scale dimension of an ANZR1 wheel with a radius of 460 mm was considered in the present study and very different results are expected. The simulated rail length was set to 1,500 mm and the height was set to 200 mm, respectively. An elasticplastic material model was applied to describe both the wheel and rail using the physical and mechanical properties listed in Table 1. As for the WEL, a comparison between the simulation results based on elastic-plastic material model and perfectly elastic material model revealed nearly no difference. Therefore, an elastic material model was mainly applied. As can be seen from Fig. 2, non-uniform meshes were conducted in the simulations. Regions close to the surface in wheel and rail and around the WEL had  higher density of meshes with very fine elements. The minimum element size was about 40 μm × 40 μm and the element size was progressively increased towards the far field domain. The total number of nodes and elements was 320,156 and 315,340 for the wheel, and 468,928 and 461,579 for the rail, respectively. The number of nodes and elements in the WEL was dependent on the geometry and summarized in Table 2.
For plastic deformation induced WELs, a deformed pearlite transition zone can be found beneath the WELs [19]. By contrast, no obvious transition region is present between the thermal induced WEL and base matrix [19,29]. In the present study, only the thermal induced WEL was simulated and no transition region was therefore assumed. Both ends of the rail were fixed in the rolling direction while the bottom of the rail was constrained in both the rolling and normal directions. Perfect interface bonding between the WEL and rail substrate was assumed. A constant vertical loading force of 13,000 N was applied in the wheel which corresponded to a maximum Hertzian contact pressure of 1.2 GPa. The Hertzian contact pressure is termed as "loading pressure" thereafter in this paper. The wheel was also assigned a rotation velocity ω equal to 70 rad/s according to a train operation speed of 110 km/h. When the wheel runs over the rail, a master-slave surface to surface contact scheme was used. Real multi-wheel passages instead of simplified multi-passage of moving pressure method, were considered, namely lifting the wheel when the wheel reaches to the rail end and then moving the wheel back to its initial position and then loading again. The data present in the paper are the results after six wheel passages. As listed in Table 3, four loading forces corresponding to loading pressures between 0.8 and 1.8 GPa, three friction coefficients between 0.1 and 0.5, five slip ratios between 0.5% to 4% were investigated. Table 3 shows a summary of simulation conditions investigated in this study.

Fatigue analysis
It has already been shown in Refs. [1-4, 6, 10, 14, 25, 26, 34, 37, 42, 47-50] that the rail is subjected to time- dependent, multi-axial, and mixed mode cyclic loading due to the repeated passages of the wheels, which can lead to the fatigue problem of rails. Therefore, fatigue analysis is essential to avoid accidents such as rail fracture and to improve the service life of wheel-rail system. There are various models for RCF life prediction [50] and they can be divided into the following groups: (1) equivalent strain approaches, (2) critical plane approaches, (3) energy and energy-density based models, (4) combined energy-density based and critical plane models, and (5) empirical models. This section provides a basic theoretical introduction of the critical plane method, which is a widely used fatigue evaluation method because of its capability to estimate life mostly within ±3 factors of life for smooth specimens [48,50]. This theory was firstly proposed by Kandil et al. [49] in the following form: where   max in Eq. (5)  This theory was later improved by Fatemi and Socie [47] as shown in Eq. (6), where the normal strain range   n was replaced by the maximum normal stress on the critical plane and influences of the mean stress and material hardening were considered to redefine the specific fatigue parameter, FP, as where  max n is the maximum normal stress on the critical plane,  y is the tensile yield strength of the material, and k is a material constant determined from axial and torsional fatigue experiments. Equation (6) can be rewritten as Eq. (7) by correlating FP with fatigue life: where f N indicates the number of cycles to failure, Due to the initial unsteady contact behavior introduced by an explicit time integration scheme, normal contact stress oscillates significantly during the early stage when the wheel moves on the rail surface. The oscillation decays gradually and then normal contact stress becomes stable for further analysis with a magnitude of about 1,050.6 MPa for a SC steel rail and about 1,140.1 MPa for a HH steel rail. However, a sharp rise of normal contact stress can be seen from Fig. 3 when the wheel passes the WEL in both cases. This can be attributed to the higher hardness of the WEL than the rail matrix as shown in Table 1, which has also been reported in Refs. [14, 19-23, 32-37, 39, 51]. Magnitude of maximum normal contact stress difference between the wheel passages of rail matrix and the WEL region is denoted as S (Fig. 3(a)). A large influence of rail material on S has been found. It is about 272.4 MPa when the rail material is SC steel, while it decreases to about 151 MPa when the rail material changes into a HH steel. According to this study, it has been found that a harder rail material results into a larger wheelrail normal contact stress in steady contact region but a smaller S between rail matrix and WEL. Influence of frictional condition on S has also been studied for a SC steel rail under loading pressure of 1.2 GPa and Sr of 0.5% as shown in Fig. 3(b), which shows a good linear relationship between S and μ. S is about 200 MPa when μ is 0.1. With increasing the μ, S increases quickly and it reaches to about 350 MPa for μ of 0.5. Figure 3(c) shows the influence of loading pressure on normal contact stress at a constant μ of 0.3 and Sr of 0.5%. As can be seen, a higher loading pressure leads to a more severe contact stress variation. S is about 220 MPa under a loading pressure of 0.8 GPa, and it then increases up to about 370 MPa when the loading pressure equals 1.8 GPa. The variation of contact stress was determined by the rail/WEL constitutive properties. When the wheel contacted with rail, the base matrix was plastically deformed, the stress was relatively low. While contacting with WEL, normal contact stress increased due to its elastic essence. Therefore, the ratio of contact patch to the WEL area has no significant influence on the variation of normal contact stress.

Stress and strain distributions around the WEL
The stress and strain distribution fields around the WEL have been studied for a SC steel rail under the loading pressure of 1.2 GPa, μ of 0.3, and Sr of 0.5%. The simulated WEL has a geometry with a/b = 5, which is corresponding to a WEL with thickness of 200 μm and length of 2 mm. As can be seen from Fig. 4(a), the maximum of effective stress (von-Mises stress) is about 726.5 MPa when the wheel passes the rail matrix. On the other hand, passage of the WEL leads to an increment of about 69% and the maximum of effective stress reaches to about 1,228.1 MPa, as seen in Fig. 4(b). Similar stress localization has also been reported by Seo et al. [35]. However, it should be mentioned that the stresses between this study and their report are not comparable. There are several reasons to interpret the differences. First, the current study shows the result after achieving a steady contact condition, while the wheel only rotated 30° in Ref. [35] and it might not have reached a steady contact region. According to Fig. 3(a), there is a significant difference between the oscillation region and the steady contact region. Second, the simulated wheel with a radius of 50 mm in Ref. [35] was very small and not describing its practical geometry, which affects the stress and strain states significantly according to the Hertzian contact theory [7]. Compared with the current study, results of contact patch size and the stress influence depth were smaller in Ref. [35] due to the application of a down-scaled wheel. In current analysis, the calculated semi-axis of contact patch was 6.6 mm but it was just 4.4 mm in Ref. [35] under the same loading pressure. Besides, the depth of peak stress occurred in rail  | https://mc03.manuscriptcentral.com/friction subsurface for current study and Ref. [35] was 0.29 and 0.17 mm, respectively, which showed a large difference. In addition, the loading, mesh condition, and a/b also differs between the two studies. The results suggested that simulating the realistic wheel geometry is very important to conduct the stress and strain analysis in wheel-rail contact.
In order to study the stress and strain evolution histories at rail surface around a typical WEL with 2 mm in length and 200 μm in thickness, seven nodes at the rail surface with different distance (d) as shown in Fig. 5(a) have been selected, where d indicates the horizontal distance from the selected node to the center of the WEL surface. Sign '+' means along the X direction while sign '-' means opposite to the X direction. The wheel will roll and contact with the rail matrix at a position of d = -2 mm, then interact with at leading edge of the WEL at d = -1 mm, then contact with the WEL at position d = -0.5, 0, and 0.5 mm, followed by rolling over the trailing edge of WEL at d = 1.0 mm and rail matrix at d = 2.0 mm. Figure 5(a) shows that effective stress at node d = -2 mm increases quickly up to 726.5 MPa when contacting with wheel and then decreases to about 350 MPa after passage of the wheel. The stress evolution at node d = 2 mm is almost the same with that of d = -2 mm, which suggests that all nodes at rail matrix surface undergo similar stress evolution in steady contact condition. By contrast, the effective stress at the WEL surface is much larger during contact with the wheel. It is interesting to observe that the node at the leading edge of the WEL (d = -1 mm) has the largest maximum effective stress which is about twice larger than at the rail matrix surface. The maximum stress decreases slightly from 1,228.3 MPa at leading edge of the WEL (d = -1 mm) to 1,056.6 MPa  at the trailing edge of the WEL (d = 1 mm). However, the residual effective stress after passage of the wheel has an opposite tendency and increases gradually from negative to positive d positions as shown in Fig. 5(a). Figure 5(b) shows the evolution of effective strain at those seven selected positions. Different from the stress, the WEL surface undergoes much lower strain compared to rail matrix. The corresponding maximum effective strain is about 0.005 for WEL and 0.022 for rail matrix, respectively. In addition, no large difference in the maximum strain at the WEL surface is observed. Significant differences in stress and strain between rail matrix and WEL should be attributed to their material property, which has been considered as elastic-plastic for rail and perfectly elastic for WEL, respectively. Another simulation assuming the WEL as elasticplastic has also been conducted. The results are almost the same with Fig. 5, because the WEL has very large yield stress as described in Table 1 and wheel-WEL contact is within the elastic deformation region. This also confirms that taking the WEL as an elastic material is reasonable in the simulations.
It should be noted that the pre-existing residual stresses generated by thermal expansion differences between the WEL and rail matrix during cyclically heating and cooling or by microstructure changes during phase transformations are not considered in the current work. According to Ahlström [52], compressive radial stresses were generated within the martensitic layer while tensile stresses were formed beneath it. Effects of these pre-existing residual stresses on the rolling contact behavior of WEL are still not clear, and they will be investigated in the future work.
In order to understand the stress gradient in the rail along its depth direction, the residual stresses at seven selected positions have been compared in Fig. 6. It should be mentioned that depth less than 1 mm has only been studied here because the most significant stress variation is located in this region according to Fig. 4, while the stress decreases slowly to zero with further increasing the depth. As can be seen from Fig. 6(a) for rail matrix (d = -2 mm and d = 2 mm), the residual effective stress increases slightly from about 370 MPa at surface to 420 MPa at a depth larger   (Fig. 6(b)) and  y (Fig. 6(c)) and in residual shear stress  xy (Fig. 6(d)). It can be found that residual stress  x plays as a dominant role and is always in compression, while  y and  xy are very small.
Al-Juboori et al. [19] have analysed the microstructure changes along the rail depth. There are three different microstructural features: (1) heavily deformed fine grains with high dislocation densities in the WEL, (2) mixture of fractured cementites and dislocations in transition region, and (3) pearlite lamellar structures in the undeformed region. However, the stress states of those microstructures are not reported in Ref. [19]. In fact, those three positions for microstructural observations in Ref. [19] are corresponding to the three residual stress distribution regions as marked in Fig. 6. Severe deformation happened in region-I leading to the formation of nano grains in the WEL region, which should have very large compressive residual stress as shown in Fig. 6. Even though the current simulations do not consider the plastic deformation in the WEL region, very large stress over 1.1 GPa still exists as shown in Fig. 5(a). In a recent study, Arechabaleta et al. [53] have reported the dislocation evolutions in a low-alloy and interstitial-free steels in the pre-yield range of a tensile test. Their theory can explain the formation of high dislocation density in region-I after a large number of wheel-rail contact cycles. By comparison, the mixed microstructures in region-II and lamellar structures in region-III have much lower residual stresses. Therefore, those microstructure gradients and dislocation density evolutions from the rail surface can be well understood by examining their residual stress variations. In our future study, a dislocation based model will be developed to understand the deformation mechanism of those different microstructures.

Influence factors on stress variation between rail matrix and WEL
From Figs. 5 and 6, it has been found that stress changes significantly when the wheel moves from rail matrix to the WEL surface. It has also been reported in many studies that cracks were mainly initiated from the leading edge of the WEL [19,35,36]. Therefore, the stress variation between rail matrix and the WEL at this position and the influences of loading pressure, μ, a/b, and Sr during contact are investigated in this section. Figure 7(a) shows the distribution of effective stresses around the leading edge of the WEL during wheel-rail contact under four loading pressures varying from 0.8 to 1.8 GPa, with a μ of 0.3, a/b of 5, and Sr of 0.5%. Under each loading pressure, the maximum effective stress is located at the leading edge of the WEL, and the stress in the WEL is higher than in the rail matrix. The result in this study has similar tendency with Ref. [36], which had much smaller magnitude in stress due to the unrealistic wheel geometry in their simulation. With increasing the loading pressure, the effective stress increases in both the WEL and rail matrix. However, it increases more quickly in the WEL than in rail matrix. The maximum effective stress in the WEL reaches to about 1,900 MPa while it is only about 800 MPa in rail matrix. The stress variation at the leading edge of the WEL has been summarized in Fig. 7(b) as a function of loading pressure, and a good linear relationship can be observed. The effective stress difference increases significantly from 150 MPa under loading pressure of 0.8 GPa to about 1,100 MPa under loading pressure of 1.8 GPa. It should be noted that loading pressure of 1.8 GPa is not a practical service condition and it is examined here only for a systematical understanding.
Previous studies [37][38][39] reported how friction or lubricated condition affected the microstructure and wear at rail surface. However, they did not consider existence of the WEL. In order to understand the influence of friction on the stress distribution in the WEL and rail matrix, three cases with different μ varying from 0.1 to 0.5, with the constant loading pressure (1.2 GPa), WEL geometry (a/b = 5), and Sr (0.5%) were simulated. The corresponding effective stress variations at the leading edge of the WEL are shown in Fig. 8, which reveals larger effective stress in both the WEL and rail matrix when μ increases. By comparison, the WEL is more sensitive to friction than the rail matrix. For a given μ, very small changes in effective stress in the rail matrix with the distance from the leading edge of the WEL can be observed. However, there is a stress decrement in the WEL with distance from the leading of the WEL. From Fig. 8(b), the difference in effective stress at the leading edge of the WEL decreases slightly from 525 to 475 MPa when μ increases from 0.1 to 0.5. Figure 9 shows the simulated effective stresses for different WEL geometries. In rail matrix shown in Fig. 9(a), the average effective stress decreases very quickly from 840 to 690 MPa when a/b increases from 2 to 5, and the corresponding α decreases from 53.2° to 22.6° and  decreases from 26.6° to 11.3° according to Table 2. The parameters α and  show the same feature. With further increasing a/b or decreasing α, decrease of effective stress becomes much slower and it is about 610 MPa for a/b = 12.5 (or α = 9.2°). By contrast, there is an opposite influence in the WEL, where a/b leading to increase of effective stress has been observed. Figure 9(b) shows significant variation of the effective stress at leading edge of the WEL as a function of a/b. As can be seen, two parts can be distinguished. In the first part, variation of effective stress increases quickly with a/b when it is less than 5. Then, variation of effective stress rises much slower when a/b is above 5 and a good linear relationship can be found. In another word, effective stress decreases slowly with α when it is smaller than 22.6° but decreases much faster with α when it is larger than 22.6°. When a/b equals 12.5 or α equals 9.2°, the maximum effective stress variation of about 645 MPa is obtained.   Figure 10(a) shows the distributions of effective stress around the leading edge of the WEL for five different slip ratios from 0.5% to 3.5%, with loading pressure of 1.2 GPa, μ of 0.3, and a/b of 5. It has been found that Sr has very small influence and the effective stress in both the WEL and rail matrix increases very slightly when Sr increases. From Fig. 10(b), we can observe similar effective stress variation for different slip ratios. It only changes from 500 to 515 MPa when Sr increases from 0.5% to 3.5%. According to Ma et al. [38], Sr less than 3.83% had very small influence on wear of a rail material and there was no clear difference in the angle and depth of cracks under those slip ratios.

Results of fatigue analysis
It is widely accepted that rolling contact fatigue is a major cause of rail failure [15,16]. It is therefore fatigue analysis of WEL is conducted in this section.
In Fig. 11(a), N f has been predicted at eight selected positions around the WEL, where position P1 is located in rail matrix, positions P2-P6 are located at the WEL surface, and positions P7 and P8 are located at WELsubstrate interface as marked. It can be found that the WEL surface (P2-P6) has very low fatigue life less than 3 × 10 6 cycles for a SC steel rail. The leading edge of the WEL (P2) has the shortest life and it is only about 2.2 × 10 6 . Figure 11(a) shows similar tendency to Ref. [33]. Fatigue life at the leading edge of the WEL, P2, is the lowest, followed by the middle of WEL at P4, and the trailing edge of WEL at P6. By contrast, rail matrix (P1) and WEL-substrate (P7 and P8) have much longer fatigue life than WEL surface. It is worth noting that the fatigue lives at P1 and P8 are about 10 times longer than at P2. The predicted results are consistent with previous rolling contact fatigue  experimental testing [35], where a lot of cracks were always found at the WEL surface. Furthermore, the metallurgical examination of ex-serviced rail indicated that cracks were mainly initiated from the leading edge of the WEL [19]. Very similar results have also been observed when the rail material changes into a HH steel. The fatigue life was affected by many parameters, including the mechanical properties of rail materials, fatigue response of rail materials, and stress state at wheel and rail interface, etc. It is interesting to report that an increase in rail strength leads to a decrease in fatigue life. Such a phenomenon can be explained by Eq. (7), which reveals that the fatigue life of a rail has an inverse relationship with its yield strength. The influence of steel grade on rail fatigue life has been observed and reported in the railway network. In Ref. [54], the occurrence of rail squats, which are typical rolling contact fatigue defects, were statistically compared between the different rail grades (HH and SC). A tendency was found that there was more squat development in HH rail than in SC rail, implying that a high strength rail would contribute to earlier fatigue initiation and decrease of fatigue life. This field observation is consistent with the current study. In addition to the steel grade, the formation of very high strength WEL on the rail surface also contributes to the reduction of fatigue life as shown in Fig. 11.
FP has also been calculated as shown in Fig. 11(b) according to Eq. (6). As can be seen from comparison with Fig. 11(a), a smaller magnitude of FP indicates a larger fatigue life. For a SC steel rail, P2 has the maximum FP about 0.043 while P8 has the minimum FP about It's 0.006. The magnitude of FP at rail matrix is about 0.017. A slightly larger FP has been found at all positions for a rail with a higher strength.
In addition to the rail material property, influence of loading pressure, μ, a/b, and Sr on the FP around the WEL has also been investigated as shown in Fig. 12, respectively. There are two features in Fig. 12(a). First, FPs at WEL surface increase very quickly when loading pressure varies from 0.8 to 1.2 GPa, and then increases much slowly with further increasing loading pressure. FPs at WEL surface under 1.8 GPa are about twice of those under 0.8 GPa, which indicates that increasing loads significantly reduces the fatigue life of the WEL. Second, there is a much smaller influence at WEL-substrate interface and rail matrix than at WEL surface, and FPs only increase slightly with loading pressure. Figure 12(b) shows the variation of FP as a function of μ. The results reveal much larger FPs at WEL surface than the other positions, but all positions have the similar tendency in FP when μ changes. It means that poor lubricated condition promotes the fatigue failure of both the WEL and rail matrix, which should be avoided in practical service. Figure 12(c) shows the FP evolutions as a function of a/b. FP decreases quickly first when a/b is less than 5 in the WEL except position P8. Then, influence of the WEL geometry can be neglected when a/b is larger than 5 and very limited increase in FP can be observed. We can also see that FPs at positions P1 and P8 keep nearly constant for all cases considered. Based on these results, a critical WEL geometry, namely a/b = 5,   Table 2. The simulations also suggest that the a/b has a complicated influence on the fatigue lives of both rail matrix and WEL. Only considering the length or thickness of the WEL is not adequate to describing its geometry influence.
In Fig. 12(d), Sr leads to gradual increase of FP at all studied positions. It is also visible that fatigue performance at WEL surface is more sensitive to Sr than at rail matrix and WEL-substrate. Wang et al. [39] reported a significant drop of fatigue life when Sr increased from 0 to 0.3%, then gradual increase of fatigue life when Sr increased from 0.3% to 1%, and then very slight drop of fatigue life again when Sr increased from 1% to 10% based on a rolling contact fatigue testing of a rail material under wet conditions. However, rail matrix (P1) in Fig. 12(d) does not show the similar phenomena to their observations. There are two main reasons causing the difference. The first one is the wheel/rail geometry and existence of the WEL, which affect the contact stress/strain distributions significantly. The second reason is the frictional condition and testing operating parameters. μ changes significantly with Sr in Ref. [39], but the constant μ is applied in this study. Large influence of Sr on the μ has also been observed by Ma et al. [38]. Comparison between Refs. [38,39] suggests that it is very difficult to accurately focus on only one influence factor in the experiments, but such a problem can be solved with the help of finite element simulation and a better understanding of each particular influence factor can be obtained by combining the experiments and simulations. As can be seen, wheel-rail rolling contact fatigue behavior is very complicated and there are many influences, such as existence of the WEL, rail material property, loading pressure, μ, Sr, and so on. Among all these influence factors, loading pressure and the a/b play the dominant role. In next study, influence of the WEL on crack initiation and development, as well as its microstructural features during wheel-rail contact will be specifically designed and investigated.

Conclusions
A systematical study on influences of rail material, loading pressure, μ, a/b, and Sr during wheel-rail rolling contact was carried out considering the practical wheel-rail contact and geometry. The following conclusions can be drawn: 1) There is a sharp rise in normal contact stress when the wheel passes from rail matrix to the WEL. The difference of normal contact stress between rail matrix and WEL in steady contact condition increased gradually with both μ and loading pressure, but was not affected by the ratio of contact path to the WEL length.
2) It has been found that the residual stress along rolling direction was dominant and it was in compression with a maximum magnitude up to 600-700 MPa. There was much higher residual stress in WEL than in rail matrix.
3) Investigation of stress gradient from the WEL surface revealed three distinguished regions within a depth of 1 mm, which was consistent with the reported microstructure evolutions. The results indicated that the nanostructures in region-I had the largest residual stress, which dropped gradually in the mixed microstructure region (or region-II). There was no obvious change in region-III where the microstructure was characterized as the pearlite lamellae. Influence of the pre-existing residual stresses within the WEL on its rolling contact behaviour has not been considered in the current model, but they will be studied using a more comprehensive model in the future. 4) Loading pressure, μ, and a/b had a large influence on stress distributions in both rail matrix and WEL around the leading edge of the WEL. By contrast, Sr only had a minor influence. The stress difference between the WEL and rail matrix was found to increase significantly with loading pressure and a/b, decrease gradually with μ, and keep almost constant with Sr. 5) Fatigue analysis suggested that WEL surface had much lower fatigue life compared to the rail matrix and WEL-substrate interface under all the investigated conditions. It has also been found that fatigue life of the WEL could be reduced by increasing the rail strength. In addition, a critical WEL geometry with a/b = 5 (or α = 22.6°) was determined. Fatigue life of the WEL increased with a/b when it was less than 5, while slightly decreased with a/b when it was over 5.
In the current study, only thermal induced WEL has been investigated based on a 2D finite element model without considering the initial residual stress and the existence of transition region between WEL and matrix. Their influence on the rolling contact behavior of the WELs needs to be investigated using a more accurate 3D model in the future.