The empirical or semi-empirical nature of Witczak and Hirsch models recommends their validation for any database that differs from those used in the calibration process. A number of studies in the open literature present comparisons between these models and actual dynamic modulus data. Results from such studies do not always converge, with some of them reporting sufficiently accurate and robust estimates of ∣E*∣ for use in mechanistic-empirical pavement performance prediction and design, while others report the opposite. Examples of the former are provided by Schwartz (8), who evaluated Witczak η-based model and obtained an R² of 97% in the logarithmic scale, and Robbins (9), who reported Witczak η-based and Hirsch models predictions to generally follow the equality line when compared to measured data. Examples of unacceptable predictions for Witczak η-based model are provided by Birgisson et al. (10), who indicated the need to apply a correction factor that was between 0.668 and 2.740, depending on the method used for measuring binder viscosity, and Shen and Yu (11), who obtained modulus estimates that were between 3 and 7 times measured moduli for high temperatures and between 25% and 50% of measured moduli for low temperatures. Shen and Yu indicated that Witczak G*-based model improved predictions for high and low temperatures, while the opposite happened for intermediate temperatures. They also indicated that Hirsch model over predicts modulus at high temperatures. This over prediction has been also reported by Schwartz (8) for Witczak η-based model and by Al-Khateeb et al. (12) for both Witczak η-based and Hirsch models.

The above mentioned contradictory results, along with the fact that such a comparison effort has not yet been conducted for Spanish mixtures, requires caution when either of Witczak models or the Hirsch one is used for mixtures and testing protocols typically used in Spain. It should be mentioned that a number of the experimental studies presented above propose alternative simplified predictive equations as a result of Witczak and Hirsch models mispredictions (11, 12). Another problem is related to binder rheology characterization which, either in terms of η or G*, is not typically available in routine mixture design or specifications in Spain. The need exists to develop simplified predictive equations or a comprehensive database for such binder mechanical properties.

The objective of the experimental study presented herein is the evaluation of the applicability of Witczak and Hirsch models for asphalt concrete mixtures and testing protocols typically used in Spain. As a result of this evaluation, it is also considered the validation of an alternative simplified model that can be used to predict dynamic modulus of asphalt concrete from routine mixture design and specification parameters.

Dynamic modulus experimental data for the present research come from eight asphalt concrete mixtures (four different gradations and two asphalt contents) whose formulations are presented in Table 1. They represent typical mixtures used as base and binder courses in Spain. The number 20 in mixtures denomination denotes the nominal maximum aggregate size in mm. The letter “G” refers to a particular aggregate gradation that is only used as a base material, under a binder plus a wearing course. The letter “S” refers to another gradation that is typically used in Spain for binder or wearing courses. B60-70 is an unmodified asphalt binder whose present denomination, according to European normative EN 12591, is 50/70, which denotes the penetration interval (mm/10) at 25 °C. BM-3c (PMB 45/80-65, according EN 14023) is highly modified with SBS polymer, and is typically used in binder and wearing courses under heavy traffic. B13-22 (15/25, according EN 13924) is a hard-unmodified bitumen, which is used in the fabrication of high modulus asphalt concrete, “mezcla de alto módulo” (MAM) in Spanish.

Asphalt mixtures were produced in plant, and reheated in the laboratory in order to prepare gyratory specimens (φ100-mm; 180-mm height). Six specimens were prepared for each of the eight mixtures, and they were tested in stress-controlled cyclic compression in order to determine the complex modulus, E*. E* is a complex number whose modulus, ∣E*∣, is the dynamic modulus and whose argument, ϕ, is referred to as phase angle, both being the so called viscoelastic properties under harmonic loading. Minimum testing temperature was –5 °C for all mixtures, while maximum temperature was adapted to binder rheology (50, 55, and 60 °C, respectively, for binders B60-70, BM-3c, and B13-22). A range of frequencies between 0.1 and 20 Hz were applied at temperature increments of no more than 15 °C. No rest period was introduced between two consecutive frequencies. Load level was adjusted for each frequency to produce ±35 µɛ axial deformation, so that linear range was not exceeded. European standard EN 12697-26:2004 “Test Methods for Hot Mix Asphalt - Part 26: Stiffness”, annex D (cyclic tension-compression), was followed except for the type of loading, which was cyclic compression in the present study. The compressive cyclic loading used in this research is also specified by AASHTO T 342-11, whose methodology is very similar to that used in producing databases for Witczak and Hirsch models calibration. The only remarkable difference between the methodology used in the present experimental study and the AASHTO normative is the rest period between two consecutive frequencies, with the later specifying between 2 and 30 minutes. Coring of φ100-mm testing specimens from φ150-mm gyratory specimens is another prescription of AASHTO T 342 which is not specified by European normative and was not followed in this research.

