# Development of a Numerical Chain to Optimize Railway Axles with Respect to Fatigue Damage

Sofiane Saad, Vincent Magnier, Philippe Dufrenoy, Eric Charkaluk, François Demilly

► To cite this version:

Sofiane Saad, Vincent Magnier, Philippe Dufrenoy, Eric Charkaluk, François Demilly. Development of a Numerical Chain to Optimize Railway Axles with Respect to Fatigue Damage. Key Engineering Materials, 2014, 611-612, pp.1683-1693. 10.4028/www.scientific.net/KEM.611-612.1683. hal-00998495

HAL Id: hal-00998495

https://hal.science/hal-00998495

Submitted on 15 Jul 2021

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# Development of a numerical chain to optimize railway axles with respect to fatigue damage

Sofiane SAAD<sup>1,2,a</sup>, Vincent MAGNIER<sup>1,2,b</sup>, Philippe DUFRENOY<sup>1,2,c</sup>,
Eric CHARKALUK<sup>2,3,d</sup>, François DEMILLY<sup>4,c</sup>

<sup>1</sup>Univ. Lille Nord de France, F-59000 Lille, France

<sup>2</sup>Univ. Lille1, Laboratoire de Mecanique de Lille, F-59650 Villeneuve d'Ascq, France,

<sup>3</sup>CNRS, UMR 8107, F-59650 Villeneuve d'Ascq, France,

<sup>4</sup>Valdunes SAS, F-59125 Trith Saint Leger, France

<sup>a</sup>sofiane.saad@polytech-lille.fr <sup>b</sup>vincent.magnier@polytech-lille.fr

<sup>c</sup>philippe.dufrenoy@polytech-lille.fr <sup>d</sup>eric.charkaluk@ec-lille.fr <sup>e</sup>Francois.DEMILLY@valdunes.com

**Keywords:** Dynamic FE simulation, Hot upsetting, Assembly operation, Fatigue load, Residual stress.

**Abstract.** In today's competitive business environment, it has become increasingly important to reduce manufacturing and raw materials cost. For this purpose, an innovative process of design and manufacturing railway axles is developed. It is based on forging hollow axles which allows a significant reduction in steel consumption. In this work, we tried to analyze how these modifications induced by this new process and design impact the service behavior and particularly the durability face to cyclic loadings that can lead to fatigue failure. In the present study, a numerical chain has been developed going from the simulation of the manufacturing process up to the analysis in fatigue. In the first step, the forging process is modeled in order to predict the residual stress field and the initial plastic strain. From this initial condition, the assembly operation of the wheel on the axle is simulated before the redistribution of stresses and strains under cyclic load. The final objective is to obtain the cyclic loading paths, in order to provide the data needed for the analysis of fatigue.

## Introduction

Railway axles were investigated at the very beginning of fatigue research and design. Regarding to safety, this part is among the most important components in railway vehicles (Fig. 1).

Fortunately accident are rather seldom. According to [1] and [2] the service failures of these components is very low today. The Railway Safety Performance Report of the European Railway Agency (ERA) of 2011 [3] shows that the number of broken axles in the European Union between 2006 and 2009 was in total 329, corresponding to approximately one fracture per two billions of axle kilometers for Europe.

In order to maintain this safety record and reduce the manufacturing process costs, a new upset axles forging process is studied which allows a significant reduction in steel consumption.

The present paper deals with the development of a numerical chain which firstly allows to model this new manufacturing process, then the press fitted process of the wheel on the axles and finally the application of a rotary bending load on the structure for fatigue analysis.

## Numerical chain background

The advantage of such a numerical chain (Fig. 2), is to keep throughout the three finite elements simulations the residual stress and plastic strain resulting from the different stages of manufacturing and assembly operation. The aim of this numerical chain is to integrate process stress and plastic strain into the fatigue design of the axle.
<figure><img src="image_2.png" /><figcaption>Fig. 1: Schematic drawing of the wheel/axle assembly</figcaption></figure>

<caption>Fig. 2: Different stages of the numerical chain</caption>

Bibliographic study shows that a particular focus is placed upon damage tolerance analysis, which uses the typical tools of fracture mechanism, fatigue growth and statistics [4, 5, 6, 7]. The purpose of this work is to control the railway axles lifetime as well as the required inspection interval.

