Analysis and modeling of unsteady effects within the wheel-rail contact

Sergio Diaz de la Fuente [ D i r e c c i ó n d e l a c o m p a ñ í a ] ANALYSIS AND MODELING OF UNSTEADY EFFECTS WITHIN THE WHEEL-RAIL CONTACT PROJECT SPONSORS: Master’s thesis to obtain the academic degree of Master in Mechanical Engineering Supervisor: Dietmaier Peter, Ao.Univ.-Prof. Dipl.-Ing. Dr.techn. External Advisors: Six Klaus, Dr. –Meierhofer Alexander, DI Graz University of Technology, 2012 Analysis and modeling of unsteady effects within the wheel-rail contact Master’s thesis 2 ACKNOWLEDGMENTS This Master’s Thesis is intended for mechanical engineers, physics, and, generally for researchers on the field of railways and contact mechanics. I wrote under the supervision of Professor Dr. P. Dietmaier of Graz University of Technology (Austria), Dr.techn. K. Six and DI A. Meierhofer, both of the Virtual Vehicle Research and Test Center (ViF). This thesis is part of the research project “D03T06DynKraS-Dynamic-Wheel-Rail Contact Forces” of the Virtual Vehicle Research and Test Center (ViF) carried out in cooperation with Siemens AG Österreich, voestalpine Schienen GmbH, LBFoster (Portec Rail/Kelsan Technologies), ÖBB Infrastruktur AG and the Institute of Applied Mechanics at Graz University of Technology. The financial support of “COMET K2 Forschungsförderungs-Programm” of the Austrian Federal Ministry of Transport, Innovation and Technology (BMVIT), the Federal Ministry of Economy, Family and Youth (BMWFJ), the Austrian Research Promotion Agency (FFG), the State of Styria and the Styrian Business Promotion Agency (SFG) is greatly appreciated. These lines are the result of half a year research at ViF, during my stay in Graz, doing an Erasmus Exchange Programme. When Dr.techn. K. Six sent me the project offer, I did not doubt to accept it because a Railways project was what I was looking for before going abroad. Since I started my degree in Barcelona, I was interested in mechanics and especially in Railways maybe because I think they are the present and the future of transport. Living and working in Austria was really exciting because it was my first experience abroad and although the beginning was hard because the theory was new to me, I learnt a lot of concepts these months and, above all, my mind changed. For these and more reasons, my first word of thanks goes to the person who took care of me at the company those months (and also later), A. Meierhofer. He helped me when I had problems with the thesis and he showed a great deal of patience. I think I still owe him a beer. The second word of thanks goes to K. Six and Prof. Dietmaier for giving me the chance to work there. Also many thanks to all the ViF team because it was a pleasure to work with all of them. Thanks to my home university, Universitat Politecnica de Catalunya, and to the Erasmus Programme, for allowing me to study abroad. I especially would like to say thank you to my parents for their financial and motivational support since I was a child. Also thank you to my family and friends. And finally, many thanks to all the friends I met these months in Austria, I also learnt a lot from each one of them. S. Diaz Barcelona, 15 September 2012 Analysis and modeling of unsteady effects within the wheel-rail contact Master’s thesis 3 ABSTRACT The traction characteristic of the wheel-rail contact is of high importance, e.g., for vehicle dynamics, and for the transmission of forces during traction and braking. This characteristic is influenced by fast changes of mechanical parameters, e.g., creepages, normal loads and coefficient of friction during the movement of the train. CONTACT, which is based on Kalker’s Theory, is able to describe these unsteady effects for low creepages (there are adhesion and slip at the same time) but it is unsuitable for Multi Body System (MBS) simulations due to its high computational effort. Therefore, a linearised model, called Unsteady Brush Model (UBM), is presented in this thesis that can be used in combination with MBS software. The results of the UBM agree qualitatively with the results of CONTACT. Analysis and modeling of unsteady effects within the wheel-rail contact Master’s thesis 4 CONTENTS ACKNOWLEDGMENTS ........................................................................................ 2 ABSTRACT .............................................................................................................. 3 LIST OF SYMBOLS ................................................................................................ 6 1 INTRODUCTION ......................................................................................... 10 2 THE CONTACT PROBLEM ....................................................................... 12 2.1 The half space ............................................................................................. 15 2.2 The normal contact problem (Hertz theory) ........................................... 16 2.3 The tangential contact problem ................................................................ 23 2.3.1 Definition of Creepages ........................................................................ 25 2.3.2 Traction bound ...................................................................................... 27 2.3.3 The traction coefficient ......................................................................... 30 2.3.4 Area of Stick and area of slip ................................................................ 31 2.3.5 Kalker’s Exact Theory: CONTACT ..................................................... 32 2.3.6 Linear theory: FASTSIM ...................................................................... 33 3 THE UNSTEADY BRUSH MODEL ........................................................... 37 3.1 Derivation for a moving coordinate system ............................................. 37 3.1.1 Constraints ............................................................................................ 40 3.2 The Unsteady Brush Algorithm (UBA). ................................................... 42 4 RESULTS ....................................................................................................... 46 Analysis and modeling of unsteady effects within the wheel-rail contact Master’s thesis 5 4.1 Parameters of the study ............................................................................. 46 4.1.1 Geometry .............................................................................................. 46 4.1.2 Cinematic and dynamic parameters ...................................................... 49 4.1.2.2 Creepage as a sine function (sine mode) ................................................................. 52 4.2 Results of the unsteady effects from CONTACT .................................... 54 4.2.1 Influence of the frequency on the traction coefficient in sine mode..... 54 4.2.2 Coefficient of friction ........................................................................... 60 4.2.3 Normal load .......................................................................................... 62 4.2.4 Contact patch’s geometry ..................................................................... 69 4.2.5 Average value of the sine function ....................................................... 71 4.3 UBA vs. FASTSIM. Steady solutions ....................................................... 75 4.3.1 Comparison of the results ..................................................................... 76 4.4 UBA vs. CONTACT. Unsteady solutions ................................................. 78 4.4.1 Comparison of the results ..................................................................... 78 5 SUMMARY&CONCLUSIONS ................................................................... 81 APPENDIX A .......................................................................................................... 83 APPENDIX B .......................................................................................................... 84 LIST OF FIGURES ................................................................................................ 86 REFERENCES ....................................................................................................... 88 Analysis and modeling of unsteady effects within the wheel-rail contact Master’s thesis 6 LIST OF SYMBOLS Formulas and units A a … semi-axe of the contact ellipse in x direction [m] … Yaw angle [-] B b … semi-axe of the contact ellipse in y direction [m] C cx … longitudinal creepage [-] cy … lateral creepage [-] cs … spin creepage [m-1] C0 … average value of the creepage in sine mode [-] Ci … constant of Hertz’s surface [-] Cp … amplitude of the creepage in sine mode [-] D d … distance between two particles [m] DQ … displacement of each particle in one time step [m] d(x,y) … vertical distance between two surfaces [m] … sum of displacements [m] … displacement of Pi [m] … time delay [ms] … Dirac Delta function [-] Δt … time step [s] Δx … discretization in longitudinal direction [m] Δy … discretization in lateral direction [m] Analysis and modeling of unsteady effects within the wheel-rail contact Master’s thesis 7 E e … eccentricity of the ellipse [-] E … combined Young’s modulus [MPa] E(e) … complete elliptic integral [-] … error [-] F f … Traction coefficient [-] F … Frequency [Hz] FN … total normal load [kN] … angle of the tangential stress components [-] G G … combined shear modulus [MPa] Gi … shear modulus of the body i [MPa] … contact angle [-] H H … Heaviside step function [-] K K … difference parameter [-] Ki … flexibility parameters [N/mm3] K(e) … complete elliptic integral [-] L L … flexibilities matrix [mm3/N] Li … FASTSIM flexibilities [mm3/N] M mx … number of elements in x direction [-] my … number of elements in y direction [-] Analysis and modeling of unsteady effects within the wheel-rail contact Master’s thesis 8 … coefficient of friction [-] N … combined Poisson’s ratio [-] … Poisson’s ratio of the body i [-] P pn … normal stress distribution [MPa] pn0 … maximum normal stress [MPa] … average pressure [MPa] R … rail radius in x direction [m] … rail radius in y direction [m] … wheel radius in x direction [m] … wheel radius in y direction [m] S S(x,y) … contact surface [m2] Sij … discreet contact surface [m2] T t … time [s] Ts … Period [s] Ti … total tangential force in i direction [kN] … tangential stress [MPa] U uij … displacement [m] V V … velocity [m/s] Analysis and modeling of unsteady effects within the wheel-rail contact Master’s thesis 9 Y yij … surface particle [-] … potential function of Boussinesq [-] Z Ωi … angular speed of body i [rad/s] Abbreviations MBS … Multibody System TB … traction bound UBA … Unsteady Brush Alogrithm UBM … Unsteady Brush Model Analysis and modeling of unsteady effects within the wheel-rail contact Master’s thesis 16 2.2 The normal contact problem (Hertz theory) In 1882, Hertz was the first one to solve the contact problem for an elastic ball on a flat rigid surface. It can be used for solving the normal contact problem in case of the wheel-rail contact (Johnson, 1987, pp. 84-106). Considering the elastic bodies 1 and 2, and mounted on axles (Figure 2.1), the Hertzian theory can be applied using the following assumptions: - Elastic material behavior - Semi-infinite spaces - Large curvature radii compared to the contact size. - Two bodies can be approximated with functions of second order. According to the assumptions, the two surfaces in contact can be described as follows: ( 2.2) Ai , Bi and Ci with i = 1,2 are constants near the contact point O and are related to the curvatures. Using a suitable choice of axes, the coefficients Ci with i = 1,2 can be zero, and the equation yields: , for i = 1,2. ( 2.3) Analysis and modeling of unsteady effects within the wheel-rail contact Master’s thesis 17 The curvatures of the bodies are defined by the four principal radii (Figure 2.4), so there will be a point O where the normal distance between them is minimal. Figure 2.4. Hertzian contact. For the wheel: ( 2.4) And for the rail: ( 2.5) Analysis and modeling of unsteady effects within the wheel-rail contact Master’s thesis 18 If the rail is considered straight, Rrx is infinite. Then, the vertical relative distance d(x,y) can be defined as follows: For the principal radii of the two bodies, A and B can be calculated as follows: ( 2.6) Before deformation, two surface particles yij are considered, with i = 1,2 and j = x,y,z. These particles satisfy Eq.( 2.2). Two distant points of the bodies (P1 and P2) move towards O parallel to z-axis (displacements ). By analogy, O moves from body 1, and from body 2 (see Figure 2.5). Then, both bodies contact each other, and due to the contact pressure, both surfaces are deformed. The elastic displacements relative to the distant points Pi are denoted by . This is also shown in Figure 2.5. Analysis and modeling of unsteady effects within the wheel-rail contact Master’s thesis 19 Figure 2.5. Graphic representation of the deformation of two bodies loaded. The particles yij are shown in red at the deformed state. If surface points yij are in contact then: ( 2.7) Writing , Eq. ( 2.7) and according to Eq. ( 2.6), the expression yields: Analysis and modeling of unsteady effects within the wheel-rail contact Master’s thesis 20 ( 2.8) where is the elastic displacement of the contact surface (sum of displacements of surface bodies 1 and 2). The contour of the contact area is assumed to be an ellipse with semi-axes a (x direction), and b (y direction). Hertz concluded the normal contact problem is analogous to the electro-static potential problem. By analogy, the normal stress distribution profile is given by ( 2.9) where pn is the normal stress distribution, and pn0 the maximum normal stress. Using the potential functions of Boussinesq for a general point in the body as a classical approach (Routh, 1908), the normal displacement of the surface is ( 2.10) with (2.11) Analysis and modeling of unsteady effects within the wheel-rail contact Master’s thesis 21 The surface displacement within the contact area satisfies the elastic displacements Eq. ( 2.8). Then, it can be written as follows ( 2.12) where M, L and are: ( 2.13) ( 2.14) ( 2.15) Where e = (1 - b2/a2)1/2, b<a is the eccentricity of the ellipse, and E, K are complete elliptic integrals ( 2.16): ( 2.16) To solve the integrals shown in Eq. ( 2.16) there are different numerical tables: K and E are tabulated in (Abramowitz & Stegun, 1964) with great precision and also in (Jahnke & Emde, 1943) with a precision of about 4 decimals. Analysis and modeling of unsteady effects within the wheel-rail contact Master’s thesis 22 The normal stress distribution is semi-ellipsoidal and, knowing the ellipsoid volume , the total normal load is: ( 2.17) And the average pressure is: ( 2.18) The size of the contact ellipse is defined by the semi-axes a and b, considering a>b the relation is: ( 2.19) Analysis and modeling of unsteady effects within the wheel-rail contact Master’s thesis 23 2.3 The tangential contact problem The tangential contact problem can be solved independently from the normal contact problem as mentioned in Subchapter 2.1. The contact forces for each wheel are: - FN for the normal force - Tx for the longitudinal creep force - Ty for the lateral creep force in the contact plane Figure 2.6. Wheelset geometry and creep forces. In a similar way to Figure 2.5, in the tangential problem there is also bodies’ deformation. There are normal displacements and also tangential displacements due to Analysis and modeling of unsteady effects within the wheel-rail contact Master’s thesis 24 the relative velocity between both bodies. These displacements and deformations are shown in Figure 2.7 : Figure 2.7. Description of the normal and tangential contact. The displacements of the particles are not only in zaxis. Analysis and modeling of unsteady effects within the wheel-rail contact Master’s thesis 25 About the tangential contact problem, Carter was the first to give an adequate solution to the force relative to the longitudinal creepage ( 2.20). Then, (Rocard, 1954) studied the relationship between the yaw angle and the guiding force, in the lateral direction ( 2.24). In the 1960s, Johnson and Kalker gave an expression of the creepage introducing variable coefficients depending on the g ratio of the contact ellipse. This expression is the most used nowadays. 2.3.1 Definition of Creepages Tangential contact problem is influenced by the relative velocity between the two bodies: the creepages. There are three different creepages which are defined as follows - Longitudinal Creepage ( : (dimensionless) ( 2.20) - Lateral Creepage ( : (dimensionless) ( 2.21) - Spin Creepage ( : (m-1) ( 2.22) V1 and V2 are the absolute velocities of the particles in bodies 1 and 2. The difference of the velocities is divided by the reference velocity Vref. In Eq. ( 2.22), are the angular velocities of the bodies, and are projected on the normal to the contact (Figure 2.8). Analysis and modeling of unsteady effects within the wheel-rail contact Master’s thesis 32 2.3.5 Kalker’s Exact Theory: CONTACT The first exact rolling code was published by (Kalker, 1967). To implement this theory, a linear programming method was created for the two-dimensional case. This method was created to treat cases of quasiidentical and non-quasiidentical cylinders with parallel axes (Kalker J. , 1971). This method is generalized to the threedimensional case (Kalker & Goedings, 1972) and implemented in a program called DUVOROL (Kalker J. , 1979). In DUVOROL only steady state rolling of quasiidentical bodies is considered. Later, (Kalker J. , 1983) implemented another program, based on (Fichera, 1964) and (Duvaut & Lions, 1972) methods. This program was called CONTACT and it is able to solve: - The normal contact problem. - The static shift problem in its incremental form. - Non-steady state rolling. - Steady state rolling. Also, this program is able to calculate local displacements, displacement gradients, linearised strains, and stresses. It is based on the generalization of Galin’s Theorem (Kalker J. , 1990, p. 94). Analysis and modeling of unsteady effects within the wheel-rail contact Master’s thesis 33 2.3.6 Linear theory: FASTSIM The expressions from Kalker’s exact theory can be linearised. (Kalker J. , 1990, pp. 103-135) To solve the tangential contact problem, first the simplified theory uses Hertz’s theory (Figure 2.13). The resulting normal stress distribution is elliptical which might lead to numerical problems, because of infinite gradients at the leading and trailing edges that make impossible to get full sliding. Then, the elliptical distribution is transformed into a parabolic one (Kalker J. , 1990) ( 2.32) ( 2.33) The parabolic pressure distribution is a reasonable approximation of the Hertz’s elliptical distribution when the ratio 2a/uz > 4 and Poisson’s ratio is <0.45. This distribution is used to solve the tangential contact problem. The two different normal stress distributions are shown in Figure 2.12 for the same normal loads: Analysis and modeling of unsteady effects within the wheel-rail contact Master’s thesis 34 Hertz’s considers deformation outside the contact zone (Kalker J. , 1990, pp. 107111). Instead, in the simplified theory, displacements outside the patch are zero (Figure 2.13.) Figure 2.13. Hertz’s Theory (a) and the Simplified theory (b) On the other hand, the tangential problem in the simplified theory can only be solved if there is quasiidentity, because then the normal displacements of both bodies are the same (Kalker J. , 1990, pp. 112-116). If there is not quasiidentity: uz1 ≠ uz2, Figure 2.12. Different distributions of the normal stresses. (a) Elliptical; (b) Parabolic A B Analysis and modeling of unsteady effects within the wheel-rail contact Master’s thesis 35 which means displacements need to be compensated by the normal traction, and this effect is not accounted in the simplified theory. Aside of this condition, the simplified theory uses the following assumptions: - The normal contact is assumed to be Hertzian - Normal stress distribution is parabolic - The surface is assumed to be brush like, i.e., made independent bristles - The stress calculation starts from the leading edge, element to element. The tangential stresses become from the linear theory (adhesion is assumed over the contact patch): (2.34) where ai is the distance from the y-axis to the element i, and yi is the distance from the x-axis to element i. The flexibilities of the bristles Li (i = 1,2,3) are given by (Kalker J. , 1990, pp. 123-125): (2.35) where c11, c22 and c23 are the Kalker’s Coefficients (Kalker J. , 1990, p. 125) , G is the Shear Modulus, a and b are the semi-axes of the ellipse. These are given as a function of b/a ratio of the contact ellipse, and of the Poisson’s ratio in Appendix A. But there is also a polynomial fit for a Poisson’s ratio of 0.27 (typical for the steels), from b/a = 1/25 to b/a = 25: Analysis and modeling of unsteady effects within the wheel-rail contact Master’s thesis 36 (2.36) Figure 2.14. Comparison of the tangential stress distribution of both theories. The results of the simplified theory are an approximation of the exact theory, and this fact is shown in Figure 2.14. As shown in Eq.(2.34), the tangential stresses are linear in x direction. The algorithm based on this simplified theory is called FASTSIM (Kalker J. , 1990), and it is not able to calculate unsteady effects. It was designed to calculate pressures and forces in a steady state. Despite that, FASTSIM is useful for calculating the solutions in a faster way, and the results obtained can be compared with the steady results from CONTACT. Analysis and modeling of unsteady effects within the wheel-rail contact Master’s thesis 37 3 THE UNSTEADY BRUSH MODEL While CONTACT gives the exact solution, the computational effort is very high. FASTSIM, on the other hand, gives approximate solutions with acceptable accuracy but only for steady state. Hence, a method is desired that is able to calculate Non-Steady solutions with the minimum time requirement possible: The Unsteady Brush Model (UBM). This method is not the only one developed for this purpose (Alonso & Giménez, 2008). 3.1 Derivation for a moving coordinate system The UBM solves a specific case of the Discrete Transport Equation of the matter in a continuum for the area of adhesion (see Subchapter 2.3.4). There, the creepage equals the strain. (Meierhofer, Hardwick, Lewis, Six, & Dietmaier, 2012). For a moving coordinate system, the principal equation is (3.1) Now, it is assumed that the tangential stress (x,y,t) is proportional to the tangential displacement (x,y,t): (3.2) Analysis and modeling of unsteady effects within the wheel-rail contact Master’s thesis 38 Here, L denotes the flexibility matrix (see Eq. (2.35)). Then, Eq.(3.2) becomes (3.3) If longitudinal, lateral and spin creepage are taken into account, Eq(3.3) can be rewritten as: (3.4) Analogue to Kalker’s derivation of the simplified theory (see chapter 2.3.6), the flexibilities Li are assumed to be different in each direction. The values for the flexibilities are taken from Kalker’s Simplified Theory (Appendix A). The inverse of the flexibility is the parameter Ki with i = x,y,s. (3.5) Eq. (3.4) yields: (3.6) Analysis and modeling of unsteady effects within the wheel-rail contact Master’s thesis 39 This equation and the contact area, assumed to be Hertzian, will now be discretized on a mx×my grid using Euler’s method with constant grid spacing Δx in the x direction and Δy in the y direction (Figure 3.1) Figure 3.1. Discretization of the contact ellipse in a 22x22 grid (3.7) Then Eq(3.6) yields: (3.8) Analysis and modeling of unsteady effects within the wheel-rail contact Master’s thesis 40 Isolating , it yields: (3.9) Where x denotes the direction, i,j the position on the grid and t the time. For the Y direction the equation is: (3.10) 3.1.1 Constraints The UBM expression is a Partial Differential Equation (PDE) and it needs a boundary condition and an initial condition to be solved. Both conditions have a physical sense. The contact patch is a closed surface defined by a function S(x,y) and discretized into mx×my elements. Analysis and modeling of unsteady effects within the wheel-rail contact Master’s thesis 41 Outside the contact patch there are no forces acting between the bodies and, at the boundary of the surface, tangential stresses are supposed to be zero. Also the traction bound constraint (see subchapter 2.3.2) is discretized in ( 3.11) The initial condition can be defined as follows: - Initialization of contact: All the tangential stresses from the previous state are zero. This is the case for pure rolling. - The calculation starts from a given tangential stress distribution . It is very useful when the input information is, e.g., a known distribution from a previous steady state. Analysis and modeling of unsteady effects within the wheel-rail contact Master’s thesis 48 The wheel and the rail are made of steel, whose material parameters are: Table 4.1. Default material parameters Parameter Symbol Value [units] Average Shear modulus G 82000 [GPa] Poisson’s coefficient 0.28 [-] Young’s modulus E 210000 [MPa] The contact patch is an ellipse generated by Hertz’s theory, using radii (see Subchapter 2.2), which are taken from the wheel profile S1002 and rail profile UIC60 (see 4.1.1): The resulting semi-axes of the contact ellipse (see chapter 2.2) are: a = 8.56 mm b = 4.08 mm Analysis and modeling of unsteady effects within the wheel-rail contact Master’s thesis 49 Figure 4.2. Description of the discreet ellipse, the semi-axes and the elements. 4.1.2 Cinematic and dynamic parameters Four other important parameters are the velocity of the wheelset, the coefficient of friction, the normal load and the creepages. - Velocity V: The default value used is V = 1 m/s - Coefficient of friction µ: The contact between both bodies in standard conditions (RH<10%, T = 20ºC), with a standard viscosity, generates a coefficient µ = 0.4, which can be used for investigations. - Normal load FN: A default value can be used for investigations is FN = 100 kN. - Creepages: In this thesis, spin creepage and lateral creepage are neglected (see Subchapter 2.3.1). So, only the longitudinal creepage is considered. If the Analysis and modeling of unsteady effects within the wheel-rail contact Master’s thesis 50 creepage is considered static in the simulation, a fixed value of cx = 0.005 is used. However, in most cases, the creepage has been varied over time. Therefore, two different input functions have been used: A step function that simulates a fast change of the contact conditions. At the beginning, the system was in steady rolling. Then, the creepage is changed instantly to a different value. The simulation lasts until steady state is reached again. A sine function that simulates changing conditions in the contact. Here, the system was initially in pure rolling, i.e., all tangential stresses were zero. Then the creepage was changed according to a sine function. 4.1.2.1 Creepage as a step function (step mode) In this case, the input function is a Heaviside step function, which is often used in the mathematics of control theory and signal processing to represent a signal that switches on at a certain instant of time and then it keeps constant. It is named after the mathematic Oliver Heaviside. The Heaviside function is the integral of the Dirac delta function: H′ = δ. This is written as: (4.1) Analysis and modeling of unsteady effects within the wheel-rail contact Master’s thesis 51 Another definition is: (4.2) Thus, in the case that the creepage is given as a step function, it is defined as (4.3) The behavior of this creepage over time is shown in Figure 4.3. Figure 4.3. Step function used as a default input function for the creepage. Analysis and modeling of unsteady effects within the wheel-rail contact Master’s thesis 52 4.1.2.2 Creepage as a sine function (sine mode) In this case, the creepage oscillates according to ( 4.4) Where C0 is the average value, Cp is the amplitude; F is the frequency. When the creepage is described by a sine function, its default values are C0 = 0.005; Cp = 0.005 and F = 10 Hz. Figure 4.4 shows the creepage over time for this case. Figure 4.4. Sine function used as a default input function for the creepage in sine mode. Analysis and modeling of unsteady effects within the wheel-rail contact Master’s thesis 53 To investigate the influence on the traction coefficient, three different frequencies were chosen: 10 Hz, 30 Hz, and 60 Hz. It is assumed that there is a deviation of the traction coefficient from the steady state. To investigate this, it is supposed that these frequencies are high enough to show these deviations. The time dependency of the resulting creepages for these three frequencies is shown in Figure 4.5. Figure 4.5. Sine functions with three different frequencies used as input functions for the creepage in sine mode. Analysis and modeling of unsteady effects within the wheel-rail contact Master’s thesis 54 4.2 Results of the unsteady effects from CONTACT This part contains the results of the initial investigations done with the program CONTACT (see Subchapter 2.3.5). In the following, the influence of parameters (coefficient of friction, normal loads, geometry, etc.) on the time dependency of the traction coefficient is presented. All the calculations were done using the two different input functions for the creepage described in Subchapter 4.1.2. 4.2.1 Influence of the frequency on the traction coefficient in sine mode First, a sine function of F = 10 Hz was used to describe the creepage (see Subchapter 4.1.2.2). In this case, the time dependency of the traction coefficient is shown in Figure 4.6. A friction coefficient µ = 1.0 was selected in order to avoid full sliding in case of high creepages. The time delay between the maximums of both signals is visible. Analysis and modeling of unsteady effects within the wheel-rail contact Master’s thesis 55 Figure 4.6. The creepage and the corresponding traction coefficient. There is a time delay between both curves ( This figure is for a µ = 1.0. The traction coefficient is delayed with respect to the creepage ( ). This time delay was also studied for higher frequencies (30 Hz and 60 Hz): Figure 4.7 shows the transient behavior of traction coefficient for a creepage that is changed with a frequency of 30 Hz. The amplitude of the traction coefficient (blue) is lower than for 10 Hz. The time delay is 4.29 ms. This delay was calculated on the second period to get a better accuracy because in the first period the traction coefficient starts from the steady state and needs time to stabilize. Analysis and modeling of unsteady effects within the wheel-rail contact Master’s thesis 56 Figure 4.7. The creepage and the corresponding traction coefficient for a frequency of 30 Hz. The results are similar to 10 Hz (see Figure 4.6), but the amplitude of the traction coefficient (blue) is lower than for 10 Hz. The time delay is round 4 ms. This delay was calculated on the second period to get a better accuracy because in the first period the traction coefficient starts from the steady state and needs time to stabilize. Analysis and modeling of unsteady effects within the wheel-rail contact Master’s thesis 57 Figure 4.8. The creepage and the corresponding traction coefficient for a frequency of 60 Hz The traction coefficient is attracted by the creepage, with a time delay of 3.68 ms. Three full periods are shown in order to investigate the stability. The traction coefficient takes 1.5 periods to be stable. The amplitude of the traction coefficient is lower than for 30 Hz (See Figure 4.7). Figure 4.8. shows the results for 60 Hz. The traction coefficient follows the creepage, with a time delay of 3.68 ms. Although the delay is similar to the two previous frequencies, due to the different frequencies it is necessary to divide the delay by the period. Table 4.2. Delay ratio for different frequencies. Frequency F [Hz] Period Ts [s] Time delay [ms] Delay ratio 10 0.1 4 0.040 30 0.0333 4.29 0.128 60 0.0167 3.68 0.220 Analysis and modeling of unsteady effects within the wheel-rail contact Master’s thesis 64 Figure 4.13. Tangential stresses using a step function for the creepage, for FN = 150 kN (a) and FN = 100 kN (b). In both plots, the colors correspond to the same time steps. The traction bound is higher and the contact ellipse is larger in (a) due to the normal load increases. This increase is not high enough to see differences on the geometry. Peaks in (a) are not caused by the numerical method (see Figure 4.15). In (b) the slip area is larger, there is more percent of the contact area that is in slip for t = 8Δt, which means the traction coefficient is higher. Figure 4.13 shows a comparison of the local tangential stress distributions using the same input function for the creepage but different normal loads FN = 150kN in (a) and FN = 100kN in (b). There are differences in: The traction bounds: A higher normal load leads to a higher value of the maximum of the normal stress distribution pn0 (Eq ( 2.33)) The area of contact: A higher normal load generates a larger ellipse. In the plots, the variation is too small to be noticed. But in Figure 4.14, the normal load is 10 times higher, so the contact area is much bigger. Analysis and modeling of unsteady effects within the wheel-rail contact Master’s thesis 65 Figure 4.14. Tangential stresses using a step function for the creepage, for FN = 1000 kN. This plot shows that the contact area increases with the normal load. The average value of stresses: For the second time step (red line in (b)), the curve reaches the traction bound (at x≈ 5.75mm) instead of the same curve in Figure 4.13(a) (x ≈ 6.5 mm). This is because of the Traction bound constraint. The result is a higher tangential force in (a). The shape of the stresses: In both cases, there are some peaks. To check if they are caused by the numerical iterations, a simulation using more grid points (a 100x100 grid) is done (Figure 4.15). There are some peaks too, which means they are not caused by the discretization. Analysis and modeling of unsteady effects within the wheel-rail contact Master’s thesis 66 Figure 4.15. Tangential stresses using a step function for the creepage, FN = 100 kN and a 100x100 grid. There are some peaks, which mean they are not caused by the discretization. Adhesion-slip areas: In Figure 4.13(a), the traction bound is smaller, so the slip area is bigger than in Figure 4.13(b). This directly affects the traction coefficient as is shown in Figure 4.16. Analysis and modeling of unsteady effects within the wheel-rail contact Master’s thesis 67 Figure 4.16. Traction coefficient for FN = 100 kN (blue) and FN = 150 kN (red), using a creepage step function. The steady states are shown in dotted lines. The initial gradients are different because of the normal loads. For higher normal loads, the traction coefficient takes a longer time to reach the steady state. In contrast to the case of the coefficients of friction (Figure 4.10), in this situation, the initial gradients are different. This effect is caused by the definition of the traction coefficient: for higher normal stresses, the traction coefficient decreases, and vice versa. (see Eq.( 2.31)) Despite having a lower steady state traction coefficient, the red curve takes longer to reach the steady state value, because of the smaller gradient caused by a lower normal load (Eq.( 2.31). Then, for higher normal loads, the traction coefficient takes a longer time to reach the steady state. When describing the time dependency of the creepage as a sine function, the effect observed is a hysteresis phenomenon. Again, the normal loads FN are varied. This time, six different normal loads were chosen. The results are shown in Figure 4.17: Analysis and modeling of unsteady effects within the wheel-rail contact Master’s thesis 68 Figure 4.17. Hysteresis effects using a sine function of 10 Hz for describing the creepage for six different normal loads. The initial condition chosen for all the curves is the steady state condition for a cx = 0.005. For higher normal loads, the hysteresis area increases and the hysteresis is also present for higher creepages: There is a longer time delay. There, the area of the hysteresis increases for higher normal loads. To investigate the influence of the frequency on the hysteresis behavior, two more simulations were done. One for F = 30 Hz and one for F = 60 Hz (Figure 4.18). There, the hysteresis shape changes according to the frequency, in a similar way to Figure 4.9. Analysis and modeling of unsteady effects within the wheel-rail contact Master’s thesis 69 Figure 4.18. Hysteresis effects using a sine function, for different frequencies and FN = 150 kN. The dotted line shows the steady state by FASTSIM. An increasing frequency generates a deviation from the steady state. For low frequencies the area becomes bigger but from F = 30 Hz to 60 Hz, the area decreases. 4.2.4 Contact patch’s geometry The geometry of the contact ellipse is also important, because the area of the contact is related to the stresses. As was shown in Subchapter 2.2, the normal load generates the final geometry of the ellipse. A circular contact patch was created that had the same area as the default one, which is elliptical: = Analysis and modeling of unsteady effects within the wheel-rail contact Master’s thesis 70 The resulting tangential stress distributions caused by the elliptical and the circular patch were compared and are shown in Figure 4.19. Figure 4.19. The stresses distributions of a circular patch (a) and the elliptical patch (b). The steps are the same in both situations. For the first two time steps (dark blue and light blue), x are higher in (a), but the slide area is bigger in (b), with respect to the semi-axe a. In steady state, the traction coefficient is higher for an elliptical patch. For initial time steps, the gradient of the traction coefficient is higher in the case of a circular patch. This is shown in Figure 4.20. Analysis and modeling of unsteady effects within the wheel-rail contact Master’s thesis 71 Figure 4.20. The traction coefficient evolution using circular and elliptical contact patches. The initial gradient is higher for the circular patch but its steady state value is lower. 4.2.5 Average value of the sine function This case of study is only made for creepages described by a sine function. The default input function presented in Subchapter 4.1.2.2 was defined in a way that the creepage was always positive. To investigate the effect of negative creepages, the parameter C0 is varied. The creepages become negative, if the average value C0 is smaller than 0.005 and the amplitude keeps constant. The values used are: And the resulting functions are shown in Figure 4.21. Analysis and modeling of unsteady effects within the wheel-rail contact Master’s thesis 72 Figure 4.21. Creepage input functions used. The red curve is the default sine function. Then, parameter C0 decreases until the sine function is symmetric with respect to the creepage (X-axis). Figure 4.22. The tangential stresses are influenced by the alternate and symmetric sine function (C0 = 0) and the distributions show different behavior from states 1 to 4. In 1, the x are at their maximum. In 2, the creepage decreases but the particles close to the trailing edge have memory of the previous creepage and are higher than in Analysis and modeling of unsteady effects within the wheel-rail contact Master’s thesis 73 case of steady state (This generates the time delay). In 3, the creepage (sine mode) is in the negative semiperiod because particles close to the leading edge presents negative stresses and particles far from the leading edge are still positive. In 4, all of the stresses are negative. These creepages generate positive and negative tangential stresses. Figure 4.22. shows this influence: There are 4 time steps presented which are described below. 1. are maximum. The tangential stresses distribution is the highest possible and there is not full sliding. 2. The creepage decreases but the particles close to the trailing edge have memory of the previous creepages and are higher than in case of steady state (This generates the time delay). 3. For negative creepage, the leading edge presents negative stresses. In steady state all the stresses should be negative, but on the right side they are still positive. 4. The stresses are all negative. The traction coefficient and the hysteresis are shown in Figure 4.23. Focusing on the curve corresponding to the sine function with C0 = 0 (blue), the hysteresis shows symmetry. And if Figure 4.22 and Figure 4.23 are compared, the conclusion drawn is that for different average values of the sine function, the hysteresis is different. Analysis and modeling of unsteady effects within the wheel-rail contact Master’s thesis 80 Figure 4.27 describes the evolution of the traction coefficient over time. The traction coefficient calculated with the UBA is higher than with CONTACT, but they both converge to the same steady state value. Next, UBA and CONTACT are compared using a sine function (see 4.1.2.2) with a frequency F = 10 Hz. The results are shown in Figure 4.28. The hysteresis generated by the UBA is smaller than the one generated by CONTACT. Figure 4.28. The hysteresis generated by UBA is inside the hysteresis generated by CONTACT. In conclusion, the results of the UBA agree qualitatively and quantitatively, but there are some differences taking into account all the uncertain wheel rail modeling. The differences are smaller than measurable values, and this approach shows effects that agree the results of CONTACT. Analysis and modeling of unsteady effects within the wheel-rail contact Master’s thesis 81 5 SUMMARY&CONCLUSIONS In the present work, unsteady effects in the wheel/rail contact have been investigated. The main effect is the hysteresis that appears when the longitudinal creepage changes fast and repeatedly. This effect is caused by a time delay between the creepage and the stresses, because the sum of the displacements is different for every time step and for every particle. These situations are easily feasible in reality, due to vehicle dynamics. These unsteady effects were first studied using CONTACT, which is based on Kalker’s Exact Theory. While its computational effort is very high, this algorithm is able to simulate unsteady effects for different conditions. Different parameters were varied and their influence on the resulting traction coefficient was investigated. This was done for two different modes: In the step mode, the time dependency of the creepage was given by a step function. And in sine mode, a sine function was used to describe the transient behavior of the creepage. In conclusion, time delays between the creepage and the traction coefficient was observed, which in turn meant that hysteresis in the traction characteristic were visible. The factors studied that lead to larger time delays are: - The increasing of the coefficient of friction µ. - The increasing of the normal loads FN. - The shape of the contact geometry. - The increasing of the frequency in case of a creepage in sine mode. Analysis and modeling of unsteady effects within the wheel-rail contact Master’s thesis 82 Next, a simplified model was created that is able to reproduce unsteady effects, called UBA. First, the UBA model was compared to an algorithm called FASTSIM, which is only able to calculate steady state solutions. The maximum relative error was below 1%. Then, the unsteady behavior of the UBA was compared to CONTACT. The UBA agrees qualitatively and quantitatively to CONTACT, but there are some differences taking into account all the uncertain wheel rail modeling. But its main advantage is that the computational effort of the UBA is much lower and may be applied in MultiBody System simulations in the future. Analysis and modeling of unsteady effects within the wheel-rail contact Master’s thesis 83 APPENDIX A TABLE A1. The Kalker’s coefficients of the linear theory of rolling contact for elliptic contact areas. ( = ln (16/g2)) C11 C22 C23= -C32 C33 g=a/b υ=0 0.25 0.5 υ=0 0.25 0.5 υ=0 0.25 0.5 υ=0 0.25 0.5 0.0 0.1 2.51 3.31 4.85 2.51 2.52 2.53 0.334 0.473 0.731 6.42 8.28 11.7 0.2 2.59 3.37 4.81 2.59 2.63 2.66 0.483 0.603 0.809 3.46 4.27 5.66 0.3 2.68 3.44 4.80 2.68 2.75 2.81 0.607 0.715 0.889 2.49 2.96 3.72 0.4 2.78 3.53 4.82 2.78 2.88 2.98 0.720 0.823 0.977 2.02 2.32 2.77 0.5 2.88 3.62 4.83 2.88 3.01 3.14 0.827 0.929 1.07 1.74 1.93 2.22 0.6 2.98 3.72 4.91 2.98 3.14 3.31 0.930 1.03 1.18 1.56 1.68 1.86 0.7 3.09 3.81 4.97 3.09 3.28 3.48 1.03 1.14 1.29 1.43 1.50 1.60 0.8 3.19 3.91 5.05 3.19 3.41 3.65 1.13 1.25 1.40 1.34 1.37 1.42 0.9 3.29 4.01 5.12 3.29 3.54 3.82 1.23 1.36 1.51 1.27 1.27 1.27 g=b/a 1.0 3.40 4.12 5.20 3.40 3.67 3.98 1.33 1.47 1.63 1.21 1.19 1.16 0.9 3.51 4.22 5.30 3.51 3.81 4.16 1.44 1.59 1.77 1.16 1.11 1.06 0.8 3.65 4.36 5.42 3.65 3.99 4.39 1.58 1.75 1.94 1.10 1.04 0.954 0.7 3.82 4.54 4.58 3.82 4.21 4.67 1.76 1.95 2.18 1.05 0.965 0.852 0.6 4.06 4.78 5.80 4.06 4.50 5.04 2.01 2.23 2.50 1.01 0.892 0.751 0.5 4.37 5.10 6.11 4.37 4.90 5.56 2.35 2.62 2.96 0.958 0.819 0.650 0.4 4.84 5.57 6.57 4.84 5.48 6.31 2.88 3.24 3.70 0.912 0.747 0.549 0.3 5.57 6.34 7.34 5.57 6.40 7.51 3.79 4.32 5.01 0.868 0.674 0.446 0.2 6.96 7.78 8.82 6.96 8.14 9.79 5.72 6.63 7.89 0.828 0.601 0.341 0.1 10.7 11.7 12.9 10.7 12.8 16.0 12.2 14.6 18.0 0.795 0.526 0.228 0.0 Analysis and modeling of unsteady effects within the wheel-rail contact Master’s thesis 84 APPENDIX B This appendix includes the auxiliary functions used by the UBM. Getcoord Function Getcoord (from “To get a coordinate”) is a function that has the following definition: [X,Y]= getcoord(mx,my,a,b) X,Y are arrays of mx×my coordinate values. The input parameters a, b are the semiaxes of the ellipse. Two variables x, y are defined: ( B.1) This function calculates the matrices element by element according to these equations: ( B.2) These equations describe the position of each element and where the calculations are referenced to. In this case, it is the center of the element (center of the rectangle). Then, the output matrices X, Y have this appearance: Analysis and modeling of unsteady effects within the wheel-rail contact Master’s thesis 85 ( B.3) And always the following equalities are satisfied: These equalities cause the origin is in the center of the ellipse. GetCXX Function GetCXX needs only three input parameters:  (Poisson’s ratio) and semiaxes a and b. It calculates the three Kalker’s coefficients explained in chapter 2.4.6. These parameters are very important to compare the steady-state results to FASTSIM and check the accuracy of the model. The output is matrix L which has got the following structure: In three-dimensional case, matrix L is always 3×3 and symmetric. Tractionbound is a subroutine that calculates the equation of the TB. Its definition is: TB= tractionbound(mx,my,a,b,X,Y,p0) The output TB is a my×mx array that contains all the values. p0 parameter is the maximum of the Hertz distribution (chapter 2.3) Analysis and modeling of unsteady effects within the wheel-rail contact Master’s thesis 86 LIST OF FIGURES Figure 2.1. Two elastic bodies rolling on each other................................................ 12 Figure 2.2. Representation of a moving coordinate system ..................................... 13 Figure 2.3. Representation of a fixed coordinate system. ........................................ 14 Figure 2.4. Hertzian contact. .................................................................................... 17 Figure 2.5. Graphic representation of the deformation of two bodies loaded.. ........ 19 Figure 2.6. Wheelset geometry and creep forces. ..................................................... 23 Figure 2.7. Description of the normal and tangential contact. ................................. 24 Figure 2.8. Geometry of the wheel in y-z plane. ...................................................... 26 Figure 2.9. A Hertzian distribution to the normal problem (a) and the. ................... 28 Figure 2.10. Application of the traction bound constraint. In (a) the module .......... 29 Figure 2.13. Hertz’s Theory (a) and the Simplified theory (b) ................................. 34 Figure 2.12. Different distributions of the normal stresses. (a) ................................ 34 Figure 2.14. Comparison of the tangential stress distribution of both theories. ...... 36 Figure 3.2. Unsteady Brush Algorithm (UBA) with the UBM program. ................. 42 Figure 4.1. a) Rail profile UIC60 b) Wheel profile S1002 ....................................... 47 Figure 4.2. Description of the discreet ellipse, the semi-axes and the elements. ..... 49 Figure 4.6. The creepage and the corresponding traction coefficient.. .................... 55 Figure 4.7. The creepage and the corresponding traction coefficient for. ................ 56 Figure 4.8. The creepage and the corresponding traction coefficient for. ................ 57 Figure 4.9. Frequency dependency of the hysteresis for a coefficient ..................... 59 Figure 4.10. Time dependency of the traction coefficient for different. .................. 60 Figure 4.11. Traction coefficient for different coefficient of friction,. .................... 61 Figure 4.12. Displacement distributions in X-direction for different time. .............. 62 Figure 4.13. Tangential stresses using a step function for the creepage,. ................ 64 Analysis and modeling of unsteady effects within the wheel-rail contact Master’s thesis 87 Figure 4.14. Tangential stresses using a step function for the creepage,. ................ 65 Figure 4.15. Tangential stresses using a step function for the creepage,. ................ 66 Figure 4.16. Traction coefficient for FN=100 kN (blue) and FN=150 kN (red),. ..... 67 Figure 4.17. Hysteresis effects using a sine function of 10 Hz for........................... 68 Figure 4.18. Hysteresis effects using a sine function, for different frequencies. ..... 69 Figure 4.19. The stresses distributions of a circular patch (a) and the elliptical. ..... 70 Figure 4.20. The traction coefficient evolution using circular and elliptical. .......... 71 Figure 4.21. Creepage input functions used. The red curve is the. .......................... 72 Figure 4.22. The tangential stresses are influenced by the alternate. ....................... 72 Figure 4.23. The hysteresis effect of the traction coefficient for. ............................ 74 Figure 4.24. The traction coefficient obtained by UBA and compared to ............... 76 Figure 4.25. Absolute error of the traction coefficient using FASTSIM. ................ 77 Figure 4.26. The stresses for a creepage of cx=0.001 calculated for. ....................... 79 Figure 4.28. The hysteresis generated by UBA is. ................................................... 80 Analysis and modeling of unsteady effects within the wheel-rail contact Master’s thesis 88 REFERENCES Abramowitz, M., & Stegun, I. (1964). Handbook of Mathematical Functions. New York: Dover Publications. Alonso, A., & Giménez, J. (2008). Non-steady state modelling of wheel-rail contact problem for dynamic simulation of railway vehicles. En Vehicle System Dynamics (págs. 179-196). Taylor & Francis. Ayasse, J.-B., & Hugues, C. (2006). Wheel-rail contact. En S. Iwnicki, Handbook of Railway Vehicle Dynamics. Taylor & Francis Group. Boussinesq, J. (1885). Application des Potentiels à l'Etude de l'Equilibe et du Mouvement des Solides Elastiques. Paris: Gauthier-Villars. Cerruti, V. (1882). Mem. fis. mat. Roma: Accademia dei Lincei. Duvaut, G., & Lions, J.-L. (1972). Les inéquations en Mécanique et en Physique. Paris: Dunod. (2006). EN 13715:2006. Brussels: European Commitee for Standardization. Fichera. (1964). Problemi elastostatici con vincoli unilaterale: il probleme di Signorini con ambigue condizioni al contorno. Memorie Accademie Nazionale dei Lincei. Gladwell, G. (1980). Contact problems in the classical theory of elasticity. Alphen a/d Rijn: Sijthoff and Noordhoff. Goldstein, H. (2001). Classical Mechanics (3rd ed.). Addison-Wesley. Gurtin, M. E. (1981). An Introduction to Continuum Mechanics. New York: Academic Press. Analysis and modeling of unsteady effects within the wheel-rail contact Master’s thesis 89 Jahnke, E., & Emde, F. (1943). Tables of Functions. New York: Dover Publications. Johnson, K. (1987). Contact mechanics. Cambridge: Cambridge University Press. Kalker, J. (1967). On the rolling contact of two elastic bodies in the presence of dry friction. Thesis Delft. Kalker, J. (1971). A minimum principle of the law of dry friction with application to elastic cylinders in rolling contact. Journal of Applied Mechanics. Kalker, J. (1979). The computation of three-dimensional rolling contact with dry friction. International Journal of Numerical Methods in Engineering. Kalker, J. (1983). Two algorithms for the contact problem in elastostatics. En E. J. Kalousek, R. V. Dukkipati, & G. Glaswell, Contact Mechanics and the Wear of Wheel-Rail Systems. University of Waterloo Press. Kalker, J. (1990). Three-Dimensional Elastic Bodies in Rolling Contact. In J. Kalker. Kalker, J., & Goedings, H. (1972). A programme for three-dimensional steady state rolling. Faculty TWI, Section TA: Report TU Delft. Logston, C. F., & Itami, G. S. (1980). Locomotive friction creep studies. En Journal of Engineering for Industry 102 (págs. 275-281). Love, A. (1952). A Treatise on the Mathematical Theory of Elasticity. Cambridge. Meierhofer, A., Hardwick, C., Lewis, R., Six, K., & Dietmaier, P. (2012). Third Body Layer - Experimental results and a model describing its influence on the traction coefficient. En Inproceedings: 9thn International Conference on Contact Mechanics and Wear of Rail/Wheel Systems.