Abstract: Backcalculating layer moduli, using the deflections measured by the FWD (Falling Weight Deflectometer) device, is the basis for structurally evaluating pavements when considering maintenance and rehabilitation options. Many algorithms have been developed to perform the backcalculation. Most of them are based on simplified assumptions such as the elastic behavior of the pavement and static load. These algorithms use the peak values of loads and deflections recorded by the FWD sensors. The dynamic analysis, to the contrary, accounts for the dynamic nature of the load and factors such as material inertia and damping. These algorithms use as inputs both load and deflection time histories measured by the FWD. In this research, three software programs, MODULUS 6.0, ELMOD 5.0 and EVERCALC 5.0, were used to do the elasto-static back-analysis and the DBSID program was employed for the dynamic backcalculation. The FWD data gathered from different test sites, including the Zanjan-Tabriz, Eivanekey-Garmsar and Garmsar-Semnan freeways and also Rafsanjan airport have been used in this evaluation. The results have been compared, the performance of each program has been evaluated and the best software for the sites under study is suggested.
INTRODUCTION
The modulus of elasticity is the major parameter in the remaining life analysis and in the design of new pavements and overlays. This critical value is computed through the backcalculation method, which uses the FWD raw data (Tawfiq, 2003). Most of the common backcalculation algorithms, which use the peak values of load and deflection at each sensor during every drop, assume that the FWD load is applied statically. The pavement system is modeled as a layered elastic system with linear or nonlinear (stress-dependent) materials. A forward analysis subroutine is used to calculate the theoretical deflections under the known load. Attempts are made to converge these calculated deflections to the measured ones. In other words, the moduli are predicted by minimizing the error between the measured and calculated deflections in every step of iteration. But in reality, the FWD applies an impulse load to the pavement that generates body waves and surface waves, which travel at finite velocities and are recorded at different times by the geophones (Lytton, 1989). This full time history data of both load and deflections can also be measured during the FWD test. This additional data can be used to backcalculating other parameters of pavement materials.
BACKCALCULATION OF LAYER MODULI, STATIC AND DYNAMIC FEATURES
Classical pavement-subgrade system analysis is based on multi-layer linear elastic analysis, infinite pavement layer dimensions and a semi-infinite subgrade. Linear elastic response assumptions are valid when there is no discontinuity, for example at joints in rigid pavements, or cracks. The assumed static loading conditions are also a simplification of reality (William, 1999).
The duration of loading runs from 20-65 msec depending on the type of FWD device. The FWD provides excellent structural information if it is fully utilized. From a review of the sensor peak values for each measurement location, it becomes clear that points farther from the center of load plate attain their peak values later than locations closer to the plate. This time differences known as phase difference of deflections and is a measure of the velocity of propagation of the shockwaves (Fig. 1). These types of time history data contain useful information that if well utilized would contribute to improved accuracy of the structural pavement evaluation. The use of only peak values discards this potentially very useful information (Matsui and Hachia, 2006).
Static Feature
Objective Function
The objective of most backcalculation programs is to determine a set of
moduli that will minimize an error term between the computed deflection and
the measure deflection (Chua, 1989). Because random errors
are unavoidable, the number of measuring deflection sensors must be greater
than the number of unknowns (for example, moduli in the linear elastic case
or material parameters for the nonlinear case) to be backcalculated (Uzan
et al., 1989).
Since, the accuracy of the sensors, which is a major source of error, is expressed in relative terms, the objective function to minimize should therefore be expressed in relative terms. The objective function is (Uzan et al., 1989):
| (1) |
| ε2 | = | Squared error |
| wim | = | Measured deflection at sensor i |
| wic | = | Computed deflection at sensor i |
| S | = | The number of sensors |
| Fig. 1: | Measured FWD data (Matsui and Hachia, 2006) |
Backcalculation Methods
The computed deflection wic can be expressed as follows
(Uzan et al., 1989):
wic = fi(Xj) |
(2) |
Xj is the unknown variables, including modulus of elasticity, poisson’s ratio, thickness of layer, pressure, contact area,… (j varies from 1 to n).
Generally, there are two groups of backcalculation programs. The first is the iterative approach in which a forward calculation scheme (typically a linear elastic program) is used within the iterative process (Chua, 1989). In other words, any solution to Eq. 1 calls for a solution to Eq. 2. The number of calls depends on the minimization algorithm used. So, the pattern search technique requires a hundreds of calls of the deflection computation program. A simplified description of the iterative process used for adjusting the modulus values is shown in Fig. 2. This illustration is for one deflection and one layer. For multiple layers and deflections, the solution is obtained by developing a set of equations that define the slope and intercept for each deflection and each unknown layer modulus as follows (Van Cauwelaert et al., 1989; Tawfiq, 2003):
Log (deflectionj) = Aij+
Sji (log Ei) |
(3) |
Where: |
| A | = | Intercept |
| S | = | Slope |
| j | = | 1, 2,…, ND (ND = No. of deflections) |
| i | = | 1, 2,…, NL (NL = No. of layers with unknown moduli) |
The data base method uses forward calculation scheme to build a data base from which regression equations are either formulated to determine the layer moduli or are used within interpolation techniques to compute the deflections, thus avoiding the use of a forward calculation scheme in the iterative process (Uzan et al., 1989).
In the case of linear elasticity and a circular contact area, Eq. 2 can be written as (Uzan et al., 1989):
| (4) |
| Fig. 2: | Relationship between deflection and modulus (Van Cauwelaert et al., 1989) |
| Fig. 3: | Deflection of elastic layers (Tawfiq, 2003) |
Where: |
| p | = | Pressure |
| a | = | Radius of contact area |
| Ei | = | Modulus of elasticity of the ith layer on the subgrade |
| Esg | = | Subgrade modulus of elasticity |
This equation is similar to the well known equation for the deflection of two layer systems, as shown in Fig. 3 (Tawfiq, 2003):
| • |
| • |
F2 is a function of E1/E2 and h1/a . The value of fi, for sensor i, which is a function of modular ratios, can be obtained from the layered system program and used as a data base (Tawfiq, 2003).
Backcalculation Programs
Most of the common backcalculation programs with elasto-static basis
use an elastic layer program for forward analysis and an error optimization
algorithm. In this study, the programs MODULUS 6.0, ELMOD 5.0 and EVERCALC
5.0 have been used to perform traditional back analysis.
Description of MODULUS6.0
This program was developed at Texas Transportation Institute (TTI) for the
Texas DOT and uses the data-base method for back-analysis. It uses WESLEA program
as a forward calculation subroutine. WESLEA is based on the multilayer linear
elasto-static theory that is traditionally used for the purposes of flexible
pavement analysis (William, 1999). The program uses WESLEA
to generate a data base of deflection bowls by assuming different modular ratios.
A pattern search technique is then used to determine the set of layers moduli
that produce a deflection basin that fits the measured one (William,
1999).
The problem with the program is that the units are in in-lb units only and that it only handles 7 sensors.
Up to four unknown layers, with or without bedrock are allowed. It is enough for the user to define the material type and the thickness for every layer. The program suggests the moduli ranges and poisson’s ratios and also the depth to bedrock and includes them in the analysis. The subgrade depth can be altered by the user if site specific information is available (William, 1999).
Description of ELMOD5.0
ELMOD is an acronym for Evaluation of Layer Moduli and Overlay Design. This
program was developed by Dynatest Consulting Inc. and accepts DYNATEST-FWD format.
The theory of Odemark-Boussinesqe transformed section is the basis of forward
analysis. Besides, it uses an iterative method for backcalculation (William,
1999). There are 3 modes of backcalculations available (Elmod
5 Quick Start Manual, 2006):
| • | Radius of curvature: The outer geophone readings are used to determine non-linear characteristics of subgrade and the inner geophones to determine the upper pavement layer moduli. The stiffness of remaining layers is then calculated based on the overall pavement response to the applied load. This approach analysis up to 4 layers |
| • | Deflection basin fit: The methodology starts with a set of estimated moduli for the pavement structure. The theoretical deflection bowl for this pavement structure is calculated. The error between the measured and calculated deflections is then assessed. The moduli in the structure are increased/decreased by a small amount (typical 10%) and if the error in either of these bowls is less than the original bowl, this is taken to be a better solution. Up to 5 layers can be estimated by this mode |
| • | FEM/LET/MET: With this option backcalculation may be carried out either with Finite Element Method (FEM), Linear Elastic Theory (LET) or the Method of Equivalent thickness. FEM makes use of a modified version of an axial symmetric finite element program, originally developed by Wilson at University of California by Duncan et al. (1968). LET makes use of WESLEA and MET is similar to the method used in the option Basin Fit, but with a simpler use of adjustment factors |
Temperature correction is also engaged. Moduli range and poisson’s ratio (equal to 0.35 for all layers in every case) is determined by the software internally (Elmod 5 Quick Start Manual, 2006).
The major assumption is that the layers are homogeneous, isotropic and linear elastic, except the subgrade which can have a nonlinear definition (Elmod 5 Quick Start Manual, 2006):
| (5) |
| σ1 | = | Major principal stress from external loading |
| pa | = | Reference stress, often taken equal to atmospheric pressure (0.1 MPa) |
C and n are constants. n varies from 0 (linear elastic material) to -0.5.
Description of EVERCALC5.0
EVERCALC5.0 was developed by Dr. Joe Mahoney at the University of Washington for Washington DOT. This program uses the WESLEA computer code for the forward analysis and a modified Gauss-Newton algorithm for solution optimization. It can handle up to 5 layers, 10 sensors and 12 drops per station (William, 1999).
The program terminates when one or more of the following conditions are satisfied (Everseries User’s Guide, 2005):
| • | Deflection tolerance: |
| (6) |
wim is the measured deflection and wic is the calculated deflection at sensor i. S is the number of sensors and n is the number of layers
| Moduli tolerance: |
| (7) |
Eki and E(k+1)i are the i-th layer moduli at the kth and (k+1)th iteration and m is the number of layers with unknown moduli
| Number of iterations has reached the maximum number of iterations. At every iteration a maximum of (m+1) calls to WESLEA is made, where m is the number of layers with unknown moduli |
Dynamic Feature
System Identification Methods
System identification methods were developed by electrical engineers who
were interested in determining the characteristics of a filter by using the
input and output signals of the filter and an assumed model of the filter. The
characteristics of the filter are changed systematically using a search technique
until the model produces an output that is acceptably close to that of the filter.
The procedure is in fact exactly analogous to what is being done in backcalculating
the moduli of pavements (Lytton, 1989).
An impulse load is applied to the pavement at a remote distance and the surface motion is sensed as signals received by two sensors, which are separated by a distance, x. The pavement between the two sensors is regarded as the filter, the characteristics of which are to be determined. The unknowns in each layer are its modulus, percent damping, thickness and unit weight (Lytton, 1989).
Description of DBSID (Dynamic Backcalculation with System Identification Method): The Fourier transform is commonly used and widely implemented as an integral transform in signal processing which decomposed e.g., a function of time into a series of harmonics. There are several efficient algorithms for implementing the discrete Fourier transform, commonly referred to as Fast Fourier Transform (FFT) (Westover and Guzina, 2006). The FWD load is first decomposed into its frequency components using the FFT. Then the transfer function, which defines the response of the pavement system to a steady state unit load, is evaluated. Finally, the transfer function is multiplied by the FFT of the load to determine the Fourier transform of the displacements (Fernando and Liu, 2002).
The dynamic analysis presented by DBSID is based on the forward model implemented in the FWD-DYN program. Based on research by Uzan et al. (1989), the developers of DBSID followed a time-domain fitting approach in the dynamic backcalculation. This approach was readily implemented using the FWD-DYN program, which predicts the displacement history for each FWD sensor, given the material properties and thickness of each pavement layer. FWD-DYN uses Fourier superposition to predict pavement response due to the impulse load by the FWD. In this method, FWD-DYN performs an inverse FFT on this Fourier transform to determine the time history of the displacements for each FWD sensor. The developers of DBSID, coupled a system identification routine to the FWD-DYN program. The Asphalt Concrete (AC) layer is modeled as a visco-elastic (or damped elastic) material and the base and subgrade layers are modeled as damped elastic materials (Fernando and Liu, 2002).
Here, point minimizing of the errors between discrete values of deflections is not the concern, but the best fit of the deflection histories.
Pavement Dynamic Material Model
Viscoelastic Material Model (User’s Manual for DBSID, 1993)
Asphalt Concrete (AC) material is a matrix of solid aggregate particles
with asphalt binder. The bonding between the asphalt binder and the aggregate
particles results in viscoelastic behavior of AC mixture. In the DBSID
procedure, the viscoelastic behavior is modeled using the following three-parameter
model that is considered realistic for loading times representative of
highway traffic speed:
D(t) = D0+ D1 tm |
(8) |
| D(t) | = | The creep compliance function, defined as D(t) = e (t)/s0 |
| e(t) | = | The total axial strain measured at time t |
| s0 | = | The applied constant stress in a creep compliance test |
| D0, D1 and m | = | Model parameters |
When t = 0, D(0) = D0, from the three-parameter viscoelastic model. This means that D0 (D0 = 1/E0) stands for the elastic response at a very small loading time. For the asphalt concrete layer, the following parameters are used in backcalculation:
| H | = | Thickness of asphalt concrete layer |
| g | = | Unit weight of asphalt concrete |
| v | = | Poisson’s ratio |
| D0 | = | 1/E0, representing the elastic response due to the solid matrix |
| D1 | = | Creep compliance constant representing the nonlinear viscous response |
| m | = | Exponent for nonlinear time dependence |
Damped-Elastic Material (User’s Manual for DBSID, 1993)
For pavement layers that are modeled as damp-elastic, only two parameters
are needed to define the response of the material: the elastic modulus
(E) and the damping ratio (ξ). The damping ratio is input to the
program and is not backcalculated. Only the elastic modulus is backcalculated.
If the surface layer is modeled as viscoelastic, the values D0, D1 and m are backcalculated. In the other damped elastic layers the moduli are backcalculated. So, together with the depth of subgrade, 6 unknowns will be defined through backcalculation.
FIELD DATA COLLECTION
The FWD data used in the analysis was collected on the freeways and airport site listed below:
| • | Eivanekey- Garmsar, km 00.00 to 26.55 |
| • | Garmsar- Semnan, km 00.00 to 104.13 |
| • | Zanjan-Tabriz, km 205.00 to 236.00 |
| • | Rafsanjan Airport, in three different parts: Taxiway, Runway and Apron |
It is noteworthy that in the last item, the HWD test has been performed.
COMPARISON BASIS
In order to evaluate the results, a comparison basis is essential. The subgrade layers are the most complex layers in the pavement structure due to their physical nature and construction practices. The backcalculation process is based on the assumption that surface deflection at a certain offset is characteristic of the elastic modulus at a certain depth (William, 1999). Several of the backcalculation programs start their analysis from the outer geophone which determines the resilient modulus of subgrade. So, accurate assumption of subgrade is vital to the rest of the analysis.
In this study, the results of soil mechanic tests for subgrade layers have been used as a ground-truth to compare with the moduli values derived from the FWD data and the ones suggested by the laboratory tests.
INPUTS TO THE PROGRAMS
General Inputs for All 4 Programs
Thickness of Layers
All of the pavement structures under study have been made up of an
asphalt surface, on a crushed stone base. These two layers rest on a subgrade
whose type is defined by mechanical tests. The thickness of each layer
at each station has also been measured.
Moduli Range and Poisson’s Ratio
In the backcalculation process (Table 1) from SHRP backcalculation
report was used to define the surface and base layers moduli ranges. For the
subgrade, the measured CBR values were used to estimate the resilient modulus
using the relationships given in Table 2 and Fig.
4 from Iran Highway Asphalt Paving Code (2000).
Based on the assembled data the following ranges were thought reasonable for the modulus of subgrade for each site:
| • | Zanjan- Tabriz: 70 to 105 MPa |
| • | Eivanekey- Garmsar and Garmsra- Semnan: 70 to 140 MPa |
| • | Rafsanjan Airport: 170 to 210 MPa |
Additional Inputs to DBSID
Unit Weight (γ)
The value of 2250 kg m-3 was selected for all asphalt concrete
layers. For base and subgrade Table 3 was used for determining
γ, if it was not clearly mentioned in the soil mechanic test results. Lastly,
the following values for the unit weight were selected:
| • | Zanjan-Tabriz: 1600 kg m-3 for the subgrade, 1902 kg m-3 for the base layer |
| • | Eivanekey-Garmsar and Garmsra-Semnan: 1920 kg m-3 for both base and subgrade layers |
| • | Rafsanjan Airport: 1920 kg m-3 for both base and subgrade layers |
| Table 1: | Default moduli ranges and Poisson’s Ratio values (William, 1999) |
| Table 2: | Relationship between CBR and resilient modulus for the CBR values less than 25 |
| Fig. 4: | Relationship between CBR and resilient modulus for the CBR values greater than 25 |
| Table 3: | Table of typical m values and unit weight for base course and subgrade materials (Lytton, 1989) |
Damping Ratio
From the sensitivity analysis performed it was found that the damping ratio
does not significantly effect the deflection time history a lot within the range
of this parameter (Fernando and Liu, 1993). The suggested
default value of 4% was fixed in the DBSID program.
After the completion of the input data, the programs were run and the outputs were reviewed for reasonableness.
DISCUSSION
The elasto-static programs terminate when the convergence condition between deflections is satisfied. The Root Mean Squared Errors (RMSE’s) higher than 5% were not accepted. Reviewing the results is always necessary because every predicted modulus may not be acceptable. For example, the moduli which are equal to the lower or upper limit of the acceptable range are difficult to justify given the known quality of the materials.
In EVERCALC program, it was found that moduli convergence usually happened earlier than deflection convergence. So, the program was rerun by changing the moduli range to achieve the best deflection fit.
The DBSID program tries to create the best fit between measured and computed deflection time histories. Figure 5a and b present the results from this study. The output files contain the RMSE between all deflection values (during the 60 msec load duration) and also the errors between the peak values of calculated and measured deflections at each sensor. The surface moduli are most dependent on the deflection measured at the sensor in the middle of the load plate (W1 at 0 inches), so in the final evaluation, not only the RMSE have been considered, but also the peak error at the W1 sensor.
| Fig. 5: | Measured and calculated deflections at one sample station (output of DBSID) (a) Low RMS and peak errors and (b) Low RMS but high peak errors |
In most of the DBSID output files from the different sections on the Zanjan-Tabriz highway both RMS and peak errors were high. The time histories of the measured deflections also had irregular shapes. The cause of this was not determined in this study but it could be related to the following reasons:
| • | Equipment Problems: FWD measurements were not taken with good accuracy |
| • | The pavement was distressed at the time of testing and the surface and base condition may have impacted sensor readings |
These output files were excluded from the final conclusion.
CLASSIFICATION OF RESULTS
After the termination of every run and omitting unacceptable values of moduli (the moduli which are not in the defined range) for each deflection input file, the results were classified in two ways:
| • | Scatter charts such as the one shown in Fig. 6 were developed to monitor the variation of modulus in consecutive stations |
| • | Bar charts such as the one shown in Fig. 7 compare the average moduli of different layers |
| Fig. 6: | Moduli changes in different stations. This example result is due to Apron number 3 of Rafsanjan airport (a) Surface (b) Base and (c) Subgard |
| Fig. 7: | Comparison of the average moduli. This example result is due to Apron number 3 of Rafsanjan airport (a) Surface and (b) Base and Subgrade |
Note that the average values in Fig. 7 are obtained from scattered values of Fig. 6.
CONCLUSIONS
FWD data collected from four different test sites were analyzed by three elasto-static backcalculation programs: MODULUS6.0, ELMOD5.0 and EVERCALC5.0 and a dynamic backcalculation program: DBSID. Qualitative comparison of the predicted moduli are presented in the following two sections:
Comparison of Three Elasto-Static Programs
| • | The moduli change in consecutive stations and the average moduli of surface layer, from the outputs of all the three programs show good consistency in all three cases |
| • | ELMOD5.0 overestimates the resilient moduli of the subgrade in some cases |
| • | MODULUS6.0 and EVERCALC5.0 show good consistency in results for all three layers |
| • | As the depth to bedrock is changeable in MODULUS6.0, if the reasonable range of the subgrade moduli is clear, the user can alter the depth to bedrock until the desired values of subgrade moduli are obtained. This is also possible in EVERCALC5.0, but the process is faster in MODULUS6.0 |
| • | Coefficients of Variance (CV), which is computed by dividing the variance by the mean value, for the predicted subgrade moduli of MODULUS6.0 are lower than that of the other programs (not greater than 19%) |
| • | As the backcalculation process often starts from predicting the subgrade modulus, a proper estimate of this value is critical for the rest of the process. The subgrade moduli of MODULUS 6.0 were in the range anticipated from known site conditions and from the CBR data described earlier |
Comparison of the Elasto-Static Analysis of MODULUS6.0 with the Dynamic Analysis of DBSID
| • | Subgrade moduli predicted by DBSID are almost 2.5 times greater than subgrade moduli predicted by MODULUS 6.0. The ratio of the depth of stiff layer from DBSID to MODULUS 6.0 is almost the same value (2.5) |
| • | Base moduli calculated by MODULUS6.0 are almost 1.3 times greater than that of DBSID. So, the base moduli of DBSID are changeable to the base moduli of DBSID by a factor of 0.77 |
| • | DBSID’s surface moduli are often equal or greater (up to 20%) than that of MODULUS 6.0 |
| • | Dynamic backcalculation is a time consuming task, analyzing each station by DBSID takes about 10 min by normal PCS |
So, by considering all the mentioned points, MODULUS6.0 is suggested as the most appropriate program for the cases under study.
ACKNOWLEDGMENT
The authors are grateful to Dr. Emanuel Fernando of Texas A and M University for providing the DBSID program.