Application of semi-analytical finite element method to analyze the bearing capacity of asphalt pavements under moving loads

Abstract

To facilitate long term infrastructure asset management systems, it is necessary to determine the bearing capacity of pavements. Currently it is common to conduct such measurements in a stationary manner, however the evaluation with stationary loading does not correspond to reality a tendency towards continuous and high speed measurements in recent years can be observed. The computational program SAFEM was developed with the objective of evaluating the dynamic response of asphalt under moving loads and is based on a semi-analytic element method. In this research project SAFEM is compared to commercial finite element software ABAQUS and field measurements to verify the computational accuracy. The computational accuracy of SAFEM was found to be high enough to be viable whilst boasting a computational time far shorter than ABAQUS. Thus, SAFEM appears to be a feasible approach to determine the dynamic response of pavements under dynamic loads and is a useful tool for infrastructure administrations to analyze the pavement bearing capacity.

Keywords

semi-analytical finite element method / bearing capacity / asphalt pavements / moving loads / dynamic response

Introduction

Infrastructure management and maintenance concepts are of upmost importance to ensure a long life whilst maintaining sufficient quality and aid the decision processes linked to rehabilitation strategies. The material performance as well as the remaining service life can be determined by analyzing the bearing capacity of asphalt pavements [¹,²]. Currently deflection measurements are the only dependable means of determining the bearing capacity of asphalt pavements [³,⁴]. Obtaining such measurements with conventional devices such as the falling weight deflectometer (FWD) necessitates the interruption of traffic. To make the process more economically feasible and to conduct the measurement along the full length high speed deflection measurement devices have been developed in recent years. These include the quest/dynatest rolling weight deflectometer (RWD), Swedish road deflection tester, Texas rolling dynamic deflectometer and traffic speed deflectometer (TSD) [⁵], which all operate at velocities between 40 km/h and 80 km/h [⁶].

An inverse mapping of the theoretical response model must be determined to characterize the bearing capacity of asphalt pavements which can be obtained by measuring deflections; this is referred to as the pavement backcalculation technique [⁷]. The processes in the backcalculation can generally be subdivided into forward calculations and backward calculations. Effectively the forward direction computes deflections in the pavement for given loading scenarios whilst the backward direction is used to adjust the calculation by comparing the computed values with measured deflections; parameter identification is used to determine the actual material parameters [⁷]. At present the forward calculation process is limited to approaches such as finite element (FE) or elastic layered theory and stationary deflection measurement methods. The more realistic approach of using dynamic moving loads requires significantly more computational resources resulting in low time efficiency; thus, further development is still needed.

In this research project the forward calculation, i.e., dynamic displacements of pavement structures are calculated under moving loads by applying a semi-analytical finite element method (SAFEM) [⁸–¹³]. The pavement structure problem was reduced to a 2D plane strain problem by exploiting spatial conditions of a typical pavement structure (Fig. 1) such as the fact that the geometric and material properties do not vary in one coordinate direction (z-direction in this case); however, loading conditions vary strongly in said direction. The computational time can be reduced significantly under the assumption that a Fourier series can be used to describe the displacement in the out-of-plane direction (y-direction in this case). By making use of the orthogonal properties such a class of problems can be reduced to a series of 2D FE-meshes [⁸–¹³].

The mathematical basis of SAFEM will be outlined, followed by the result verification by means of commercial FE-software ABAQUS and results obtained from field measurements. Finally, the findings are summarized and conclusions are drawn.

Description of semi-analytical finite element method

A Fourier series with z ranging from zero to a is used to represent the general shape functions defining the displacement variations in SAFEM (Fig. 2) [⁸–¹²]:

(1)u={\displaystyle \sum _{k=\mathrm{1}}^{\mathrm{6}}{N}_k(x,y,z)}{u}_k={\displaystyle \sum _{l=\mathrm{1}}^L{\displaystyle \sum _{k=\mathrm{1}}^{\mathrm{6}}\left\{\right[\overline{N_k}(x,y)\cos {\displaystyle \frac{l\pi z}{a}}]+[\overline{\stackrel{‾}{N_k}}(x,y)\sin {\displaystyle \frac{l\pi z}{a}}\left]\right\}{u}_k^l}}

