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).
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.,
||Measured deflection at sensor i
||Computed deflection at sensor i
||The number of sensors
The computed deflection wic can be expressed as follows
(Uzan et al., 1989):
Xj is the unknown variables, including modulus of elasticity,
poissons 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,
Log (deflectionj) = Aij+
Sji (log Ei)
, ND (ND = No. of deflections)
, 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):
||Radius of contact area
||Modulus of elasticity of the ith layer on the subgrade
||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,
for flexible pavement
for rigid pavement
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,
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,
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 poissons 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,
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 poissons 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):
||Major principal stress from external loading
||Reference stress, often taken equal to atmospheric pressure (0.1
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,
The program terminates when one or more of the following conditions are satisfied
(Everseries Users Guide, 2005):
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
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
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,
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 (Users 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:
||The creep compliance function, defined as D(t) = e (t)/s0
||The total axial strain measured at time t
||The applied constant stress in a creep compliance test
|D0, D1 and m
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:
||Thickness of asphalt concrete layer
||Unit weight of asphalt concrete
||1/E0, representing the elastic response due to the solid
||Creep compliance constant representing the nonlinear viscous response
||Exponent for nonlinear time dependence
Damped-Elastic Material (Users 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
It is noteworthy that in the last item, the HWD test has been performed.
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 Poissons 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
||Default moduli ranges and Poissons Ratio values (William,
|| Relationship between CBR and resilient modulus for the CBR
values less than 25
|| Relationship between CBR and resilient modulus for the CBR
values greater than 25
||Table of typical m values and unit weight for base course
and subgrade materials (Lytton, 1989)
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.
The elasto-static programs terminate when the convergence condition between
deflections is satisfied. The Root Mean Squared Errors (RMSEs)
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
||Measured and calculated deflections at one sample station
(output of DBSID) (a) Low RMS and peak errors and (b) Low RMS but high peak
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
||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
||Bar charts such as the one shown in Fig. 7 compare
the average moduli of different layers
||Moduli changes in different stations. This example result
is due to Apron number 3 of Rafsanjan airport (a) Surface (b) Base and (c)
||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.
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
||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
||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
||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
||DBSIDs 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.
The authors are grateful to Dr. Emanuel Fernando of Texas A and M University
for providing the DBSID program.