As a first step, the measured dynamic moduli were compared to Witczak and Hirsch models predictions. This comparison involves 352 data points (44, on average, per mix), where each individual point is the mean of the six replicate specimens. Binder rheological properties to use in the models were determined according to AASHTO TP5-98 “Test Method for Determining the Rheological Properties of Asphalt Binder Using a Dynamic Shear Rheometer (DSR)”. The DSR device applies a torque to a thin film of asphalt binder placed between two parallel plates, while measuring the angular strain used to calculate the complex shear modulus, G*. A range of frequencies between 0.1 and 30 Hz were applied at temperatures between 4 and 82 °C, at 6 °C increments. It would have been desirable that the lower bound of the testing interval had been close to or below –5 °C, which was the lower bound of the corresponding temperature interval for asphalt mixture testing. Unfortunately, measuring binder stiffness at such low temperatures is beyond the loading capacity of the DSR used in this research.

An example of dynamic moduli measured for one of the mixtures is presented in Figure 1. Values reflect temperature and frequency susceptibility of asphalt materials. Time-temperature superposition principle can be also deduced from the figure. According to this principle, which is applicable to thermorheologically simple materials like asphalt, a change in temperature is equivalent to a compression/extension of the time scale. This compression/extension means a shift in the logarithmic scale. Consequently, the sets of moduli measured for each temperature can be shifted in the log scale in order to produce a unique set, as reflected in Figure 1. This set defines a curve known as the “dynamic modulus master curve”, and the frequency resulting after the shifting process is known as the “reduced frequency”. The set was fitted by a sigmoidal function in this particular study. The same sigmoidal function, according to equation [1], is used by both CalME and the MEPDG. log ( ∣ E * ∣ ) = δ + α 1 + e β + γ ⋅ log ( f red )

where, ∣E*∣ is dynamic modulus

  f _red is reduced frequency: log (f _red) = log(f) +log (a_T)

  f is frequency (Hz)

 a_T is time-temperature shift factor; 20 °C is the reference temperature, for which a_T =1

δ, α, ϒ, & γ are model parameters

note: log refers to decimal logarithm

Each of the parameters of the sigmoidal function has a distinct physical meaning. δ represents the minimum dynamic modulus of the mixture (10^δ), that is, the lower asymptote of the master curve. δ + α represents the maximum modulus (10^{δ + α}), which is the higher asymptote. ϒ is a location parameter, which determines the horizontal position of the curvature sign change point, while γ represents the maximum slope of the master curve. Sigmoidal function is an analytical curve frequently used in mix/binder viscoelastic characterization within the frame of mechanistic-empirical pavement design. Its robust nature in curve fitting experimental data (both in frequency and time domain) with some noisy makes it also a good choice in pre-smoothing experimental data before fitting a more computational efficient analytical curve (i.e. Prony series) to perform stress-strain analysis of pavement layers using the Finite Element Method (13).

Binder dynamic shear modulus data were also fitted by a master curve, after the application of the corresponding shift factors. An example is presented in Figure 2 for the B60-70 binder, where three sets of points can be observed. Each set corresponds to a different binder aging condition, which is mainly related to oxidation and volatile loss. The “original” condition refers to the binder before mixture fabrication; the “RTFO” condition refers to the binder residue after the artificially accelerated aging process of the Rolling Thin-Film Oven test (AASHTO T 240), conceived to reproduce aging during mixture fabrication; finally, the “PAV” condition refers to the residue after the Pressurized Aging Vessel test (AASHTO R 28), conceived to reproduce long-term aging (5–10 years) under service conditions. Binder mechanical properties used in mixture dynamic modulus predictive models must reflect actual binder aging condition in the mixture. In the present study, mixtures were tested after plant fabrication, so the binder aging condition is the one simulated by the RTFO test. CAM model (14), reflected in equations [2], was used to fit the ∣G*∣ master curve. The same model parameters were used for the three aging conditions, although a factor “k” was introduced that multiplies reduced frequency for RTFO and PAV data. This approach is based on aging effect being equivalent to a frequency shift in the logarithmic scale. The validity of this approach is based on the excellent fit that was achieved for original and PAV data while using a unique set of CAM model parameters, as observed in Figure 2. The same agreement was achieved for the rest of the binders. | G * | = G g · [ 1 + ( f c / f r e d )   log ( 2 ) R ] - W · R log ( 2 ) δ = 90 · W [ 1 + ( f r e d / f c )   log ( 2 ) R ]