With l is the term of the Fourier series; L is the number of terms considered; \overline{N_k} & \overline{\overline{N_k}} is the shape functions of a node in the XY plane(equal to displacement approximations in 2D problems).

Variation of loads along the z-axis is given by a loading function correspondingly [⁹]:

(2)f={\displaystyle \sum _{l=\mathrm{1}}^L\left\{\right[\bar{p}(x,y)\cos {\displaystyle \frac{l\pi z}{a}}]+[\overline{\stackrel{‾}{p}}(x,y)\sin {\displaystyle \frac{l\pi z}{a}}\left]\right\}{u}_k^l}

With \bar{p}(x,y) and \overline{\stackrel{‾}{p}}(x,y) is the external load function.

The pavement degrees of freedom are fixed at the boundaries z= 0 and z= a whilst unrestricted motion in the z-direction.

The displacement function in u, v and w can be written as such [¹⁰]:

(3)U=\left\{\begin{array}{l}u\\ v\\ w\end{array}\right\}={\displaystyle \sum _{l=\mathrm{1}}^L{\displaystyle \sum _{k=\mathrm{1}}^{\mathrm{6}}{N}_k\left[\begin{array}{lll}\sin {\displaystyle \frac{l\pi z}{a}}& \mathrm{0}& \mathrm{0}\\ \mathrm{0}& \sin {\displaystyle \frac{l\pi z}{a}}& \mathrm{0}\\ \mathrm{0}& \mathrm{0}& \cos {\displaystyle \frac{l\pi z}{a}}\end{array}\right]}}\left\{\begin{array}{l}{u}_k^l\\ v_k^l\\ w_k^l\end{array}\right\}={\displaystyle \sum _{l=\mathrm{1}}^L{N}^l\cdot U^l}

With u_k^l, v_k^l and w_k^l is the displacements of the node at the term of the Fourier series along x-, y- and z-directions, respectively.

The loading function for the pavement analysis can be simplified as [⁹]:

(4)f={\displaystyle \sum _{l=\mathrm{1}}^{L}p(x,y)\sin {\displaystyle \frac{l\pi z}{a}}}={\displaystyle \sum _{l=\mathrm{1}}^L{\{p\}}^l}

(5)p(x,y)={\displaystyle \sum _{t=\mathrm{1}}^n\left({\displaystyle \frac{\mathrm{2}{P}_t}{l\pi }}\right)}\left[\cos {\displaystyle \frac{l\pi }{a}}{Z}_{t\mathrm{1}}-\cos {\displaystyle \frac{l\pi }{a}}{Z}_{t\mathrm{2}}\right]

With P_t is tire loading pressure; Z_t1 is z coordinate of tire load beginning; Z_t€2 is z coordinate of tire load ending.

Beginning with the displacement of a single element the resulting strain and stress can be determined with geometric and physical equations; resulting in the strain-displacement matrix B_k^l as follows:

(6)B_k^l=\left[\begin{array}{l}\begin{array}{lll}{\displaystyle \frac{\partial N_k}{\partial x}}\sin {\displaystyle \frac{l\pi z}{a}}& \mathrm{0}& \mathrm{0}\\ \mathrm{0}& {\displaystyle \frac{\partial N_k}{\partial y}}\sin {\displaystyle \frac{l\pi z}{a}}& \mathrm{0}\\ \mathrm{0}& \mathrm{0}& -{\displaystyle \frac{l\pi }{a}}{N}_k\sin {\displaystyle \frac{l\pi z}{a}}\end{array}\\ \begin{array}{lll}{\displaystyle \frac{\partial N_k}{\partial y}}\sin {\displaystyle \frac{l\pi z}{a}}& {\displaystyle \frac{\partial N_k}{\partial x}}\sin {\displaystyle \frac{l\pi z}{a}}& \mathrm{0}\\ \mathrm{0}& {\displaystyle \frac{l\pi }{a}}{N}_k\cos {\displaystyle \frac{l\pi z}{a}}& {\displaystyle \frac{\partial N_k}{\partial y}}\cos {\displaystyle \frac{l\pi z}{a}}\\ {\displaystyle \frac{l\pi }{a}}{N}_k\cos {\displaystyle \frac{l\pi z}{a}}& \mathrm{0}& {\displaystyle \frac{\partial N_k}{\partial x}}\cos {\displaystyle \frac{l\pi z}{a}}\end{array}\end{array}\right]\mathrm{.}