So far, there is no agreed method for collecting and estimating the fatigue phenomenon on the axle. Many works like [8, 9], have studied the impact of residual stress resulting from the assembly process of the wheel/axle on the fatigue behavior. Among the weakness of this finite element simulation they does not take into account residual manufacturing stress. Other study [10] has developed a numerical chain to integrate only the hardening process effect on the wheel, in order to study the rail/wheel fatigue contact.

In our case, ABAQUS commercial finite element modeling software has been used to develop the numerical chain to 1 : 3 scale. In the first step of the numerical chain we simulate the forging process to obtain the final axisymmetric geometry of the axle. Then we proceed by transferring results to simulate the assembly process of the wheel on the axle. The final simulation is related to the application of fatigue load on the structure.

Hot forging process

The objective here is to provide a model which generates from a pipe the final geometry of the axle in one step. This process of hot upset forging is illustrated in Fig. 3, is developed by Barriol & Dallière Industries in order to reduce the steel consumption and manufacturing cost.

Several studies have focused on simulating the forging process of railway wheel [11] and axles[12]. In this paper we present the axisymmetric finite elements model to upsetting the axle and the different employed approach to reduce the computation time.

<figure><img src="image_3.png" /><figcaption>Fig. 3: Schematic drawing of the upset axle forging</figcaption></figure>
**Material law.** In this section we present the EA1N (XC38) material flow which is widely adopted for manufacturing railway axles [13]. A high strain, strain rate and temperature in upsetting process have strong influence on material's flow. To get better analysis results a material model that combine all these aspects is necessary. So in this simulation, the Thermo-elastic/viscoplastic Johnson and Cook [14] law (Eq. 1) is opted. The parameters used in the FE simulation are detailed in Tab. 1.

$$ \sigma = [A + B\bar{\epsilon}^n] \left[ 1 + C \ln \left( \frac{\bar{\epsilon}}{\dot{\epsilon}_0} \right) \right] \left[ 1 - \left( \frac{T - T_r}{T_{fus} - T_r} \right)^m \right] \quad (1) $$

Where $\sigma$ is the equivalent stress, $A$ is the initial yield strength (MPa) of the material at room temperature, $B$ is the hardening modulus, $n$ is the work-hardening exponent, $C$ is the coefficient dependent on the strain rate, $m$ is the thermal softening coefficient and $\bar{\epsilon}$ represents the equivalent plastic strain. The strain rate $\dot{\epsilon}$ is normalized with a reference strain rate $\dot{\epsilon}_0$. $T_r$ and $T_{fus}$ represent the room temperature and melting temperature, respectively.

<table><caption>Table 1: Workpiece Johnson & Cook parameters [15]</caption><thead><tr><th>A(MPa)</th><th>B(MPa)</th><th>n</th><th>c</th><th>m</th></tr></thead><tbody><tr><td>490</td><td>600</td><td>0.21</td><td>0.015</td><td>0.6</td></tr></tbody></table>

**Numerical methods.** An explicit scheme [16] was adopted to solve the coupled mechanical and thermal problem (ABAQUS/Explicit). The advantages of these method is to allow the simulation of the dynamic phenomena present in the forging operation and it presents a better contact manager with simple algorithms. In addition to that, this method allows the use of mass-scaling approach to reduce the simulation computing time. It is based on the artificial increasing of the material density $\rho$. As a result an increase of the time increment $\Delta t = \frac{L_e}{C_d}$ and thus reduces the computation time. $C_d$ is the speed of wave propagation $C_d = \sqrt{\frac{E}{\rho}}$ and $L_e$ is the length of the smallest element edge.

The excessively use of such approach leads to wrong results known as Hourglass problem. Many works [17, 18, 19] prove that to have a good calculation accuracy we must have the kinematic energy less than 5% to the internal energy of the FE simulation. In addition to that we should maintain an artificial energy less than 10% of the internal energy.

In order to avoid the mesh distortion, the Arbitrary Lagrangien-Eulerian (ALE) [20] approach is used. This approach combines the advantages of the Lagrangian formulation which is summarized in a good definition of boundaries conditions. As the advantages of the Eulerian formulation which permits to introduce a large distortions. So this method is presented as the most appropriate formulation for modeling large deformation problems and the study of coupled problems [21]. Concerning the modeling of friction, a Coulomb law was implemented with a constant coefficient.