where, ∣G*∣ is dynamic shear modulus

 δ is phase angle

 G_g is binder glassy modulus (1 GPa)

  f _red is reduced frequency

R, w, & f _c are model parameters

Time-temperature shift factors obtained for G20 mixture and B60-70 binder are presented, as an example, in Figure 3. Agreement between both experimental sources was expected, and it happened for all temperatures except for the highest one. No definitive explanation was found for this disagreement, that was also observed for the rest of the mixtures except for those using B13-22, the hardest binder. This led authors to believe that the disagreement obtained for the highest temperature might be related to viscoplastic behavior, which is not captured by dynamic modulus. It is important to bear in mind that dynamic modulus can be considered a mechanical property of the asphalt mixture for medium and low temperatures, when linear viscoelasticity hypotheses are reasonably fulfilled. For high temperatures, non-linear viscoplastic behavior is expected, which means that the ratio of peak cyclic stress to peak cyclic strain will not be a constant anymore. Just as an example, Huang et al. (15) showed that confining pressure had a minimal effect on asphalt mixtures tested at 10 °C (linear behavior), while the opposite happened at 54 °C. At the last temperature, dynamic modulus values measured with confining pressures of 100 and 200 kPa were significantly greater than values measured without confining pressure, which indicates that dynamic modulus was not a material property.

Binder DSR minimum testing temperature was 5 °C, so low-temperature shift factors were extrapolated using the well-known Williams-Landel-Ferry equation reflected in Figure 3. The –5 °C factors obtained for the binders this way were in good agreement with values obtained for the asphalt mixtures. This result, together with the fact that all binders tend to reach a limiting modulus of 1 GPa for low temperatures, seem to overcome the DSR testing limitation.

Witczak first predictive equation is based on binder viscosity. This variable was determined by using the approach proposed by Bonaquist et al. (16), presented in equation [3], which is based on the classical Cox-Merz rule. It should be mentioned that binder behaves like a Newtonian fluid only at very high temperatures (17). For pavement service temperatures, result from equation [3] must be considered only an indication of binder stiffness rather than an intrinsic mechanical property. η = | G * | ω · ( 1 sin ( δ ) ) 3.63922 + 0.13137. ω - 0.0009 ω 2

where, η is binder viscosity in Pa·s

 ∣G*∣ is dynamic shear modulus in Pa

 ω is angular frequency (rad/s)

 δ is phase angle

The approach for the formulation of the new model is based on the relationships existing between parameters of the sigmoidal master curve (equation [1]) and mixture volumetric and gradation properties as well as binder rheological variables. There is a general agreement among researchers that δ and α parameters of the sigmoidal function are determined by volumetric and gradation properties of the mixtures rather than by binder rheology (8). In the upper limit of the reduced frequency (low temperatures), all binders become glassy materials with a modulus around 1 GPa, i.e., binder type does not make a difference. In the lower limit (high temperatures), binder stiffness becomes so low that aggregate skeleton determines asphalt mixture stiffness. Nonetheless, the efforts to relate δ and α parameters to mixture volumetric and gradation variables have not been very successful to date. The only variable for which consensus exists in the way it influences dynamic modulus is VMA, the voids in the mineral aggregate skeleton (6, 8, 12), which is included in Hirsch model but not in Witczak models. Unfortunately, the range of this variable was very narrow in this study, and such influence could not be directly determined. In order to overcome this limitation, results from this study were combined with those obtained by Shen and Yu (11) from seven mixtures typically used in Washington State. The sum of δ and α (log of maximum mixture modulus) is presented in Figure 5 for both experimental sources. The best fit equation presented in the figure can be used to extrapolate results from this study to mixtures with higher VMA. Extrapolation should not go much further than the actual VMA range tested in this research, which should be enough for most AC mixtures used in Spain. From Figure 5, Hirsch model understates the influence of VMA on dynamic modulus, which is in line with results presented above.

