Residual fatigue life evaluation of rail at squats seeds using 3D explicit finite
element analysis
Xiangyun Deng*, Maysam Naeimi, Zili Li, Zhiwei Qian, Rolf Dollevoet Delft University of Technology, Delft, The Netherlands
Abstract: A modeling procedure to predict the residual fatigue life of rail at squats seeds is developed
in this article. Two models are involved: a 3D explicit Finite Element (FE) model to compute the stress and strain at squats in rail, and the J-S fatigue damage model to determine the residual fatigue life on the basis of the computed stress and strain. In the FE model dynamic effects of wheel-rail system under rolling contact is taken into account. Bilinear isotropic elastic-plastic material properties are adopted to represent the hardening of wheel and rail. Squats are subject to multiple loading cycles. The geometry of the squat is varied in the simulation corresponding to a growing squat at different ages. It is found that small squats lead to fatigue failure while severe ones lead to ratcheting failure.
Keywords:
1 Introduction
! " ! # $ %&' ( ) ! ( ( ( ( %*' + ( , + ( %- .' $ / " / # " # ( / 0 "*1&& *1&*#%2 3' / , ! *PhD Candidate. Section of Railway Engineering, Faculty of Civil Engineering and Geosciences, Delft
University of Technology. Stevinweg 1, 2628 CN, Delft, The Netherlands. X.Deng@tudelft.nl
Full Bibliographic Reference:
X. Deng, M. Naeimi, Z. Li, Z. Qian, R. Dollevoet, Residual fatigue life evaluation of rail at squats seeds using 3D explicit finite element analysis, Proceedings of the International Conference on Ageing of Materials & Structures, Delft, The Netherlands, May 2014.
The methods of fatigue analysis in this study are not new and they are previously employed by other researchers for different structures such as pressure vessels, gas turbines, aircraft, vehicles and components such as rotating disks, axles and crankshafts and even railway track components. However, the application of these techniques together with the solution of frictional rolling contact with a 3D finite element method and dynamic simulation is the novel.
Prediction of the fatigue life of rail material under the impact loading of the squats seeds in rail material is another pioneering adventure of this study. Based on numerical simulations and field observations, it is shown by Li et al. [7] that squat development is closely related to the dynamic contact force which is excited by the squat’s seeds and which is determined by the local eigen characteristics of the vehicle-track system. That work gives an appropriate understanding about the squat growth process and the wave pattern that often follows squats. The corresponding estimation of rail life time regarding various stage of squat growth has not yet been studied. This study is an attempt to acquire an estimative judgment about the effect of squats seeds on rail fatigue life behavior.
2
Fatigue analysis of wheel-rail material
In rolling contact fatigue, the causes of cracks initiation can be either ratcheting or low cycle fatigue[8]. To investigate the fatigue mechanism, it is necessary to recognize the loading path in rolling contact. In the wheel-rail rolling contact problem, the rail is subjected to a non-proportional multiaxial stress state, which results in the variation of the principal stress and maximum shear stress-strain directions during a passage of the wheel [9]. Therefore, to predict the life of rail, it calls for a multiaxial stress criterion including non-proportional loading. In this study, a model based on the energy density and a critical plane approach is used to predict the life of rail at squats. In addition, a well-recognized criterion for ratcheting failure in rail material proposed by Kapoor is used to predict the life of ratcheting.
Jiang and Sehitoglu [10, 11] proposed a multiaxial low cycle fatigue criterion for RCF phenomenon based on critical plane approach. In this criterion, it is postulated that both normal and shear components of stress and strain, on the critical plane, contribute to the damage of the material. The model is expressed as following equation:
=∆ 〈 〉 + . ∆ . ∆ (1)
where ∆ is the normal strain range, is the maximum normal stress, ∆ is the shear strain range, ∆ is the shear stress range, J is a material-dependent constant and 〈 〉 denotes the McCauley bracket 〈 〉 = (| | + )/2.
All the stress and strain quantities in Eq. 1 are on the critical plane where the fatigue parameter FP is the maximum (FPmax). Through a tensor rotation for the stress and strain, the
maximum FP and the critical plane are determined by surveying all the possible planes at a material point. The first term in Eq.1 considers the mean stress effect. The proposed multiaxial fatigue model has the correct form to capture the synergism between the shear and normal stress components. The relationship between fatigue parameter and crack initiation life is described by the following equation:
( − ) . = (2)
where Nf is the crack initiation life corresponding to fatigue parameter FP, FP0, m and C are
material fatigue properties obtained by best fitting base line experimental data. Fatigue damage is assumed to accumulate linearly. When the fatigue parameter FP is equal to or smaller than FP0, no fatigue damage is predicted.
According to Kapoor (1994) [12], if a material displays a constant ratcheting rate, the ratcheting rate parameter is considered constant in the present work as well. In [12], the equivalent ratcheting plastic strain per cycle can be used as the ratcheting strain, it is expressed as:
∆ = (∆ ) + (∆ /√3) (3)
where ∆ and ∆ are the average ratcheting normal and shear strain per cycle in three dimensional stress states. The equivalent ratcheting strain is calculated on the plane with the largest shear strain accumulation.
The ratcheting life can be estimated by
= "
∆ # (4)
where $ is the critical strain for failure by ratcheting.
The above formulations allow identifying the parameters that control the failure rates. It is suggested that the failure mechanisms of fatigue and ratcheting were independent and competitive so that the life of the component was governed by whichever would be expected to cause failure in the shorter number of cycles. If a material has a low ductility, low cycle fatigue failure could occur in finite number of cycles in the form of crack initiation and propagation. For a ductile material, extensive ratcheting strain can cause the extrusion of thin slivers at the surface. In either case, the combination of plastic ratcheting and fatigue can be used as a measure of failure prediction under these conditions [10, 11].
The rail steel grade that is generally being used in the Netherlands is the typical pearlitic steel, normal grade 900A. The material constants which are used in the J-S failure model can be generated by normal fatigue tests like tension-compression and torsion fatigue tests. In the current study, however, these constants are obtained from [13] for the equivalent pearlitic steel as: J= 0.32, m=2.50, C=1500000, a constant value of $ = 7.1 is assumed for the critical strain of rail material as it is proposed in [14], the value of FP0 is assumed to be 0.2 .
3
Numerical simulations
3.1 FE modelling of vehicle-track system
A three-dimensional finite element model of vehicle-track system with the implementation of elastic-plastic material and transient dynamic simulation was used to study the states of stresses and strains in rail material as it is schematically shown in Figure 1. The wheel is at its starting position. The wheel rotates with an initial speed to reach a distance which is long enough to relax the dynamic effects that are induced by the sudden loading of the wheel.
The current model of vehicle-track interaction system is previously employed by Zhao et al. (2012) [4] to study wheel-rail impact and the dynamic forces at discrete supports of rails. That numerical model is further developed in the current study to understand the detailed mechanism of rail fatigue and to obtain the RCF crack initiation life in rail material. The 3D FE model has been validated by Zhao and Li (2011) [5] in the normal and the tangential contact solutions against Hertz theory and Kalker’s program CONTACT. Applying some field observations with measurement of axle box acceleration, the model on the vehicle-track interaction has also been validated by field monitoring tests [7] when it was applied to squats.
The relevant structures of the vehicle and the track are both considered in the model together with their actual sizes. The primary suspension of the vehicle is considered in the simulation. In order to simulate the high-frequency dynamic behavior of wheel-set and track, the detailed flexibilities of the vehicle system is disregarded and the sprung mass of the vehicle is lumped into Mc that is connected to the wheel set through the primary suspension Kc and Cc (see Figure 1
and 2). The track system is a model of typical ballasted railway track, in which the supports of the rail are composed of the fastenings, the sleepers, and the ballast. Only a half wheel set and a half straight track are modeled in view of the symmetry of the system.
Figure 1 Schematic geometry of the vehicle-track FE model in longitudinal (left) and lateral (right) directions.
Figure 2 Left: 3D representation of the vehicle –track FE model. Right: Finest mesh pattern in the solution zone.
The FE model (Fig.2) is analyzed using explicit time integration scheme. The wheel–rail interaction is calculated using a surface-to-surface contact algorithm. The parameters of the FE model are listed in Table 1, in which the variables for defining nonlinearity of the wheel-rail materials are also given. In this work, the bilinear isotropic hardening model is employed for the materials of wheel and the rail. The stiffness and damping parameters in this table are taken from [4].
A constant friction coefficient of 0.6 is assumed for the simulations. The tangential loading of the contact is modeled by applying a traction force corresponding to the traction coefficient. Typical velocity of the trains in this work is considered 140 km/h based on the maximum operating speed in the Dutch railway network. Two wheel passages are simulated in the present work.
Table 1 The values of parameters used in the numerical simulations
Parameters (unit) Values Parameters (unit) Values
Static wheel load, Mc (kN) 116.8 Stiffness of primary suspension, Kc (kN/m) 880
Wheel weight (kg) 900 Damping of primary suspension, Cc (N.s/m) 4000
Rail weight per length(Kg/m) 54.42 Young’s modulus of wheel-rail material, Er (GP) 210
Sleeper mass Ms (kg) 280 Poisson’s ratio of wheel-rail material, νr 0.3
Friction coefficient 0.35 Density of wheel-rail material, ρr (kg/m3) 7800
Traction coefficient 0.15 Yield stress- work hardened rail (GP) 1.12 Rolling speed (km/h) 140 Yield stress- work hardened welds (GP) 0.99 Stiffness of ballast, Kb (kN/m) 45000 Tangent modulus of elastic-plastic rail (GP) 21
Damping of ballast, Cb (N.s/m) 32000 Young’s modulus of concrete material, Ec (GP) 38.4
Stiffness of rail pad, Kp (kN/m) 1300000 Poisson’s ratio of concrete sleeper material, νc 0.2
Damping of rail pad, Cp (N.s/m) 45000 Density of sleeper material, ρc (kg/m3) 2520
K c , Cc M c 600mm K p , Cp K b , Cb M s Rail Wheel 300mm 700mm Fastening Sleeper Ballast Solution zone, location of squat
Initial wheel location End of simulation 450mm K c , Cc M 1435/2 mm K p , Cp K b , Cb M s Rail Wheelset Fastening Sleepe Ballast S y m m e tr ic a l b o u n d a ry c o n d it io n
3.2 Modelling of squats on the rail surface
Squats are some kinds of singular rail surface defects which occur on the rail top and cause excitation in wheel–rail impact condition leading to large dynamic forces. Due to the large dynamic impact conditions, the fatigue life of material could be affected as a result of squat appearance in rail. Some examples of the recent studies about the influence of track geometry irregularities on rolling contact fatigue can be found in [4, 6, 15].
To understand the effects of dynamic forces excited by squats and find possible estimation of fatigue life of rail material, squats are modeled in the FE simulation as the geometric irregularities in the rail top. According to the observations and measurements in [7], a squat seed in its early stage is a small irregularity in the rail top surface without cracks, as shown in Figure 3(a). Such an irregularity due to, e.g. indentation, burn initially a V-shape in its longitudinal-vertical profile, and gradually develops into a W-shape. Cracks develop in later stages, usually first in the middle ridge of the W-shape, as shown in Figure 3(b). Based on these observations, four idealized defects are studied(scenarios1, 2, 3 and 4) in the present work as the geometric models of the squats either in V or W shapes at different stages, besides a smooth rail as the reference situation (scenario 0). Figure 4 shows the simulated defect profiles (including two V-shape and two W-shape defects) in the longitudinal-vertical plane and one 3D representation.
a b
Figure 3 Squats in different stages, (a) a V shape in early stage, (b) a severe W shape in later stage.
Figure 4 All simulated scenarios of rail squats. (a) Longitudinal–vertical profile at the middle of the running
band. (b) 3D representation of one W-shape defect sample.
3.3 Strategy of Calculation for the residual life
First of all, the material properties need to be appropriately considered to anticipate the fatigue life. The purpose of the present work is to predict the residual life of the rail under operation. Therefore, the object of investigation is a hardened rail material rather than a new rail. In the present work, hence, the value of 990MPa is used as the yield stress of the hardened rail material. This yield stress is induced from hardness of rail steel which was measured in rails under operation. The fatigue parameters for the new rail steel can be found in literature, ([13] for instance). Using these parameters, the whole life can be predicted. The residual life can be obtained by the whole life subtracting the current operation time which can be measured by
-60 -40 -20 0 20 40 60 80 100 120 -0.20 -0.15 -0.10 -0.05 0.00 0.05
S0-Basic scenario(smooth rail) S1-Defect_V1 S2-Defect_V2 S3-Defect_W1 S4-Defect_W2 Ve rt ic a l h e ig h t (m m ) Longitudinal position (mm) Rolling direction from left to right
0.56 0.58 0.60 0.62 0.64 0.1584 0.1588 0.1592 -0.010 -0.005 0.000 0.005 0.010 0.015 V e rt ic a l, z ( m ) Late ral, y (m ) longitudin al, x(m) Rolling direction
considering the relationship between yield stress and operation time in the field or the one obtained by plastic law. In the present work, however, only the first step is considered.
Secondly, the state of stress and strain at squats is changed with the growth of squats, which can excite impact force. It is assumed in the present work that the stress and strain are constant in each cycle for each defect model. The variation of the stress and strain is achieved by modeling different geometry of defects as squats representative for different stages.
Lastly, for the position of fatigue, according to the observations and detections, usually, squats cracks initiate in the middle ridge of the W-shape at later stages. The point of cracks initiation is called critical point in the present work. Since stress level influence the cracks initiation, the critical point usually has the maximum von Mises stress in one wheel passage at W-shape stages. In this work, therefore, the critical point in two W-shape defects is investigated to predict the rest life of rail. To predict the residual life of squats at different stages and the influence of growth of squats on the residual life, the same point in V shape defects and smooth rail are considered to predict the life. The material characteristics for all scenarios of defects are assumed to be the same.
4
Results of simulations and fatigue analysis
4.1 Stress-strain responses
Numerical analyses are performed for all the defects models (one case for smooth rail and four cases for models with squats seed) that are defined in the previous section. In Figure 5(a), it shows the maximum von Mises stress variation with time step in the second wheel passage for the same point which is the critical point determined in S3 and S4. The critical point is located at the middle ridge of the W- shape defect. This location is agreed with that of observations. It is probed then to see whether the plastic deformation occurs or not on the critical point at different stages of the squat’s growth process. In the present case, plastic deformation occurs in S3 and S4, while there is no plastic deformation in the rest three cases.
The gradual evolution of different stress components at the critical points of the rail at different time steps of a wheel passage is depicted in Figure 5(b). Hereby only the result of S4 scenario is demonstrated as a sample. This shows that the stress distribution diagrams are getting higher values when the wheel is approaching to the critical point. Also the stress alteration diagrams of the critical point are non-proportional and the material is exposed to a complex stress states.
Figure 5 Variation of stress in the element with critical condition (second cycle), (a) Variation of Von
Mises stress for all cases of defects, (b) Variation of stress components for scenario S4 as a sample. 4.2 Calculation of fatigue crack initiation life
Considering the stress and strain histories for the critical point, a tensor rotation technique is used to obtain the stresses and strains on an arbitrary material plane for the critical point in second wheel passage. According to the fatigue model, Eq. (1), the variations of fatigue parameter FP with respect to the angles '( and ' for scenario S1 and S4 as two samples are
720 740 760 780 800 0 200 400 600 800 1000 V o n M is e s s tr e s s ( M P a ) Time step S0 S1 S2 S3 S4 a 720 740 760 780 -2000 -1500 -1000 -500 0 500 S tr e s s e s ( M P a ) Time step σyy σzz σxx τyz τzx τyx b
calculated, as shown in Figure 6. '( and ' in this figure are the two characteristic angles of the normal vector of the material plane, while the third angle can be calculated by spherical law of cosines.
Figure 6 An illustration of the search of the critical plane for scenario S3 (a) and S4 (b).
By applying J-S fatigue damage model (Eq. 2), the fatigue crack initiation life Nf is predicted
on the possible crack plane determined above. Figure 7 shows the results of so called fatigue curves for crack initiation life variations with the possible critical plane for S1 and S4 scenarios. It can be apparently seen that the fatigue life parameter progressively increases with the reduction of fatigue parameter (FP). All of the individual spots in the fatigue curve are corresponding to the all possible critical plane at the critical point. The minimum value of Nf by
the way is the most critical case of fatigue crack initiation in material since it is equivalent maximum value of FP. The smallest Nf value is remarked in the abscissa to demonstrate the
number of cycles to crack initiation in rail material for two scenarios. Employing the same procedure, the fatigue life of rail in other scenarios is determined.
Figure 7 Results of fatigue parameter and fatigue life of material for two samples of FE simulations, (a) S1
scenario, (b) S4 scenario.
4.3 Calculation of ratcheting failure life
In the current work only two sequential cycles of wheel passage are simulated for the ratcheting analysis since the realistic material model of hardened rail under operation is considered. The ratcheting rate parameter is considered constant in the present work. Figure 8 shows the shear stress-strain and the normal stress-strain on the plane with the largest shear strain at each time step in two cycles for scenarios S3 and S4 in which the plastic deformation occurs at the critical point. 0 45 90 0 45 90 0 0.3 0.6 θ1 θ2 F P a 0 45 90 0 45 90 0 0.3 0.6 θ1 θ2 F P b 106 107 108 109 1010 1011 0.2 0.4 0.6 0.8 1 F a ti g u e p a ra m te r, F P
Number of cycles to fatigue, Nf
Nf Nf = 2796,028 S1 a 106 107 108 109 1010 1011 0.2 0.4 0.6 0.8 1 F a ti g u e p a ra m e te r, F P
Number of cycles to fatigue, Nf
Nf
Nf =1737,732 S4 b
Figure 8 Variation of shear strain and normal strain in the critical point for two scenarios, (a) scenario S3, (b) scenario S4.
From this figure, it can be found that the critical point experiences incremental shear strain and normal strain between first and second loading cycles. The accumulative shear strain and normal strain at the critical point is transferred to the equivalent plastic strain using Eq.3. These values are listed in Table.2. Using these incremental values and employing Eq. 4, the number of cycles to crack initiation by ratcheting (Nr) is estimated as shown in Table 2. Similar to the
normal fatigue life, the minimum value of Nr in the whole simulation stands for the ratcheting
life of rail material. It is worth noting that the ratcheting life Nr is unlimited for smooth rail and
both S1 and S2 defect scenarios since no plastic deformation occurs in these cases.
Table 2 Calculation of ratcheting failure life for all cases
Defect models Max. plastic normal strain Max. plastic shear strain Nr
S0- smooth rail 0 0 ∞
S1- defect V1 0 0 ∞
S2- defect V2 0 0 ∞
S3- defect W1 3.0E-05 2.0E-05 207,346
S4- defect W2 5.9E-05 2.2E-05 117,643
4.4 Discussion, dominant fatigue mechanism
The values of Nf and Nr are calculated for all scenarios using the same procedure. The results are
summarized in Figure 9 which shows the crack initiation life of material governed either by fatigue or ratcheting failure in all scenarios. For defect S0, S1 and S2, there are no plastic deformation occurring, therefore the ratcheting strain is mathematically zero and the ratcheting life is unlimited. In these three cases therefore, the crack initiation life of rail material is dominantly affected by fatigue. Comparing the required number of cycles for fatigue in all scenarios together with ratcheting life in S3 and S4, it can be found that the crack initiation life of material is dominantly affected by the presence of squat as the surface defects. Nf is reduced
once the larger size squat geometry is applied on rail surface. Nr for S3 and S4 is moreover
reduced with the growth of squat. The overall rate of reduction in the number of cycles has been more drastic for the fatigue life especially for transition between scenarios S3 and S4.
In contrast with the ratcheting life, a gentle decrease in fatigue life of material is drawn for S3 and S4 scenarios. Comparing the results of Nf and Nr in the two scenarios, it is worth noticing
that the value of Nf has been relatively higher than Nr. The amount of discrepancy between the
fatigue life and ratcheting failure for these two cases is the relative life (Nf/Nr) of 11 and 9.84,
respectively. In S3 and S4 cases, therefore, the crack was initiated by ratcheting mechanism rather than fatigue one.
720 750 780 -0.003 -0.002 -0.001 0.000 0.001 0.002 0.003 First cycle Second cycle S tr a in Time step First cycle Second cycle Maximum shear strain
Normal strain a 720 750 780 -0.003 -0.002 -0.001 0.000 0.001 0.002 0.003 Normal strain: First cycle Second cycle S tr a in Time step First cycle Second cycle Maximum shear strain: b
0 1.000.000 2.000.000 3.000.000 4.000.000 S0 N u m b e r o f c y c le s t o f a ti g u e Defect cases Nf Nr S1 S2 S3 S4 2814660 2796030 2617010 2287000 1737730 207,346 176,759
Figure 9 Results of Nf and Nr for all defect scenarios.
5
Concluding remarks
A dynamic finite element model was employed in this study to simulate the dynamic interaction of vehicle-track system at squats. Based on field observation and measurement, a series defects with different geometry representing a squat at different stages are modeled to predict the rest life of rail at different stages of a squat. Jiang’s fatigue model and Kapoor’s ratcheting model are then employed for prediction of fatigue life Nf and ratcheting failure life Nr in material. The
approach in this study therefore considered the initiation of cracks under combined ratcheting plastic strain and multiaxial fatigue. The results of Nf and Nr were summarized and compared for
different cases. Based on these results, following conclusions can be made.
- For smooth rail and small squats seeds, crack initiation is mainly caused by fatigue while ratcheting phenomenon does not exist. When squats seeds become typical W-shape, stress increases in their middle ridge and ratcheting strain is found.
- Predicted lives of RCF crack initiation and ratcheting failure for the scenario S3 were observed around 5.71×105 and 2.07×105 cycles of wheel passage, respectively, followed
by 4.34×105 and 2.07×105 for scenario S4. The ratcheting life for these two cases can be
theretofore considered as the life to crack initiation.
- The fatigue life of the rail material is significantly decreased by the presence of severe W-shape squats on rail surface. The main cause for this is ratcheting.
Acknowledgement
This research is supported by the Dutch Technology Foundation STW, which is part of the Netherlands Organization for Scientific Research (NWO). ProRrail is kindly acknowledged for providing part of the funding and technical supports.
References
[1]. A. Fatemi and D.F. Socie, A Critical Plane Approach to Multiaxial Fatigue Damage
Including out‐of‐Phase Loading. Fatigue & Fracture of Engineering Materials & Structures,
1988. 11(3): p. 149-165.
[2]. K. VAN and M. Maitournam, Rolling contact in railways: modelling, simulation and damage
prediction. Fatigue & Fracture of Engineering Materials & Structures, 2003. 26(10): p.
939-948.
[3]. K. Karttunen, E. Kabo, and A. Ekberg, A numerical study of the influence of lateral
Institution of Mechanical Engineers, Part F: Journal of Rail and Rapid Transit, 2012.
226(6): p. 575-586.
[4]. X. Zhao, Z. Li, and J. Liu, Wheel–rail impact and the dynamic forces at discrete supports of
rails in the presence of singular rail surface defects. Proceedings of the Institution of
Mechanical Engineers, Part F: Journal of Rail and Rapid Transit, 2012. 226(2): p. 124-139.
[5]. X. Zhao and Z. Li, The solution of frictional wheel‐rail rolling contact with a 3D transient
finite element model: Validation and error analysis. Wear, 2011. 271(1-2): p. 444-452.
[6]. X. Zhao, Z. Li, and R. Dollevoet, The vertical and the longitudinal dynamic responses of the
vehicle–track system to squat‐type short wavelength irregularity. Vehicle System
Dynamics, 2013. 51(12): p. 1918-1937.
[7]. Z. Li, R. Dollevoet, M. Molodova, and X. Zhao, Squat growth—Some observations and the
validation of numerical predictions. Wear, 2011. 271(1): p. 148-157.
[8]. Y. Jiang and H. Sehitoglu, Rolling contact stress analysis with the application of a new
plasticity model. Wear, 1996. 191(1): p. 35-44.
[9]. J.W. Ringsberg, Life prediction of rolling contact fatigue crack initiation. International Journal of fatigue, 2001. 23(7): p. 575-586.
[10]. Y. Jiang, A fatigue criterion for general multiaxial loading. Fatigue and fracture of engineering materials and structures, 2000. 23(1): p. 19-32.
[11]. Y. Jiang and H. Sehitoglu, A model for rolling contact failure. Wear, 1999. 224(1): p. 38-49. [12]. A. Kapoor, A re‐evaluation of the life to rupture of ductile metals by cyclic plastic strain.
Fatigue & fracture of engineering materials & structures, 1994. 17(2): p. 201-219.
[13]. O. Onal, D. Canadinc, H. Sehitoglu, et al., Investigation of rolling contact crack initiation in
bainitic and pearlitic rail steels. Fatigue & Fracture of Engineering Materials & Structures,
2012. 35(11): p. 985-997.
[14]. M. Akama, H. Matsuda, H. Doi, and M. Tsujie, Fatigue crack initiation life prediction of rails
using theory of critical distance and critical plane approach. Journal of Computational
Science and Technology, 2012. 6: p. 54-69.
[15]. T.S. Deng, X. Zhao, B. Wu, et al., Prediction of crack initiation of rail rolling contact fatigue. Applied Mechanics and Materials, 2013. 344: p. 75-82.