At the end of forging stage, a cooling step is applied with air convection.

**Hot upset forge results.** In this section, we present numerical results of the upsetting operation. The first results from the forging simulation gives a good correlation in terms of the outer and internal shape of the axle manufactured by BDI shown in Fig. 4. This comparison allows to show geometrical precision to predict the external and internal numerical shape of the axle. Fig. 5 represents the numerical result of the upsetting process in the end of forging step. In accordance with the forging tests, we notice that the simulation predicts an increased thickness of the axle journal and keeps the same thickness for the axle body.

After the cooling process, the residual stress in the wheel seat is represented by the Fig. 6. The stresses in the radial direction and the axial directions respectively S11 and S22 globally has a compressive fields. We notice that the upsetting process generates a compressive peaks in the two axles radius respectively -200MPa and -400MPa in the radial direction. Unlike in the axle direction there is a decrease in compressive stress at this specific location.
<figure><img src="image_4.png" /><figcaption>Fig. 4: Finale axle shape manufacturing by BDI</figcaption></figure>

<figure><img src="image_5.png" /><figcaption>Fig. 5: Axisymmetric numerical axle (PE11: Plastic stain in the radial direction)</figcaption></figure>

<figure><img src="image_6.png" /><figcaption>Fig. 6: Stress distribution after upsetting process</figcaption></figure>

## Press fitting operation

The assembly operation of the wheel on the axle is made by *Valdunes* company. This operation allows tightening of the wheel on the axle to ensure tightness between these two bodies. The principle is to provide a diameter of a portion of the axle, called wheel seat, greater than the diameter of the wheel hub (Fig. 7).

This operation is performed with a feed rate of $1mm.s^{-1}$. To ensure coaxiality between the axle and wheel a chamfer of $\alpha = 2$ degrees from the horizontal is included on the beginning of the axle (Fig. 7). This operation is validated if the force issue by the press included between the two master curves shown in (Fig. 8). This industrial validation of the assembly operation is very interesting method to confirm the numerical results.

The axisymmetric FE simulation of this operation is made with an Implicit scheme. The advantage of changing from the explicit scheme, used in the previous forging simulation, to the implicit is based on the fact that we do not model dynamic and thermal phenomena. Furthermore during the assembly operation the implicit method (ABAQUS/Standard) [22] presents more accurate results with a lower computational time. The results transfers are performed from the forging to the assembly simulation with Abaqus. The aim of such transfers is to keep the residual stress and plastic strain as an initial state of the assembly model.
<figure><img src="image_7.png" /><figcaption>Fig. 7: Assembly operation</figcaption></figure>

<figure><img src="image_8.png" /><figcaption>Fig. 8: Experimental force of the assembly</figcaption></figure>

**Material law.** Based on the fact that dynamic and thermal phenomena are neglected in this operation, an elastic-plastic material law is sufficient here. Two constitutive laws are used for this model, for the wheel and for the axle. A non-linear kinematic hardening plastic model [23] using Armstrong-Frederick's law [24] (Eq. 2) is affected to the wheel. Identifying parameters was done in this work [25] based on cyclic tests. Tests were performed for the ER7(XC48) material under reversed tensile and reversed torsion loadings. Tab. 2 summarizes the parameters used to modeling the wheel behavior with ER7 grade.

$$d\bar{X} = C d\bar{\epsilon}^{pl} - \gamma \bar{X} |d\bar{\epsilon}^{pl}| \quad (2)$$

$\bar{\epsilon}^{pl}$ is the plastic strain tensor, $\bar{X}$ is the backstress tensor, $k$: the yield stress and $C$, $\gamma$ are the hardening parameters of the Armstrong-Frederick constitutive law.

<table><caption>Table 2: Armstrong-Frederick's ER7 wheel parameters[25]</caption><thead><tr><th>E(MPa)</th><th>k(MPa)</th><th>C</th><th>γ</th></tr></thead><tbody><tr><td>208550</td><td>296</td><td>128900</td><td>276</td></tr></tbody></table>

In order to model the hardening of the axle in this simulation and the asymptotic response of the structure in the next fatigue simulation, we have chosen Lemaître and Chaboche law [26]. This behavior law combine isotropic (Eq. 3) and kinematic hardening model. Five parameters used to model the axle behavior.

$$R = R_0 + Q_{\infty}(1 - \exp(-b\epsilon_{eq}^{pl})) \quad (3)$$