The sum δ + α varied from 4.43 to 4.55 among the different mixtures of this study, as reflected in Figure 5. In this case, the master curve of each mixture was determined independently from the others. The global R² of the master curve calibration was 99.93% in the logarithmic scale, which indicates an almost perfect fit, as can be seen in Figure 1. A second calibration was conducted where δ and α were forced to be equal for all the mixtures. The global R² of the new master curve calibration was still 99.82%, indicating that a single set of these two parameters can be used for all mixtures without losing accuracy. ϒ parameter was found to depend on binder stiffness at reference temperature, whether it was formulated in terms of ∣G*∣ or η. This result was expected, as it has been pointed out by Schwartz (8) and reflected by Witczak models. But ϒ was also found to depend on binder content, whether it was formulated in terms of VFA, VA, or effective binder volume (VMA-VA). This dependency is not acknowledged by Witczak models, and may be a reason for the low sensitivity of their predictions with respect to mixture volumetric properties. Finally, Υ parameter was found to depend on β. Again, this dependency is not acknowledged by Witczak η-based model, which makes Υ= -0.31335 in any case. The pattern described above was used to formulate a new predictive model that was calibrated from dynamic modulus data measured for the eight mixtures included in Table 1. The model is presented in equations [4], and Figure 6 indicates the model's goodness of fit. R² of the model was 98.7% in the log scale (96% in arithmetic scale), and MSE was 12.6% in the arithmetic scale for the model based on G* and 12.5% for the model based on η. δ = 1 .572

where, G_ref is binder dynamic shear modulus (MPa) at 20 °C (T_ref) and 10 rad/s

 η_ref is binder steady state viscosity (cP) at 20 °C (T_ref)

VA is voids in asphalt (expressed as ratio)

note: log refers to decimal logarithm

One limitation of the suggested model was related to using only three binder types in the calibration. Specifically, temperature susceptibility was relatively reduced in the three cases, with VTS (viscosity temperature susceptibility) parameter ranging from –2.85 to –2.75, which is a very narrow interval. VTS is the slope of the well-known relationship that exists between double logarithm of viscosity (cP) and logarithm of temperature (in Rankine). In order to overcome this limitation, an additional term was calibrated on the basis of a database generated by Onofre (19) at Universidade Federal do Ceará. Two asphalt binders with seven different types of modifiers were used by Onofre in his research. This way, the VTS varied from –4.25 to –2.90. As expected, temperature susceptibility influenced the Υ parameter of the master curve, with higher temperature susceptibility corresponding to higher Υ in absolute value. This influence could be modeled by adding the term “–0.1536·(∣VTS∣–2.824)” to the Υ expression in equations [4]. This term can be used to extrapolate results from this study to mixtures with binders whose VTS parameter is below –2.85.

Binder rheological data required by the new equation is very limited when compared to Witczak and Hirsch models, since binder stiffness must be only defined for a single temperature. But still, this information will not be available in daily practice of pavement design in Spain, since consistency prescriptions for binders are based on classical penetration and ring and ball tests rather than on DSR or viscosity testing.

A simplified approach is here presented in order to estimate η_ref and VTS from available data. First, it is assumed that the ring and ball softening point corresponds to a viscosity of 13000 Poises, which is a hypothesis based on research conducted by Shell Oil and later confirmed by Mirza and Witczak (20). Secondly, equation [5] is used to estimate log(η_ref) from binder penetration. This empirical equation was calibrated (R² = 97.5%; MSE = 1.17%) on the basis of ten binders typically used in Spain, that were tested within the frame of a research project focused on binder rheology (21). Once viscosity is known for two different temperatures, the VTS parameter can be determined. log ( η ref ) = 12.91 – 2.412 ⋅ log ( p )

where, η_ref is binder steady state viscosity (cP) at 20 °C

 P is binder penetration (mm/10) at 25 °C

note: log refers to decimal logarithm

Predictions of the new equation were compared to dynamic modulus data measured for five mixtures coming from the CEDEX accelerated pavement testing facility. Mixtures had been tested between -15 and 40 °C. The simplified approach described above was used to estimate η_ref and VTS from binders penetration and ring and ball softening point. The RTFO test residue data were used, since asphalt mixtures were tested after fabrication (no additional aging had taken place). Model's goodness of fit is appreciated in Figure 7. R² of the prediction was 97.9% in the log scale (94.1% in the arithmetic scale) and MSE was 14.4% in the arithmetic scale. No bias is appreciated in the error. These statistics are similar to those obtained in the calibration, which validates the suggested predictive equation.