The element stiffness matrix can be represented through a sub-matrix {(K^{lm})}^e by means of the principle of minimum potential energy [¹⁰]:

(7){(K^{lm})}^e={\displaystyle {\iiint }_{vol}{(B^l)}^{T}D{B}^m\mathrm{d}x\mathrm{d}y\mathrm{d}z}

The force vector is typically:

(8){(F^l)}^e={\displaystyle {\iiint }_{vol}{(N^l)}^T{\{p\}}^l\mathrm{d}x\mathrm{d}y\mathrm{d}z}

With Eqs. (6), (7) the stiffness matrix of one element includes [¹⁰]:

(9)\begin{array}{l}{I}_{\mathrm{1}}={\displaystyle {\int }_{\mathrm{0}}^a\sin {\displaystyle \frac{l\pi z}{a}}\cos {\displaystyle \frac{m\pi z}{a}}}\mathrm{d}z\\ I_{\mathrm{2}}={\displaystyle {\int }_{\mathrm{0}}^a\sin {\displaystyle \frac{l\pi z}{a}}\sin {\displaystyle \frac{m\pi z}{a}}}\mathrm{d}z\\ I_{\mathrm{3}}={\displaystyle {\int }_{\mathrm{0}}^a\cos {\displaystyle \frac{l\pi z}{a}}\cos {\displaystyle \frac{m\pi z}{a}}}\mathrm{d}z\end{array}

The orthogonal properties of the integrals ensure that the following relations are true:

(10)I_{\mathrm{2}}=I_{\mathrm{3}}=\{\begin{array}{ll}{\displaystyle \frac{\mathrm{1}}{\mathrm{2}}}a,& \mathrm{f}\mathrm{o}\mathrm{r}l=m\\ \mathrm{0},& \mathrm{f}\mathrm{o}rl\ne m\end{array}

In the case that both l and m are odd or even numbers the integral I₁ becomes zero. When I₁ is zero all terms in the B_l matrix containing I₁ become zero as well, resulting in the diagonalization of the {\left(K^{lm}\right)}^e matrix. Thus, the stiffness matrix can be reduced to resemble the following:

(11)\left[\begin{array}{l}\begin{array}{llll}{K}^{\mathrm{11}}& & & \end{array}\\ \begin{array}{llll}& K^{\mathrm{22}}& & \end{array}\\ \begin{array}{llll}& & \ddots & \end{array}\\ \begin{array}{llll}& & & K^{LL}\end{array}\end{array}\right]\left\{\begin{array}{l}{U}^{\mathrm{1}}\\ U^{\mathrm{2}}\\ ⋮\\ U^L\end{array}\right\}+\left\{\begin{array}{l}{F}^{\mathrm{1}}\\ F^{\mathrm{2}}\\ ⋮\\ F^L\end{array}\right\}=\mathrm{0}

The final equation shows that system of equations is reduced to contain exactly L individual problems.

(12)K^{ll}{U}^l+F^l=\mathrm{0}

To increase the computational efficiency, it is important to decouple problems and allow for parallel computing as opposed to sequential schemes. In this case the equilibrium equations are decoupled for the harmonics of the Fourier series (Eq. (12)) [¹⁴]. The time independent FE analysis is conducted with Eq. (12). As the dynamic loading and the respective response include a temporal component all objects are four dimensional (x, y, z, t); therefore Eq. (12) must be rewritten. The time coordinates must be included in the FE algorithm. The derivation of this is found in numerous literature sources [¹⁵] and will not be described further.

Verification of the safem

Results from the commercial FE-software ABAQUS [¹⁶] and in field measurements on a test track of German Federal Highway Research Institute (BASt) were used to analyze the accuracy of SAFEM. The experiments and calculations were conducted based on a truck passing at 30 km/h due to spatial restrictions of the test track (Fig. 3).

The pavement length was set to be 4800 mm, which is 20-fold the minimum required loading length for a computational calculation. The width was set to 3750 mm corresponding to the width of the test track [¹⁸–²⁰]. The pavement layer materials and the respective material properties are given in Table 1. The first three layers are totally bound. Each contact layer among the asphalt base course, gravel base layer, frost protection layer and sub-grade was defined as being partially bound.

The geometric information of the truck from the experiment is given in Fig. 4. To simplify the model only the first wheel of the left axle was regarded. The distances of the wheels are estimated to be far enough away from the left wheel on the first axle that their influence is negligible. The loading area is assumed to be rectangular with dimensions 240 mm × 290 mm, and exhibit a loading pressure of 0.515 MPa. The center of the loading path is 1100 mm from the left boundary of the test track.

The respective meshes from SAFEM and ABAQUS are depicted in Fig. 5. The SAFEM mesh consisted of 2D 6-node triangular elements in the XY plane (Fig. 5a). The ABAQUS mesh consisted of 3D 10-node quadratic tetrahedron elements (Fig. 5b).

In both meshes the element size increases from the loading path towards the periphery. In ABAQUS the finest elements within the loading area are 24 mm × 29 mm (two right-angle sides); in SAFEM the finest elements are 24 mm (the surface side). The total number of increments is 210; the increment time is 0.00288 s as a result of the loading speed of 30 km/h. Thus, the total step time for the load to move through the pavement section is 0.6048 s. Some stages of the moving load on the pavement surface are given in Fig. 6 (ABAQUS).

The surface deflection from ABAQUS and SAFEM are compared in Fig. 7 when the loading arrived at the center (X= 2400 mm). The dynamic response of the pavement results in the deflection curve not being axis symmetrical. The maximum deflection results from ABAQUS and SAFEM are 0.221 mm and 0.212 mm. Given that ABAQUS exhibits a coefficient of determination R² of 0.999 the relative error of SAFEM of 4.21% is regarded as acceptable, taking the different meshes into account.

The results from SAFEM and ABAQUS were checked against obtained measurement data. The comparison was based on measurements by sensors embedded in the test track which measured the strain on the bottom of the base course and the vertical tensile stress on the top of the gravel base course in the moment when the load was exactly above.

The results are given in Table 2. The computational strains and stresses exceed the field measurements on the underside of the asphalt base course and are smaller on the top of the gravel base course.

Discrepancies of 20% can be tolerated as a result of the uncertainties and fluctuations of in field measurements [¹⁷–¹⁹] resulting in well coinciding strain values obtained from the computation. However, the stress values of both computational models are far beyond the 20% threshold. This may be attributed to non-linear material properties of the gravel base layer and will be the focus of future research. The two computational models exhibit very good correspondences indicating the reliability of SAFEM.

The verification results prove the reliability and accuracy of SAFEM. To assess the factor of computational efficiency the computational time of SAFEM as well as ABAQUS were run on an Intel Core Duo 3.4 GHz, 32 GB RAM; the results are given in Table 3. The different mesh algorithms and resulting element and node numbers only allow for a qualitative comparison. SAFEM appears to be 21 times more efficient than ABAQUS.

Summary and conclusion

This research project analyzes the suitability of SAFEM for predicting asphalt responses under moving loads with regard to strain and stress values as well as computational efficiency. The accuracy is assessed based on comparisons with the commercial FE-software ABAQUS as well as field measurements. The computed strain at the critical point in SAFEM corresponds well with field measurements. The computed stresses of both SAFEM as well as ABAQUS on the other hand, are far beyond a tolerable margin; this will be the focus of future research. Furthermore the analysis of computational time showed that SAFEM is far more efficient than ABAQUS.

The conclusion to be drawn is that the SAFEM exhibits the potential to reliably predict dynamic responses of asphalt pavement under moving loads and represents a theoretical tool for administrations in the scope of analyzing pavement bearing capacity under high speed measurement regimes and can be included in strategic processes. SAFEM also supports the implementation of material properties such as viscoelasticity for asphalt and nonlinear elasticity for the sub-base of the pavement. The inclusion of these enhancements will make SAFEM far more accurate at predicting the response of asphalt pavements under moving loads.

References

[1] Hou Y, Wang L, Yue P, Pauli T, Sun W. Modeling Mode I Cracking Failure in Asphalt Binder by Using Nonconserved Phase-Field Model. Journal of Materials in Civil Engineering, 2014, 26(4): 684–691

[2] Hou Y, Wang L, Pauli T, Sun W. Investigation of the Asphalt Self-healing Mechanism Using a Phase-Field Model. Journal of Materials in Civil Engineering, 2015, 27(3): 04014118

[3] Ferne B W, Langdale P, Roundm N, Faircloughm R. Development of a Calibration Procedure for the U.K. Highways Agency Traffic-Speed Deflectometer. Transportation Research Record, 2009a, 2093(1): 111–117

[4] Ferne B W, Langdale P, Roundm N, Faircloughm R. Development of the UK Highways Agency Traffic Speed Deflectometer. Bearing Capacity of Roads, Railways and Airfields, London, Champaign, IL: Taylor & Francis Group, 2009b, 409–418.

[5] Muller W B, Roberts J. Revised Approach to Assessing Traffic Speed Deflectometer Data and Field Validation of Deflection Bowl Predictions. International Journal of Pavement Engineering, 2013, 14(4): 388–402

[6] Baltzer S, Pratt D, Weligamage J, Adamsen J, Hildebrand G. Continuous Bearing Capacity Profile of 18,000 km Australian Road Network in 5 Months. In: the 24th ARRB Conference Proceedings, Melbourne, Australia, 2010

[7] Goktepe A B, Agar E, Lav A H. Advances in Backcalculating the Mechanical Properties of Flexible Pavements. Advances in Engineering Software, 2006, 37(7): 421–431

[8] Fritz J J. Flexible Pavement Response Evaluation Using the Semi-Analytical Finite Element Method. International Journal of Materials and Pavement Design, 2002, 3(2): 211–225

[9] Hu S, Hu X, Zhou F. Using Semi-Analytical Finite Element Method to Evaluate Stress Intensity Factors in Pavement Structure. Rilem International Conference on Cracking in Pavements, 2008

[10] Zienkiewicz O C, Taylor R L. The finite element method for solid and structural mechanics. 6th edition, Elsevier Butterworth-Heinemann, Oxford, 2005

[11] Liu P, Wang D, Oeser M, Chen X. Einsatz der Semi-Analytischen Finite-Elemente-Methode zur Beanspruchungszustände von Asphaltbefestigungen. Bauingenieur, 2014, 89(7/8): 333–339 (in German)

[12] Liu P, Wang D, Oeser M. The Application of Semi-Analytical Finite Element Method Coupled with Infinite Element for Analysis of Asphalt Pavement Structural Response. Journal of Traffic and Transportation Engineering, 2015, 2(1): 48–58

[13] Liu P, Wang D, Oeser M. Leistungsfähige semi-analytische Methoden zur Berechnung von Asphaltbefestigungen. Tangungsband, 3.Dresdner Asphalttage, 2013 (in German)

[14] Kim J R, Kim W D, Kim S J. Parallel Computing Using Semianalytical Finite Element Method. AIAA Journal, 1994, 32(5): 1066–1071

[15] Zienkiewicz O C, Taylor R L. The Finite Element Method, Volume 1 The Basis. 5th edition, Elsevier Butterworth-Heinemann, Oxford, 2000

[16] ABAQUS. ABAQUS analysis user’s manual. Dassault Systemes Simulia Corp., 2011

[17] Gohl S. Vergleich der gemessenen mechanischen Beanspruchungen der Modellstraße der BASt mit den Berechnungsergebnissen ausgewählter Programme. TU Dresden, Professur für Straßenbau, Diplomarbeit, 2006 (In German).

[18] Rabe R. Bau einer instrumentierten Modellstraße in Asphaltbauweise zur messtechnischen Erfassung der Beanspruchungssituation im Straßenaufbau. AP 03 342, interner Bericht, Bundesanstalt für Straßenwesen, Bergisch Gladbach, 2004 (In German)

[19] Rabe R. Messtechnische Erfassung der Beanspruchungen im Straßenaufbau infolge LKW-Überfahrten über eine Modellstraße in Asphaltbauweise‘ AP 04 342, interner Bericht, Bundesanstalt für Straßenwesen, Bergisch Gladbach, 2007 (In German)

[20] Rabe R. Angaben zum Aufbau der Modellstraße und Angabe von ausgewählten Ergebnissen und Materialkennwerten. Bergisch Gladbach, 2014 (In German)

RIGHTS & PERMISSIONS

Higher Education Press and Springer-Verlag Berlin Heidelberg