Where $R_0$ is the yield stress at zero plastic strain and $Q_{\infty}$ and $b$ are material parameters. $Q_{\infty}$ is the maximum change in the size of the yield surface, and $b$ defines the rate at which the size of the yield surface changes as plastic straining develops.

**Assembly simulation results.** From the upsetting axle mechanical state, we proceed by the press fitting of wheel. The Fig. 9 shows the evolution of the wheel advance and the redistribution of radial stress.
<figure><img src="image_9.png" /><figcaption>Fig. 9: Radial stress redistribution in the wheel seat before and after assembly operation</figcaption></figure>

These radials stresses allows a tightening between the two pieces. A pressure of about 350 MPa is obtained in the wheel seat (Fig. 10). A presence of two peaks of compressions at the beginning of the range and at the end of a value of 680 MPa and 500 MPa respectively. These two peaks explained by formation of a bead at the beginning of the assembly operation which crashing in this position. The second peak is due to a geometrical singularity. This simulation is validated by using experimental forces data (Fig. 8). We proceed by comparing the force required for the numerical assembly (Fig. 11) with experimental data.

<figure><img src="image_10.png" /><figcaption>Fig. 10: Stress redistribution after assembly</figcaption></figure>

<figure><img src="image_11.png" /><figcaption>Fig. 11: Force of the assembly operation</figcaption></figure>

This comparison shows that the numerical effort of the assembly operation is between the two master curves, so it has the same magnitude order.

Fatigue simulation

The aim of this step simulation is to get the asymptotic response of the axle. Based on this stabilized
response a fatigue criterion will be used into the axle design. For the non-axisymmetric character of the
fatigue phenomenon, an implicit quasi-three-dimensional model based on the precedent axisymmet-
ric assembly simulation is developed (Fig. 12). A Python program is developed to ensure the results
transfers in terms of stress and plastic strain tensor between the axisymmetric model to the quasi-3D
*CAXARN* elements model.
<figure><img src="image_12.png" /><figcaption>Fig. 12: Development of the fatigue simulation</figcaption></figure>

After, we apply a service loading exerted on the journal by the subroutine *DLOAD*. The same behavior law as the assembly simulation are used. The combined isotropic and nonlinear kinematic hardening plastic model presents the advantage of modeling the Elastic shakedown, Plastic shakedown and Ratcheting that can occur in the axle.

**CAXA Elements type.** The choice of such elements is based on their ability to accurately predict the response bending load. They are used to give a simulated three-dimensional response by revolving a two-dimensional surface around the centerline of symmetry. The use of *CAXA* elements increases the efficiency of the model when compared to a true three-dimensional model while still maintaining accurate results. [27, 28] showed that a *CAXA* element is capable of simulating complicated model conditions normally done utilizing three-dimensional elements. By adopting such elements in the fatigue numerical simulation we can achieve a gain of computing time by about 20% compared to the three-dimensional models.

A three dimensional response is simulated using quasi three-dimensional Fourier analysis elements (*CAXA*) available within ABAQUS. These elements use standard isoparametric interpolation in the radial - symmetry axis plane, combined with the Fourier interpolation with respect to the angle of revolution. The asymmetric deformation is assumed to be symmetric with respect to the radial - symmetry axis plane at an angle equal to 0 or π [29].

**Rotary bending modeling.** In order to properly model the service load (Fig. 13) generated by the train weight coupled with the axle rotation a *DLOAD* subroutine is developed. The *DLOAD* is typically used when a load is a complex function of time and/or position.

<figure><img src="image_13.png" /><figcaption>Fig. 13: Problem contact simplification</figcaption></figure>

<figure><img src="image_14.png" /><figcaption>Fig. 14: Hertz Cylinder-Cylinder pressure</figcaption></figure>

In our case we model a Hertz pressure distribution cylinder-cylinder [30] shown by Fig. 14. The particularity of this load that is simulated the rotation of axles hence the importance of the change in pressure position respect to the angle θ with a certain frequency.
In the subroutine *DLOAD* we implement Gaussian function corresponding to the hertz contact distribution between rolling and axle (Eq. 4). This function *P* represent a great interest to model the contact pressure distribution and the maximum pressure by manipulating the *σ* parameter of the function. On the other hand, we can change the center position of the load by varying the *m* = *w*.*t parameter (*w* is the angular speed of the axle and *t* is the time). The Fig. 15 shows the shape of the pressure applied by this subroutine which presents the extreme load case.

$$P(\theta) = \frac{1}{\sigma\sqrt{2\pi}} \exp\left(-\frac{(\theta - m)^2}{2\sigma^2}\right) \quad (4)$$

<figure><img src="image_15.png" /><figcaption>Fig. 15: *DLOAD* Pressure distribution respectively for *m* = π/2 and *m* = π</figcaption></figure>

**Fatigue load results.** In this part we interest in the distribution of plastic strain in the axle and his evolution in each cycle under rotary load applied by the *DLOAD* subroutine. Based on the shape of the stress-strain loops, outcome of the element which undergoes the greatest plastic flow, we determine the structure response under a rotating bending load. In this example, we define a pressure in the subroutine *DLOAD* reflecting a 15KN load applied to the axle journal with a frequency of 20Hz, which corresponding to the axle rotation frequency. The Fig. 16 shows the evolution of the stress-strain loops for 500 cycles simulated. We deduce that the asymptotic behavior of the structure is the plastic shakedown.

<figure><img src="image_16.png" /><figcaption>Fig. 16: Stress-Strain loop</figcaption></figure>

## Summary

A numerical chain is developed to predict the residual stress and plastic strain throughout three simulations.
**Step 1:** In the first step, the upsetting operation is simulated with explicit scheme to take into account the hot forging dynamic phenomenon. A Johnson and Cook law was adopted to describe the material flow and an Arbitrary Lagrangian Eulerian formulation to avoid excessive mesh distortion. This simulation predict the final axisymmetric geometry of the axle and particularly the stress and strain distribution.

**Step 2:** In the second step, we proceed by transferring results in terms of the forging residual and plastic strain into the new simulation based on the implicit scheme. The implicit resolution is used to simplify the assembly model and reduce the computation time. Considering this mechanical state, we proceed by simulating the axisymmetric cold assembly operation to predict the redistribution of stress and plastic strain under the wheel seat. A combined isotropic and kinematic hardening behavior is used.

**Step 3:** The final simulation is related to the fatigue load applicable on the axle journal. The final assembly mechanical state is transferred using a specific Phyton algorithm. This implicit quasi-three dimensional simulation is developed using the CAXA elements based on the Fourier interpolation. The complexity of pressure load which combine the axle rotation and the wagon weigh, lead to the use of the DLOAD subroutine, which integrate a Gaussian pressure distribution. According to the structure response: Elastic shakedown, Plastic shakedown and Ratcheting a fatigue criterion will be used.

Finally, the paper improves the dimensional axle approach [31] by taking into account the manufacturing and assembly residual stress and strain. This method is optimized using different approach to reduce the calculation time. In terms of results, the different simulations present a good correlation with the experimental data.

## References

[1] R.A. Smith and S. Hillmansen. A brief historical overview of the fatigue of railway axles. Proceedings of the Institution of Mechanical Engineers, Part F: Journal of Rail and Rapid Transit, 218(4):267-277, 2004

[2] K. Hirakawa, K. Toyama, and M. Kubota. The analysis and prevention of failure in railway axles. International Journal of Fatigue, 20(2):135-144, 1998

[3] Railway Safety Performance in the European Union. European railway agency (ERA). Technical report, 2011

[4] M. Ognjanovic, A. Simonovic, M. Ristivojevic, and T. Lazovic. Research of rail traction shafts and axles fractures towards impact of service conditions and fatigue damage accumulation. Engineering Failure Analysis, 17(7-8):1560-1571, October 2010.

[5] U. Zerbst, S. Beretta, G. Köhler, A. Lawton, M. Vormwald, H.Th. Beier, C. Klinger, I. Cerny, J. Rudlin, T. Heckel, and D. Klingbeil. Safe life and damage tolerance aspects of railway axles a review. Engineering Fracture Mechanics, 2012

[6] U. Zerbst, M. Schödel, and H.Th. Beier. Parameters affecting the damage tolerance behaviour of railway axles. Engineering Fracture Mechanics, 78(5):793-809, March 2011

[7] T. Makino, T. Kato, and K. Hirakawa. Review of the fatigue damage tolerance of high-speed railway axles in Japan. Engineering Fracture Mechanics, 78(5):810-825, March 2011

[8] A. Yameogo, A. Carrado, A. M. Marechal, S. Pommier, C. Prioul, and A. Lodini. Residual stress redistribution due to cyclic loading in a railway Wheel/Axle assembly. Journal of Neutron Research, 12:63-68, September 2004
[9] A. Yameogo. Etude expérimentale et numérique de l'amorçage et de la propagation de fissures de fretting dans un assemblage roue/essieu ferroviaire. PhD thesis, Ecole Centrale de Paris, 2004

[10] A.M.G. Langueh. Une méthode de dimensionnement des roues ferroviaires sous sollicitations de contact entraînant un phénomène de fatigue et d'écaillage. PhD thesis, université Lille 1, 2013

[11] K. Davey, B.C. Miller, and M.J. Ward. Efficient strategies for the simulation of railway wheel forming. Journal of Materials Processing Technology, 118(1-3):389-396, December 2001

[12] L.-P. Lei, J. Kim, and B.-S. Kang. Analysis and design of hydroforming process for automobile rear axle housing by FEM. International Journal of Machine Tools and Manufacture, 40(12):1691-1708, September 2000

[13] Railway applications - wheelsets and bodies - axles - production requirements, 2003

[14] G.R. Johnson and W.H. Cook. A constitutive model and data for metals subjected to large strains, high strain rates and high temperatures, volume 21. 1983

[15] T.T. Öpöz and X. Chen. Finite element simulation of chip formation. School of Computing and Engineering Researchers' Conference, pages 166-171, December 2010

[16] A. Lucia Garcia. Etude thermo-mécanique et modélisation numérique de l'emboutissage à chaud de l'Usibor 1500. PhD thesis, Ecole des mines de Paris, January 2004

[17] H.H. Choi, S.M. Hwang, YH. Kang, J. Kim, and BS. Kang. Comparison of implicit and explicit finite-element methods for the hydroforming process of an automobile lower arm. The International Journal of Advanced Manufacturing Technology, 20(6):407-413, 2002

[18] Z. Wang and T. Nakamura. Simulations of crack propagation in elastic-plastic graded materials. Mechanics of Materials, 36(7):601-622, July 2004

[19] H.H. Wisselink and J. Huétink. 3D FEM simulation of stationary metal forming processes with applications to slitting and rolling. Journal of Materials Processing Technology, 148(3):328-341, May 2004

[20] D.J. Benson. An efficient, accurate, simple ale method for nonlinear finite element programs. Computer Methods in Applied Mechanics and Engineering, 72(3):305-350, 1989

[21] R. Boman. Développement d'un formalisme Arbitraire Lagrangien Euléerien tridimensionnel en dynamique implicite. Application aux opérations de mise en forme. PhD thesis, Université de Liège, March 2010

[22] P.G. Bergan, G. Horrigmoe, B. Brakeland, and T.H. Soreide. Solution techniques for nonlinear finite element problems. International Journal for Numerical Methods in Engineering, 12(11):1677-1696, 1978

[23] W. Prager. Problèmes de plasticité théorique. Dunod edition, 1975

[24] P.J. Armstrong and C.O. Frederick. Technical Report RD/B/N731, Berkeley Nuclear Laboratories, 1966

[25] A.M.G. Langueh, J-F. Brunel, E. Charkaluk, P. Dufrenoy, and J.-B. Tritsch. Effects of sliding on rolling contact fatigue of railway wheels. Fatigue Fract Engng Mater Struct, 36:515-525, October 2012

[26] J. Lemaître and J-L. Chaboche. Mécanique des matériaux solides. Dunod 2 edition, 2001
[27] M. Wierszycki, W.K. akol, and T. Lodygowski. ABAQUS users conference, 2006

[28] B. Sukumaran, W.O. McCarron, P. Jeanjean, and H. Abouseeda. Efficient finite element techniques for limit analysis of suctioncaissons under lateral loads. Computers and Geotechnics, 24(2):89-107, March 1999

[29] ABAQUS. Abaqus 6.11 Theory Manual. 2011

[30] H. Hertz. Über die berührung fester elastischer körper. Journal für reine und angewandte Mathematik, 92:156-171, 1882

[31] Norme européenne Applications ferroviaires essieux montés et bogies-essieux-axes porteurs
méthode de conception, NF EN 13103+A1, January 2011