Research, Development, and Technology

Turner-Fairbank Highway Research Center

6300 Georgetown Pike

McLean, VA 22101-2296

Modeling of Hot-Mix Asphalt Compaction:

A Thermodynamics-Based Compressible 

Viscoelastic Model

Publication No. FHWA-HRT-10-065 December 2010
FOREWORD

This report details the activities and findings of research on mathematical modeling of compaction 

of hot-mix asphalt (HMA) pavements. The primary objective of this project was to develop a

material model in a general mechanical setting to study the mechanics associated with asphalt

mix compaction. The researchers developed a viscoelastic compressible fluidlike model, which 

was then adapted for simulating compaction of HMA employing the finite-element (FE) method. 

The model was developed within a thermodynamic framework that was capable of demonstrating

the nonlinear dissipative response associated with the asphalt mix. This work should lead to further

development of simulation models capable of capturing the mechanics of the compaction processes

both in the laboratory and in the field. 

The FE method was employed to simulate laboratory and field compaction under various loading 

and boundary conditions. The FE results agree well with the data obtained from various laboratory

and field compaction projects. The field compaction results validate the applicability of the model

to predict compaction in highway projects.

This study can be used by practitioners and researchers to design mixtures with desirable

compaction characteristics and improved performance. The material presented here is the final 

report for the project.

Jorge E. Pagán-Ortiz

Director, Office of Infrastructure

Research and Development

Notice

This document is disseminated under the sponsorship of the U.S. Department of Transportation 

in the interest of information exchange. The U.S. Government assumes no liability for the use 

of the information contained in this document. This report does not constitute a standard, 

specification, or regulation.

The U.S. Government does not endorse products or manufacturers. Trademarks or 

manufacturers’ names appear in this report only because they are considered essential to the 

objective of the document.

Quality Assurance Statement

The Federal Highway Administration (FHWA) provides high-quality information to serve 

Government, industry, and the public in a manner that promotes public understanding. Standards 

and policies are used to ensure and maximize the quality, objectivity, utility, and integrity of its 

information. FHWA periodically reviews quality issues and adjusts its programs and processes to 

ensure continuous quality improvement.
TECHNICAL REPORT DOCUMENTATION PAGE

1. Report No.

FHWA-HRT-10-065

2. Government Accession No. 3. Recipient Catalog No.

4. Title and Subtitle

Modeling of Hot-Mix Asphalt Compaction: A Thermodynamics-Based 

Compressible Viscoelastic Model

5. Report Date

December 2010

6. Performing Organization Code

C-00053

7. Author(s)

Eyad Masad, Saradhi Koneru, Tom Scarpas, Emad Kassem, and 

K.R. Rajagopal

8. Performing Organization Report No.

DTFH61-07-C-00053-F-10-01

9. Performing Organization Name and Address

Texas Transportation Institute

The Texas A&M University System

College Station, TX 77845-3135

10. Work Unit No. (TRAIS)

11. Contract or Grant No.

DTFH61-07-C-00053

12. Sponsoring Agency Name and Address

Office of Acquisition Management

Federal Highway Administration

1200 New Jersey Avenue, SE

Washington, DC 20590-0001

13. Type of Report and Period Covered

Final Report

November 2007–May 2010

14. Sponsoring Agency Code

HAAM-30

15. Supplementary Notes

This project was completed in coordination with Eric Weaver from FHWA by researchers at Texas A&M University and 

researchers at the Section of Structural Mechanics, Civil Engineering and Geosciences, Delft University of Technology, 

The Netherlands.

16. Abstract 

Compaction is the process of reducing the volume of hot-mix asphalt (HMA) by the application of external forces. As a 

result of compaction, the volume of air voids decreases, aggregate interlock increases, and interparticle friction increases.

The quality of field compaction of HMA is one of the most important elements influencing asphalt pavement performance.

Poor compaction has been associated with asphalt bleeding in hot weather, moisture damage, excessive aging and 

associated cracking, and premature permanent deformation. 

This study was conducted to develop a model within the context of a thermomechanical framework for the compaction 

of asphalt mixtures. The asphalt mixture was modeled as a nonlinear compressible material exhibiting time-dependent 

properties. A numerical scheme based on finite elements was employed to solve the equations governing compaction 

mechanisms. The material model was implemented in the Computer Aided Pavement Analysis (CAPA-3D) finite-element 

(FE) package. Due to the difficulty of conducting tests on the mixture at the compaction temperature, a procedure was 

developed to determine the model’s parameters from the analysis of the Superpave® gyratory compaction curves. A 

number of mixtures were compacted in the Superpave® gyratory compactor using an angle of 1.25 degrees in order to 

determine the model’s parameters. Consequently, the model was used to predict the compaction curves of mixtures 

compacted using a 2-degree angle of gyration. The model compared reasonably well with the compaction curves. FE 

simulations of the compaction of several pavement sections were conducted in this study. The results demonstrated the

potential of the material model to represent asphalt mixture field compaction.

The developed model is a useful tool for simulating the compaction of asphalt mixtures under laboratory and field 

conditions. In addition, it can be used to determine the influence of various material properties and mixture designs on the 

model’s parameters and mixture compactability. 

17. Key Words

Compaction, Hot-mix asphalt, Mixture design, Constitutive model, 

Finite element, Simulation, Viscoelastic, Compressible

18. Distribution Statement

No restrictions. This document is available to the

public through the National Technical Information 

Service, Alexandria, VA 22312.

19. Security Classification (of this report)

Unclassified

20. Security Classification (of this page) 

Unclassified

21. No. of pages

110

22. Price

Form DOT F 1700.7 (8-72) Reproduction of completed page authorized
ii

SI* (MODERN METRIC) CONVERSION FACTORS 

APPROXIMATE CONVERSIONS TO SI UNITS

Symbol When You Know Multiply By To Find Symbol 

LENGTH

in inches 25.4 millimeters mm 

ft feet 0.305 meters m 

yd yards 0.914 meters m 

mi miles 1.61 kilometers km 

AREA 

in2 square inches 645.2 square millimeters mm2

ft2 square feet 0.093 square meters m2

yd2 square yard 0.836 square meters m2

ac acres 0.405 hectares ha

mi2 square miles 2.59 square kilometers km2

VOLUME 

fl oz fluid ounces 29.57 milliliters mL 

gal gallons 3.785 liters L 

ft3 cubic feet 0.028 cubic meters m3 

yd3 cubic yards 0.765 cubic meters m3 

NOTE: volumes greater than 1000 L shall be shown in m3

MASS oz ounces 28.35 grams g

lb pounds 0.454 kilograms kg

T short tons (2000 lb) 0.907 megagrams (or "metric ton") Mg (or "t") 

TEMPERATURE (exact degrees) o

F Fahrenheit 5 (F-32)/9 Celsius oC 

or (F-32)/1.8 

ILLUMINATION 

fc foot-candles 10.76 lux lx

fl foot-Lamberts 3.426 candela/m2 cd/m2

FORCE and PRESSURE or STRESS 

lbf poundforce 4.45 newtons N 

lbf/in2 poundforce per square inch 6.89 kilopascals kPa 

APPROXIMATE CONVERSIONS FROM SI UNITS 

Symbol When You Know Multiply By To Find Symbol 

LENGTH mm millimeters 0.039 inches in 

m meters 3.28 feet ft 

m meters 1.09 yards yd 

km kilometers 0.621 miles mi 

AREA mm2 square millimeters 0.0016 square inches in2 

m2 square meters 10.764 square feet ft2 

m2 square meters 1.195 square yards yd2 

ha hectares 2.47 acres ac

km2 square kilometers 0.386 square miles mi2 

VOLUME 

mL milliliters 0.034 fluid ounces fl oz 

L liters 0.264 gallons gal 

m3 cubic meters 35.314 cubic feet ft3 

m3 cubic meters 1.307 cubic yards yd3 

MASS 

g grams 0.035 ounces oz

kg kilograms 2.202 pounds lb

Mg (or "t") megagrams (or "metric ton") 1.103 short tons (2000 lb) T 

TEMPERATURE (exact degrees) o

C Celsius 1.8C+32 Fahrenheit oF 

ILLUMINATION 

lx lux 0.0929 foot-candles fc 

cd/m2 candela/m2 0.2919 foot-Lamberts fl

FORCE and PRESSURE or STRESS 

N newtons 0.225 poundforce lbf 

kPa kilopascals 0.145 poundforce per square inch lbf/in2

*SI is the symbol for the International System of Units. Appropriate rounding should be made to comply with Section 4 of ASTM E380. 

(Revised March 2003)
iii

TABLE OF CONTENTS

CHAPTER 1. INTRODUCTION.................................................................................................1

PROBLEM STATEMENT .....................................................................................................1

OBJECTIVES OF RESEARCH PROJECT.........................................................................1

SCOPE AND LIMITATIONS................................................................................................2

OUTLINE OF RESEARCH REPORT..................................................................................2

NOTATIONS ...........................................................................................................................3

CHAPTER 2. REVIEW OF THE LITERATURE.....................................................................5

HMA COMPACTION.............................................................................................................5

FACTORS AFFECTING COMPACTION...........................................................................6

LABORATORY VERSUS FIELD COMPACTION............................................................6

MODELING HMA RESPONSE............................................................................................7

COMPACTION MODELS.....................................................................................................8

INTELLIGENT COMPACTION ..........................................................................................9

CHAPTER 3. A THERMOMECHANICAL MODEL FOR ASPHALT CONCRETE........11

THERMODYNAMICS OF DISSIPATIVE MATERIALS...............................................11

CONSTITUTIVE MODEL...................................................................................................12

CONSTITUTIVE EXPRESSION FOR STRESS...............................................................13

THE EVOLUTION EQUATION FOR BBκκpp(tt)......................................................................13

ANALYTICAL CALCULATIONS USING A SEMI-INVERSE APPROACH..............13

ONE-DIMENSIONAL CREEP DEFORMATION............................................................13

STRESS RELAXATION.......................................................................................................17

CHAPTER 4. NUMERICAL IMPLEMENTATION IN THE FE METHOD ......................21

VERIFICATION OF FE IMPLEMENTATION................................................................21

CONSTANT SHEAR.............................................................................................................23

CONSTANT SHEAR RATE ................................................................................................25

CHAPTER 5. FE SIMULATION OF SGC COMPACTION..................................................29

RELATIONSHIPS OF MODEL-INPUT PARAMETERS TO MATERIAL 

PROPERTIES........................................................................................................................29

Model Verification and Calibration...................................................................................30

PARAMETRIC ANALYSIS OF MODEL PARAMETERS .............................................32

SENSITIVITY OF COMPACTION TO ANGLE OF GYRATION ................................37

MODEL CALIBRATION AND VERIFICATION............................................................39

CHAPTER 6. FE SIMULATION OF FIELD COMPACTION..............................................47

INTRODUCTION..................................................................................................................47

FE MODEL ............................................................................................................................48

FE Mesh Used for Simulation............................................................................................48

Boundary Conditions for Simulating Field Compaction ...................................................50

Response to Moving Load .................................................................................................51

Implementing a Contact-Area Algorithm..........................................................................54

Normal Stresses Predicted by Model.................................................................................58

Interface Elements .............................................................................................................59
iv

Use of Impedance Layers...................................................................................................61

Effect of Base Stiffness......................................................................................................64

Compaction Equipment Characteristics.............................................................................67

CHAPTER 7. SIMULATION OF FIELD COMPACTION....................................................71

SENSITIVITY ANALYSIS ..................................................................................................71

Observations from the Sensitivity Analysis.......................................................................72

Parametric Analysis...........................................................................................................73

Correlation with Laboratory Compaction..........................................................................77

ANALYZING LONGITUDINAL JOINTS.........................................................................77

VERIFICATION OF SIMULATION MODEL..................................................................78

Highway SH-21 .................................................................................................................78

Highway US-87 .................................................................................................................82

Highway US-259 ...............................................................................................................84

Observations from Field Compaction Simulations for US-87 and US-259 ......................86

CHAPTER 8. SUMMARY AND CONCLUSIONS..................................................................91

FE IMPLEMENTATION .....................................................................................................91

SIMULATION OF GYRATORY COMPACTION ...........................................................92

SIMULATION OF FIELD COMPACTION ......................................................................92

RECOMMENDATIONS FOR FUTURE INVESTIGATIONS........................................93

REFERENCES.............................................................................................................................97
v

LIST OF FIGURES

Figure 1. Equation. Gradient and divergence operators with respect to reference configuration ...4

Figure 2. Equation. Gradient and divergence operators with respect to current configuration.......4

Figure 3. Illustration. Evolution of natural configurations associated with microstructural 

transformations resulting from the material response to an external stimulus ..............................11

Figure 4. Equation. Energy-dissipation inequality ........................................................................11

Figure 5. Equation. Reduced energy-dissipation inequality ..........................................................12

Figure 6. Equation. Helmholtz potential for mixture.....................................................................12

Figure 7. Equation. Rate of dissipation function ...........................................................................12

Figure 8. Equation. Shear-modulus function .................................................................................12

Figure 9. Equation. Viscosity function ..........................................................................................12

Figure 10. Equation. Constitutive equation for stress....................................................................13

Figure 11. Equation. Evolution equation .......................................................................................13

Figure 12. Illustration. Schematic for the one-dimensional compression problem.......................14

Figure 13. Chart. Creep solution for the material model for HMA—creep load...........................14

Figure 14. Chart. Creep solution for the material model for HMA—Almansi strain....................15

Figure 15. Chart. Total stretch calculated from the solution to the one-dimensional 

creep problem.................................................................................................................................15

Figure 16. Chart. Evolution of natural configuration calculated from the solution to the 

one-dimensional creep problem.....................................................................................................16

Figure 17. Chart. Stored energy calculated from the solution to the one-dimensional 

creep problem.................................................................................................................................16

Figure 18. Chart. Rate of dissipation calculated from the solution to the one-dimensional 

creep problem.................................................................................................................................17

Figure 19. Illustration. Schematic for the one-dimensional stress-relaxation problem.................18

Figure 20. Chart. Applied strain calculated from the solution to the one-dimensional 

stress-relaxation problem...............................................................................................................18

Figure 21. Chart. Stress relaxation calculated from the solution to the one-dimensional 

stress-relaxation problem...............................................................................................................19

Figure 22. Chart. Stored energy calculated from the solution to the one-dimensional 

stress-relaxation problem...............................................................................................................19

Figure 23. Chart. Rate of dissipation calculated from the solution to the one-dimensional 

stress-relaxation problem...............................................................................................................20

Figure 24. Chart. FE solution to the model response to applied constant stress—creep load.......22

Figure 25. Chart. FE solution to the model response to applied constant stress—total stretch.....22

Figure 26. Chart. FE solution to the model response to applied constant strain 

(compressive 0.05)—applied strain ...............................................................................................22

Figure 27. Chart. FE solution to the model response to applied constant strain 

(compressive 0.05)—normal stress response.................................................................................23

Figure 28. Illustration. Schematic of the constant shear loading applied to a unit cube ...............23

Figure 29. Chart. Shear stress (T12) observed in response to constant shear loading ....................24

Figure 30. Chart. Comparison of the first normal stress (T11–T22) response to constant 

shear loading ..................................................................................................................................24

Figure 31. Chart. Comparison of the second normal stress (T22–T33) response to constant 

shear loading ..................................................................................................................................25
vi

Figure 32. Illustration. Schematic of the constant shear rate applied to a unit cube .....................25

Figure 33. Chart. Shear stress (T12) observed in response to constant shear rate loading .............26

Figure 34. Chart. Shear stress (T12) as a function of the shear rate (using the model 

parameters in table 1).....................................................................................................................26

Figure 35. Chart. Comparison of the first normal stress (T11–T22) response to 

constant shear rate loading.............................................................................................................27

Figure 36. Chart. Comparison of the second normal stress (T22–T33) response to 

constant shear rate loading.............................................................................................................27

Figure 37. Photo. Superpave® gyratory compactor .......................................................................29

Figure 38. Photo. Static steel-wheel roller.....................................................................................30

Figure 39. Illustration. FE mesh used in modeling the SGC .........................................................32

Figure 40. Chart. Analysis of the sensitivity of compaction to μμ̂..................................................33

Figure 41. Chart. Analysis of the sensitivity of compaction to n1 .................................................34

Figure 42. Chart. Analysis of the sensitivity of compaction to ηη̂..................................................35

Figure 43. Chart. Analysis of the sensitivity of compaction to λ2 .................................................35

Figure 44. Chart. Analysis of the sensitivity of compaction to q2 .................................................36

Figure 45. Equation. Modified shear-modulus function................................................................36

Figure 46. Equation. Modified viscosity function .........................................................................36

Figure 47. Chart. Illustration of the relationship of the model’s parameters to the 

compaction process........................................................................................................................37

Figure 48. Chart. Influence of angle of gyration on the compaction curve ...................................38

Figure 49. Chart. Maximum shear stress at the top of the specimen for a gyration angle 

of 1.25 degrees...............................................................................................................................39

Figure 50. Chart. Maximum shear stress at the top of the specimen for a gyration angle 

of 2.0 degrees.................................................................................................................................39

Figure 51. Chart. Fitting of the compaction data at 1.25 degrees for project IH-35......................40

Figure 52. Chart. Fitting of the compaction data at 1.25 degrees for project US-259...................41

Figure 53. Chart. Fitting of the compaction data at 1.25 degrees for project SH-36.....................41

Figure 54. Chart. Fitting of the compaction data at 1.25 degrees for project SH-21.....................42

Figure 55. Chart. Fitting of the compaction data at 1.25 degrees for project US-87.....................42

Figure 56. Chart. Prediction of the compaction data at 2.0 degrees for project IH-35..................43

Figure 57. Chart. Prediction of the compaction data at 2.0 degrees for project US-259...............43

Figure 58. Chart. Prediction of the compaction data at 2.0 degrees for project SH-36.................44

Figure 59. Chart. Prediction of the compaction data at 2.0 degrees for project SH-21.................44

Figure 60. Chart. Prediction of the compaction data at 2.0 degrees for project US-87.................45

Figure 61. Illustration. Pavement structure typically employed for studying field compaction....48

Figure 62. Illustration. Sectional view of the FE mesh used for setting up field compaction 

simulations.....................................................................................................................................49

Figure 63. Illustration. Schematic diagram illustrating the edges of the lane that correspond 

to fixed and free edges of the mesh in figure 62............................................................................50

Figure 64. Chart. Typical displacement curve for a node under the cylindrical load for a 

cycle with forward and return passes followed by a forward pass with the load removed ...........51

Figure 65. Chart. Vertical displacement observed at point P on the path of the roller..................52

Figure 66. Chart. Roller location from observation point P ..........................................................52

Figure 67. Chart. Deflection at the node of interest when a load is applied for a short duration

and then removed...........................................................................................................................53
vii

Figure 68. Chart. Permanent deformation prediction of the model in multiple-pass loading .......53

Figure 69. Illustration. Roller contact geometry for static indentation and during motion ...........54

Figure 70. Chart. Indentation of a cylindrical roller into pavement ..............................................54

Figure 71. Chart. Change in nodal reaction forces as the load is applied over a smaller area ......55

Figure 72. Chart. Comparison of vertical displacement for the two different patterns 

of loading .......................................................................................................................................56

Figure 73. Chart. Comparison of the normal stress distribution at the pavement’s 

top surface due to different loading patterns..................................................................................57

Figure 74. Chart. Comparison of the shear-stress (XY) distribution of the pavement’s 

top surface due to different loading patterns (on the plane of symmetry).....................................57

Figure 75. Equation. Mean contact pressure..................................................................................58

Figure 76. Equation. Normal stress in elastic medium..................................................................58

Figure 77. Chart. Comparison of the normal stress predicted by the model to the normal stress 

predicted for an elastic medium (for the same mean contact pressure, Pm)...................................59

Figure 78. Chart. Comparison of pavement vertical response as the normal stiffness of the 

interface elements is lowered by an order of magnitude from 1,451 ksi (10,000 MPa) 

to 145 ksi (1,000 MPa)...................................................................................................................60

Figure 79. Chart. Comparison of pavement response in the rolling direction as the normal 

stiffness of the interface elements is lowered by an order of magnitude from 

1,450 ksi (10,000 MPa) to 145 ksi (1,000 MPa)............................................................................60

Figure 80. Chart. Comparison of pavement response to increasing the shear stiffness 

in the x-direction............................................................................................................................61

Figure 81. Chart. Comparison of pavement response to increasing the shear stiffness 

in the z-direction ............................................................................................................................61

Figure 82. Chart. Comparison of the dampening effect provided by impedance layers 

surrounding the structure laterally in the x-direction.....................................................................62

Figure 83. Chart. Comparison of the dampening effect provided by impedance layers 

surrounding the structure laterally in the y-direction.....................................................................63

Figure 84. Chart. Comparison of the horizontal dampening effect provided by impedance 

layers surrounding the top layer (pavement and old asphalt) laterally in the x-direction..............63

Figure 85. Chart. Comparison of the vertical dampening effect provided by impedance 

layers surrounding the top layer laterally in the x-direction..........................................................64

Figure 86. Chart. Comparison of the volumetric component of the viscous-evolution gradient 

for soft (72.5 ksi (500 MPa)) versus stiff (290 ksi (2,000 MPa)) bases ........................................65

Figure 87. Chart. Comparison of the volumetric component of the viscous-evolution gradient 

at two base-stiffness moduli of interest .........................................................................................65

Figure 88. Chart. Comparison of the effect on the x-displacement of a node in the roller path 

as the base stiffness is varied from 72.5 ksi (500 MPa) (soft base) to 290 ksi (2,000 MPa) 

(stiff base) ......................................................................................................................................66

Figure 89. Chart. Comparison of the effect on the y-displacement of a node in the roller path 

as the base stiffness is varied from 72.5 ksi (500 MPa) (soft base) to 290 ksi (2,000 MPa) 

(stiff base) ......................................................................................................................................66

Figure 90. Chart. Comparison of the deflection for two base-stiffness moduli of interest............67

Figure 91. Chart. Compaction over a sequence of passes as the amplitude of vibration increases...68

Figure 92. Chart. Compaction over a sequence of passes at different frequencies .......................69

Figure 93. Chart. Material response to change in dead load carried by each roller.......................69
viii

Figure 94. Chart. Field compaction response at a constant frequency over multiple passes 

on a point........................................................................................................................................70

Figure 95. Chart. Evolution of the volumetric viscous gradient with a change in values of 

individual parameters μμ̂, ηη̂, λ1, and λ2 ...........................................................................................71

Figure 96. Chart. Evolution of the volumetric viscous gradient with a change in values of 

individual parameters n1, n2, q1, and q2..........................................................................................72

Figure 97. Chart. Regions of influence of model parameters in gyratory compaction..................73

Figure 98. Chart. Evolution of the volumetric viscous gradient with a change in λ1 ....................74

Figure 99. Chart. Evolution of the volumetric viscous gradient with a change in q1 ....................74

Figure 100. Chart. Evolution of the volumetric viscous gradient with a change in n2 ..................75

Figure 101. Chart. Evolution of the volumetric viscous gradient with a change in λ2 ..................75

Figure 102. Chart. Evolution of the volumetric viscous gradient with a change in q2 ..................76

Figure 103. Chart. Evolution of the volumetric viscous gradient with a change in n1 ..................76

Figure 104. Chart. Plot representing the final compacted state of the material along the 

width of the pavement....................................................................................................................78

Figure 105. Illustration. Schematic of a roller on a material with three locations for density 

measurements.................................................................................................................................79

Figure 106. Chart. Measurements of the %AV in the asphalt mix................................................80

Figure 107. Chart. Measurements of change in %AV in the asphalt mix .....................................80

Figure 108. Chart. Measurements and modeling results of %AV at point A of the pavement 

locations shown in figure 105........................................................................................................80

Figure 109. Chart. Measurements and modeling results of %AV at point B of the pavement 

locations shown in figure 105........................................................................................................81

Figure 110. Chart. Measurements and modeling results of %AV at point C of the pavement 

locations shown in figure 105........................................................................................................81

Figure 111. Illustration. Pavement structure for the US-87 project...............................................83

Figure 112. Chart. Schematic for the rolling patterns for the US-87 project.................................83

Figure 113. Illustration. Pavement structure for the US-259 project.............................................84

Figure 114. Chart. Schematic for the rolling patterns for the US-259 project...............................85

Figure 115. Chart. Comparison of the total percent compaction from simulations with the 

general trend of the %AV measured at the end of the field compaction process for US-87.........87

Figure 116. Chart. Total percent compaction from simulations compared to the 

general trend of the %AV measured at the end of the field compaction process for US-259.......87

Figure 117. Chart. Comparison of prediction of percent compaction per roller pass for 

US-87 and US-259.........................................................................................................................88

Figure 118. Chart. Prediction of percent compaction per roller pass across the material for 

US-87 (cores taken at four locations) ............................................................................................89

Figure 119. Chart. Prediction of percent compaction per roller pass across the material for 

US-259 (cores taken at four locations) ..........................................................................................90

Figure 120. Illustration. Micromechanical response of an asphalt mix.........................................94

Figure 121. Chart. Representation of tasks involved in modeling asphalt mix compaction .........95
ix

LIST OF TABLES

Table 1. Parameters employed for the study of the material response when subject to 

one-dimensional constrained compression ....................................................................................21

Table 2. Model parameters used in the parametric study ..............................................................33

Table 3. Model parameters obtained from compaction data..........................................................40

Table 4. Typical material properties of pavement structural layers...............................................50

Table 5. Parameters used to study the typical material response during field compaction ...........51

Table 6. Typical roller-speed ranges..............................................................................................55

Table 7. Material parameters used to investigate the contact algorithm, interface effects, 

impedance element effects, and the effect of base stiffness ..........................................................57

Table 8. Parameters employed for the sensitivity study of the material........................................72

Table 9. Material properties used for the SH-21 project ...............................................................79

Table 10. Model parameters used for projects SH-21, US-87, and US-259..................................79

Table 11. Summary of mixture designs.........................................................................................82

Table 12. Summary of properties of mixture constituents.............................................................82

Table 13. Material properties used for the US-87 project .............................................................83

Table 14. Rolling pattern and %AV measured in the field for the US-87 project.........................84

Table 15. Material properties used for the US-259 project ...........................................................85

Table 16. Rolling pattern and %AV measured in the field for the US-259 project.......................86

Table 17. Change in measured %AV in the field and calculated percent compaction for 

the US-87 project...........................................................................................................................86

Table 18. Change in measured %AV in the field and calculated percent compaction for 

the US-259 project.........................................................................................................................86
x

LIST OF ABBREVIATIONS

%AV Percent air void

CAPA-3D Computer Aided Pavement Analysis

FE Finite-element

FHWA Federal Highway Administration

HMA Hot-mix asphalt

IC Intelligent compaction

PG Superpave® performance grading

SGC Superpave® gyratory compactor

VMA Void in the mineral aggregate
1

CHAPTER 1. INTRODUCTION

Compaction is a process that reduces the volume of a mixture of hot asphalt binder, aggregates, filler

materials, and air voids by the application of external forces to form a dense mass. This densification

causes an increase in the unit weight of the material and improves the aggregate interlock. The 

goals of compacting an asphalt pavement are to achieve an optimum air-void content, to provide 

a smooth riding surface, and to increase the load-bearing capacity of the material under construction.

(1)

Improper compaction generally leads to poor performance from the asphalt pavement in spite of 

all other desirable mixture-design characteristics being met. This can lead to premature, irreparable

damage to the built-in infrastructure. Such damage typically manifests in the form of rutting, 

permanent deformation, cracking, and moisture damage. At the elevated temperatures at which

compaction is typically performed, hot-mix asphalt (HMA) is a nonlinear fluidlike material. Toward

the end of the compaction process, the material slowly transforms into a nonlinear, highly viscous

viscoelastic fluid. The challenge to modeling compaction is choosing an appropriate material

model capable of representing such a transformation. 

PROBLEM STATEMENT

Compaction is an important process that has significant impact on the performance of asphalt 

pavements. Poor compaction of asphalt mixtures makes them susceptible to permeation of oxygen

and moisture. As a result, these mixtures become more prone to failure mechanisms such as 

cracking and permanent deformation. 

This study targeted the development of a theoretical and computational platform that can be 

adapted for the simulation of the typical compaction processes: laboratory compaction using a 

gyratory compactor and field compaction using a sequence of rolling passes. There are no effective

methods currently available that are capable of addressing the simulation of these two processes. 

Efforts were undertaken in this study to fill this void and to provide a reasonable starting point to 

modeling compaction of HMA. The model developed and employed is an isothermal model 

that can be adapted in the future to account for nonisothermal phenomena and is capable of

exhibiting a compressible viscoelastic fluidlike response when stimulated by external mechanical 

means. The model was implemented in the finite-element (FE) method to simulate the conditions 

most relevant to laboratory and field compaction processes. 

This study aimed to provide a better understanding of the modeling aspects of HMA compaction 

and has the potential to lead to the development of more general integrated approaches to performance

prediction and design improvement of pavements such as intelligent compaction (IC) systems. Once

an appropriate computer-simulation environment is in place, users will be able to quickly adapt it 

to provide feedback to the controller of the compaction process in a highway construction project.

OBJECTIVES OF RESEARCH PROJECT

The primary objectives of this research project are as follows:

Derive a three-dimensional constitutive model for highly compressible viscoelastic materials. 

Implement the constitutive model in finite elements. 
2

Develop a method to determine the model’s parameters. 

Utilize the model to describe the compaction of five asphalt mixtures under the forces 

applied in the Superpave® gyratory compactor (SGC).

Utilize the model implemented in finite elements to simulate the field compaction of 

asphalt mixtures.

SCOPE AND LIMITATIONS

The modeling of asphalt mixtures is rather extensive and complex. It would be too ambitious to 

attempt to solve the issues of the whole field in a single research project. Therefore, restrictions 

were applied to limit the scope of the current research. This project is restricted to deducing material

parameters for a single mixture type, HMA. A second restriction is that researchers have assumed

isothermal conditions during laboratory and field compaction processes. The third restriction is 

that simulations for rolling compaction have only considered static steel rolling. Despite the

limitations, the modeling and simulation framework developed here does work for HMA mixtures

for static steel rollers. Researchers expect that the system can be adapted to work for other 

mixture-roller combinations as long as they have comparable characteristics.

OUTLINE OF RESEARCH REPORT

The following is an outline of the remainder of this report:

Chapter 2. Review of the Literature.

o A review of compaction-modeling research is presented.

o A number of relevant studies on asphalt mixture compaction are presented. 

Chapter 3. A Thermomechanical Model for Asphalt Concrete.

o Details of the constitutive model used for capturing HMA response during compaction

(utilizing a framework of multiple natural configurations) are presented.

o The equations of the model are formulated to make them amenable to implementation 

in finite elements. 

Chapter 4. Numerical Implementation in the FE Method.

o Algebraic equations are formed from the governing differential equations obtained 

from the constitutive specification. 

o Implementation of the model in Abaqus and Computer Aided Pavement Analysis 

(CAPA-3D) is explained briefly, and algorithms to implement the model’s system of 

equations in a general FE framework are presented.(2)

o The response of the constitutive and numerical model to the application of simple 

deformations is studied and presented.
3

Chapter 5. FE Simulation of SGC Compaction.

o The mechanics involved with Superpave® gyratory compaction are described briefly. 

o The boundary conditions necessary to simulate the SGC process are identified and 

described. 

o Simulations that were performed to study the model response with changes in the 

model parameters are presented. 

o The results of a sensitivity analysis are presented to gauge the relevance of the 

parameters in influencing the compaction response.

o SGC simulations are presented as validation for the FE model by comparing experimental

data at two different angles of compaction in the SGC with the model results.

Chapter 6. FE Simulation of Field Compaction.

o The FE model used for field compaction is described. 

o The effects of compaction equipment characteristics such as drum weight, amplitude 

of applied load, and frequency of applied vibratory loads are studied. 

o Parametric analysis is performed to gauge if model response is similar to what can be 

expected based on laboratory compaction simulation. 

o Sensitivity analysis is performed to reconcile significant parameters in field 

compaction simulations with those from SGC simulations. 

o The ability to simulate longitudinal joint modeling is demonstrated. 

Chapter 7. Simulation of Field Compaction.

o Field compaction cases are simulated using the developed FE model.

o Three highway projects are studied.

Chapter 8. Summary and Conclusions.

o Observations of the material response of the model are considered when subject to 

simple deformation mechanisms.

o Relevant conclusions are drawn correlating the laboratory and field compaction results. 

o Factors that impact laboratory and field compaction processes are presented.

o Suggestions are provided for further model improvement and future work on HMA

modeling. 

NOTATIONS

The notations used in this report are similar to those used in standard continuum-mechanics 

texts. Vectors and tensors (second-order and fourth-order tensors) are represented as follows:

a is a vector.

T is a second-order tensor.
4

The standard gradient and divergence operators are employed as presented in figure 1 and figure 2.

Grad a = ∂a

∂X , Grad [a]ij = ∂aaii

∂XXjj

,

Div TT = ∂TT

∂XX , Div [TT]i = ∂TTiijj∂XXjj

,

Figure 1. Equation. Gradient and divergence operators with respect 

to reference configuration

grad a = ∂a

∂x

, grad [a]ij = ∂aaii

∂xxjj

,

div TT = ∂TT

∂x , div [TT]i = ∂TTiijj∂xxjj

Figure 2. Equation. Gradient and divergence operators with respect 

to current configuration.
5

CHAPTER 2. REVIEW OF THE LITERATURE

Asphalt concrete is a composite material made up of graded aggregate rocks bound together by 

an asphalt binder along with the presence of air voids. The nature of this mixture is dependent on 

the type of each of the constituents selected to prepare the mix. Asphalt concrete pavement, more 

commonly known as HMA pavement, consists of bound layers of a flexible pavement structure. 

For most applications, asphalt concrete is placed as HMA, which is a mixture of coarse and fine 

aggregate and asphalt binder. HMA derives its nomenclature from the fact that it is mixed, placed,

and compacted at elevated temperatures. Asphalt concrete pavement can also be placed at ambient

air temperatures, but HMA is the primary placement method for roads and interstates.

Because of its composition, asphalt concrete displays a nonlinear response even at small strain and

exhibits a different response in tension and compression, particularly for long-duration loading. 

Also, the extreme temperature sensitivity of asphalt concrete causes a change in temperature to 

have a more significant impact on the mechanical behavior than changes in loading magnitude.

(3)

The deformation resistance of asphalt concrete is mainly derived from the aggregate matrix and 

the viscous asphalt binder matrix. Because of the change in microstructure, due to either mechanical

changes such as reduction of air voids or chemical changes such as aging of the asphalt, the response

of the aggregate matrix and the asphalt binder matrix to traffic loading changes over time. Also, 

the ability of the pavement to stress-relax upon load removal changes because of continual changes

in the microstructure. This change in microstructure and changes in loading and environmental 

conditions cause phenomena such as ruts, fatigue cracking, low-temperature cracking, and 

moisture-induced damage. Distress due to rutting is caused by the accumulation of deformation 

under repeated traffic loading, which results in the development of longitudinal ruts along the 

pavement. This accumulation of deformation depends largely on the one-dimensional densification

due to air-void reduction and the flow of the asphalt-mortar matrix. A thorough review of the nature

and uses of asphalt concrete and the many attempts at modeling the material are presented by

Krishnan and Rajagopal.

(4)

The current approach to flexible pavement design relates the engineering properties of asphalt 

mixtures to pavement distress. Therefore, current models are mostly distress-prediction models 

and not constitutive models as usually understood in mechanics. In distress-prediction models, 

some measure of distress—for example, the number of cycles to failure and amount of permanent

deformation—is related to material properties, assuming asphalt mixtures to be either linear

viscoelastic or linear elastic. However, HMA behaves like a nonlinear fluidlike material at elevated

temperatures and slowly transforms into a highly viscous nonlinear viscoelastic fluid toward the 

end of the construction process.

HMA COMPACTION

Proper construction of roadways requires that pavements be laid down according to certain 

specifications governing the desirable characteristics for the material. Construction of HMA

pavements takes advantage of the direct and indirect compacting forces exerted by rollers passing

over the loose mix to produce dense layers of structurally durable material. Compaction reduces 

the volume of a mixture of hot asphalt binder, aggregates, and filler materials to form the required

dense, impervious mass. The motivation for compacting an asphalt pavement is a desire to achieve
6

an optimum air-void content, to provide a smooth riding surface, and to increase the load-bearing 

capacity of the material under construction.

(1) The densification of the mixture due to compaction 

causes an increase in the unit weight of the material and improves the aggregate interlock.

(5,6)

The literature notes extensively that improper compaction generally leads to poor performance 

from the asphalt pavement, in spite of all other desirable mixture-design characteristics being met. 

This can lead to premature irreparable damage to the built-in infrastructure, typically in the form of 

rutting, permanent deformation, cracking, and moisture damage.

FACTORS AFFECTING COMPACTION

The compaction process is influenced by many factors, including the properties of the materials 

in the mixture, environmental variables, conditions at the work site, and the method of compaction,

as detailed by the U.S. Army Corps of Engineers and summarized in some detail by Kassem.

(5,7)

The required compaction effort increases with an increase in aggregate angularity, size, and 

hardness. The grade and amount of asphalt binder also influence the compaction process. A

mixture produced with too little asphalt is stiff and usually requires more compaction effort than 

a mixture with high asphalt-binder content. The temperatures of the air, mixture, and base are 

also important factors that influence compaction, according to the U.S. Army Corps of Engineers.

(5)

In addition, the compaction effort increases with an increase in layer thickness.

LABORATORY VERSUS FIELD COMPACTION

Several studies have attempted to examine the relationship between field compaction methods, 

laboratory compaction methods, and mechanical properties. Consuegra et al. and Harvey and 

Monismith evaluated several laboratory compaction methods based on the relationship between 

the mechanical properties of laboratory specimens and field cores.(8,9) These studies provided

recommendations concerning devices that produce laboratory-compacted specimens with properties

that better relate to those of field cores. Peterson et al. evaluated the influence of changing the 

compaction parameters in the SGC on the mixture’s mechanical properties and their correlation 

with the mechanical properties of field cores.(10) The control parameters that were varied in the 

study were the angle of gyration, specimen height, gyratory compaction pressure, and temperature

of the compaction mold. Peterson et al. found the angle of gyration to be the most important 

parameter influencing mechanical properties.(10)

Masad et al. used image-analysis techniques to study the air-void distribution in the SGC.(11)

The results showed that the air-void distribution in SGC specimens is not uniform and that 

the top and bottom layers have a higher air-void content than the middle layer. Tashman et al.

compared the air-void distribution and mechanical properties of SGC-compacted specimens 

to field cores compacted using different compaction patterns.(12) That study showed that the 

compaction parameters in the SGC can be changed in order to improve the relationship between 

the internal structure (including air-void distribution) and mechanical properties of SGC 

specimens and field cores.

More detailed descriptions of the various factors influencing compaction and compaction 

methods/techniques and a summary evaluation of various studies into the compaction of HMA 

are presented in Huerne and Kassem.(13,7)
7

MODELING HMA RESPONSE

The motivation for most studies of asphalt concrete behavior is to develop models to understand the 

distress an asphalt concrete pavement is subjected to and to determine experimental variables based 

on the analysis of such distress models. Unlike asphalt, for which models that describe the behavior

of viscoelastic fluids have been developed, asphalt concrete has constitutive specifications that are

mostly related to empirical correlations for different types of distress. However, studies assuming 

asphalt concrete to be a viscoelastic material also assume that the macroscopic mechanical behavior

of asphalt concrete is viscoelastic. These studies use either a spring-dashpot analogy in the form of

a Burger’s model or some simple form of viscoelastic constitutive equation. (See references 13–28.)

As noted by Krishnan and Rajagopal, the deformation resistance of asphalt concrete is mainly

derived from the aggregate matrix and the viscous asphalt mastic.(4) Because of the change in the 

microstructure (either due to mechanical changes such as reduction of air voids or chemical 

changes such as aging of asphalt), the response of the aggregate matrix and the asphalt mastic to 

traffic loading changes with time. Also, the ability of the pavement to stress-relax upon load

removal changes as the microstructure is continuously modified. This change in microstructure 

and changes in loading and environmental conditions cause phenomena such as rutting, fatigue 

cracking, low-temperature cracking, and moisture-induced damage. For instance, distress due to 

rutting is caused by the accumulation of deformation under repeated traffic loading, resulting in 

the development of longitudinal ruts along the pavement. While most studies have not considered

the mechanism of densification (assuming that the pavement will be compacted well during

construction), the plastic flow of the asphalt concrete has been assumed to be dependent on the 

temperature, loading rate, and loading time interval.(29,30)

There are quite a few models for asphalt concrete that consider the microstructure, but as with most

of the phenomenological models, they neglect the evolution of the microstructure during the lifetime 

of the pavement or consider it, for instance, by means of some shift factors. Nijboer conducted 

microstructural modeling of asphalt concrete using the analogy of soil mechanics, postulating 

that the entire deformation resistance of bituminous mixes can be explained in terms of initial

resistance, internal friction, and viscous resistance.(31) Huschek used a three-phase system

consisting of regions characterized by viscosity, modulus of elasticity, and modulus of plasticity.(23)

Van der Poel modeled the behavior of asphalt mixes by calculating the rigidity of concentrated 

solutions of elastic spheres in an elastic medium, using a method developed for dilute dispersions 

by Frolich and Sack.

(32,33) Hills developed models for the long-term creep behavior of asphalt 

mixes by characterizing the internal structure of the mix by means of the asphalt film thickness.(34)

Cheung et al. and Deshpande and Cebon have developed models for asphalt-concrete mixes using

isolated contact models and shear box models.(35–37) Boutin and Auriault used the analogy of a 

porous medium saturated by a viscoelastic fluid to classify the macroscopic behavior of asphalt 

concrete as biphasic, elastic, or viscoelastic depending upon the ratio of the dimension of the 

pores to the macroscopic wave length.(38) Florea used the viscoplastic potential to develop an 

elastic/viscoplastic model for bituminous concrete.(39,40)

Recently, Krishnan and Rao attempted to model the air-void reduction in asphalt concrete using 

the continuum theory of mixtures and the theory of linear elastic material with voids.(41,42) A 

thermodynamic framework has recently been put into place and can be used for the constitutive 

description of asphalt concrete. This framework has a reasonably general structure within which 
8

a host of dissipative processes can be described. This framework to model asphalt concrete 

recognizes the change in the microstructure of the material through the changes in the natural 

configurations of the body. For the purpose of this study, a natural configuration is a stress-free 

configuration, with further details of the framework presented by Rajagopal and later adopted by 

Rajagopal and coworkers. (See references 4 and 43–46.)

COMPACTION MODELS

Very little research has been directed toward modeling HMA compaction and the material properties

that influence compactability. Guler et al. proposed the use of a porous elastoplastic compaction 

model using a modified Gurson-Tvergaard yield function.

(47) An incremental constitutive relation 

for the porous material was formulated for this purpose. Researchers focused on obtaining 

statistically significant parameters for this constitutive relation and a correlation between the model

parameters and mixture variables (i.e., volumetric properties and particle size). Simple linear 

models were built to predict the model parameters. The displacement field used to represent 

three-dimensional compaction was an approximation of the actual motion in an SGC. Also, the 

model was formulated assuming the small-strain theory, was time independent, and assumed 

isothermal conditions (no changes in temperature).

Huerne used a modified form of soil critical-state theory in modeling asphalt-mixture compaction.(13)

The critical-state theory describes granular material behavior by means of a closed yield locus, which

gives a boundary between stress states that cause elastic (recoverable) deformations and plastic 

(irrecoverable) deformations. Huerne’s implementation simulated void reduction by means of 

plastic volume changes. The Hveem device was used for determining the model’s parameters. 

The model was developed assuming small strain deformation and, consequently, is limited in 

modeling the high strains involved in the compaction process. Also, the model has many 

parameters that are not directly linked to mixture properties.

Krishnan and Rao developed a constitutive model for asphalt mixes using mixture theory to model

the one-dimensional compaction of asphalt mixtures under a static load.(41) This model utilizes 

the fundamental balance laws to obtain mathematical relations to describe the performance and 

characteristics of asphalt mixes. While their work places the modeling within a general framework

that takes into account the balance laws of mechanics, it ignores certain critical issues concerning 

the material response, such as the evolution of the microstructure of the material being compacted

during the process. Also, such an approach to modeling compaction of HMA is limited by the

restrictive experimental techniques available to measure the various mixture properties involved 

in the model.

In summary, previous attempts to model HMA compaction were limited in the following ways:

They did not account for changes in the material structure during compaction either 

explicitly or implicitly through appropriate choices for the total dissipation.

They were not validated using various laboratory and field loading and boundary conditions. 

They did not attempt to develop methods for the determination of the model’s parameters 

based on mixture design and material properties. 
9

Recently, efforts have been made to develop material models to describe asphalt concrete behavior

based on a thermodynamically consistent methodology by employing the framework of multiple 

natural configurations formalized by Rajagopal.

(48,49,43)

INTELLIGENT COMPACTION

An important motivation for developing a material model to be used in simulations is the concept of 

IC. Moore defines IC as a system that applies to a vibratory roller and that automatically adjusts the

energy output of the roller to neither undercompact nor overcompact the materials.(50) IC of soil and 

asphalt layers is a relatively new technology introduced in the United States for compaction 

control.

(51) IC for HMA compaction is in an even earlier stage of development. This new

compaction concept is, however, gaining attention from the asphalt-paving industry. The application 

of the model developed here can be extended to obtain correlations between pavement designs and 

model characteristics. This will enable the decision process in an IC control system.

11

CHAPTER 3. A THERMOMECHANICAL MODEL FOR ASPHALT CONCRETE

The framework this study adopts to model asphalt concrete recognizes that, corresponding to any 

current configuration of a body, several associated natural configurations are possible. Of these 

natural configurations, the body may go to one specific natural configuration when the tractions 

are removed, preferably a stress-free state (see figure 3).

Figure 3. Illustration. Evolution of natural configurations associated with microstructural

transformations resulting from the material response to an external stimulus.

The specific natural configuration that the body goes to depends on the process. For example, 

perfectly elastic materials have a unique stress-free configuration, deferring only by a rigid body 

motion from each other. During deformation, the underlying microstructural mechanism does not 

change, and the body returns to the same stress-free state on the removal of applied traction. For a

classical viscoelastic material capable of instantaneous elastic response, the stress-free state obtained

by an instantaneous elastic unloading response can be associated with a natural configuration. To 

develop a model for HMA capable of representing the transformation of the mix (i.e., a low-viscosity

fluidlike state to a high-viscosity fluidlike state), researchers began with the same assumption by 

associating the natural configuration of the material with the stress-free state obtained upon

instantaneous elastic unloading. The model developed is an isothermal, compressible, viscoelastic,

nonlinear fluidlike model. The chief constitutive assumptions for the model are presented briefly, 

and further details regarding model development have been published elsewhere.(48,49)

THERMODYNAMICS OF DISSIPATIVE MATERIALS

The material model developed has to satisfy the second law of thermodynamics in the form of 

the energy-dissipation inequality shown in figure 4.

T ∙ L − ρρψψ̇ − ρρζζθθ̇ − q ∙ gradθθ

θθ = ρρ Ξ = ξξ ≥ 0

Figure 4. Equation. Energy-dissipation inequality.
12

Splitting the entropy production into two parts—that due to thermal effects and that due to 

mechanical dissipation—and requiring both of them to be individually nonnegative leads to the 

equation in figure 5.

T ∙ L − ρρψψ̇ = ξξ̂ ≥ 0

Figure 5. Equation. Reduced energy-dissipation inequality.

CONSTITUTIVE MODEL

In the current model, temperature is assumed to be constant during laboratory compaction. In 

general, the temperature effect can be accounted for by allowing the material parameters in the 

constitutive model to be different at different temperatures.

The specific Helmholtz potential for a solid, ψ, is given in figure 6.

ψψ = μμ(IIIG)

2ρρκκpp(tt)

�IBκκpp(tt) − 3 − ln �IIIBκκpp (tt)��

Figure 6. Equation. Helmholtz potential for mixture.

The rate of dissipation due to the mechanical working is specified in figure 7.

ξξ̂ = ηη(IIIG)Dκκpp(tt) ∙ Bκκpp(tt)Dκκpp(tt)

Figure 7. Equation. Rate of dissipation function.

Since the mixture is considered a single constituent, the shear-modulus function, μ, can be thought

of as a material property that reflects the characteristics of the aggregate matrix, and the viscosity 

function, η, reflects the characteristics of the asphalt mastic. Such a specification of the rate of 

dissipation corresponds to a material with a power-law-type viscosity and is therefore used to 

model the non-Newtonian fluidlike nature of the mix. During compaction, the aggregate matrix,

which is initially in a loose form, slowly evolves into a dense aggregate matrix. This densification

is further aided by the changes in the asphalt mastic due to densification and reduction of mastic 

film thickness. The specific forms for the shear-modulus and viscosity functions chosen (figure 8

and figure 9) represent such evolution.

μμ = μμ̂(1 + λ1(IIIG)2nn1 )qq1

Figure 8. Equation. Shear-modulus function.

ηη = ηη̂(1 + λ2(IIIG)2nn2 )qq2

Figure 9. Equation. Viscosity function.

To take into consideration the highly compressible behavior during the initial stages of compaction

and the material hardening as it reaches the limits of compressibility toward the end of the 

compaction process, choices are made for μ and η in figure 8 and figure 9.

(xˆ) ( )
13

CONSTITUTIVE EXPRESSION FOR STRESS

Considering the form chosen for the Helmholtz potential and the rate-of-dissipation function, the 

constitutive form for stress is obtained in figure 10.

TT = 2ρρIIIBκκpp(tt)

∂ψψ

∂IIIBκκpp(tt)

I + 2ρρ

∂∂ψψ

∂∂IBκκpp(tt)

Bκκpp(tt)

Figure 10. Equation. Constitutive equation for stress.

THE EVOLUTION EQUATION FOR BBκκpp(tt)

Thereafter, the evolution equation for Bκκpp(tt) is obtained in figure 11.

Bκκpp(tt) = − 2

ηη �Vκκpp(tt)

−1 TVκκpp(tt) − ρρIIIG

∂ψψ

∂IIIG

I�

Figure 11. Equation. Evolution equation.

ANALYTICAL CALCULATIONS USING A SEMI-INVERSE APPROACH

For validation of the model, consideration needs to be given to three-dimensional compressive 

deformation with applied static and dynamic loads. However, to do this, the model needs to be 

implemented in an FE methodology that will enable modeling the three-dimensional behavior in 

general. Therefore, researchers chose to study the developed model through analytical calculations

aided in part by the computations. Presented in this section are the formulation for such a study and 

the calculations performed using numerical computation software MATLAB®.

(52)

The governing equations of motion for the constrained one-dimensional compression of the material

model were developed using a semi-inverse approach by prescribing the necessary deformation 

field for this case. The model’s response to constant stress and constant strain is illustrated in the 

following sections.

ONE-DIMENSIONAL CREEP DEFORMATION

The diagonal components of stress are unknown except for the prescribed creep load (see figure

12). Hence, a constant compressive stress is applied along the axial z-direction.

Considering that the specific application is a one-dimensional constraint creep, the axial stress

and other similar dimensional quantities can be normalized using the initial applied creep load. 

∇
14

Figure 12. Illustration. Schematic for the one-dimensional compression problem.

Typical plots obtained for the given material by calculating the solution to the semi-inverse creep 

problem using MATLAB® are depicted in figure 13 through figure 18. Kinematical terms are 

nondimensional quantities. The material model is subject to a sudden normal stress, and the 

deformation response shows the viscoelastic nature of the material. 

Figure 13. Chart. Creep solution for the material model for HMA—creep load.
15

Figure 14. Chart. Creep solution for the material model for HMA—Almansi strain.

Figure 15. Chart. Total stretch calculated from the solution to the one-dimensional 

creep problem.
16

Figure 16. Chart. Evolution of natural configuration calculated from the solution to the 

one-dimensional creep problem.

Figure 17. Chart. Stored energy calculated from the solution to the one-dimensional 

creep problem.
17

Figure 18. Chart. Rate of dissipation calculated from the solution to the one-dimensional 

creep problem.

Strain lags behind stress because of the material’s viscoelasticity. Upon instantaneous removal of 

the load, the material rebounds partially to attain a permanent strain. Plots of some of the relevant

field variables calculated after obtaining the solution are presented in figure 15 through figure 18.

All quantities have been nondimensionalized. The nondimensionalized stored energy has negative

values due to the nondimensionalization of shear modulus with respect to the initial applied load. 

The stored energy starts high and then, as the material evolves, dissipates (though the rate of 

dissipation decreases) to none upon unloading. The natural configuration of the material evolves 

as well during the process, and its evolution is presented in figure 16. The plot shows that the 

material returns to a natural (stress-free) configuration upon unloading; therefore, the dissipation 

occurs because of a change in the microstructure of the material. 

STRESS RELAXATION

Similarly, researchers solved for the motion of the one-dimensional deformation under a constant 

applied strain (i.e., with respect to the current configuration) and observed the corresponding 

stress response as a function of time (see figure 19). 
18

Figure 19. Illustration. Schematic for the one-dimensional stress-relaxation problem.

The calculations of the response of the model to a step-input strain results in the characteristics 

shown in figure 20 through figure 23. All quantities have been nondimensionalized.

Figure 20. Chart. Applied strain calculated from the solution to the one-dimensional 

stress-relaxation problem.
19

Figure 21. Chart. Stress relaxation calculated from the solution to the one-dimensional 

stress-relaxation problem.

Figure 22. Chart. Stored energy calculated from the solution to the one-dimensional 

stress-relaxation problem.
20

Figure 23. Chart. Rate of dissipation calculated from the solution to the one-dimensional 

stress-relaxation problem.

The applied compressive strain is 0.005, and the resultant nondimensionalized stress response 

depicted in figure 20 through figure 23 shows that the stress relaxes slowly and reaches a

stress-free state. The stress-free state is reached after the material dissipates the energy it stores to 

change its microstructure. Again, the nondimensionalized stored energy has negative values due 

to nondimensionalization.

Therefore, the material model predicts the expected response to mechanical loads for compressible,

viscoelastic material models.
21

CHAPTER 4. NUMERICAL IMPLEMENTATION IN THE FE METHOD

The modeling framework and the proposed model used for the compaction study were introduced

in chapter 3. The model was further developed and improved to account for the behavior of asphalt

concrete subjected to compaction in the field. To this end, the model was implemented in the FE

software CAPA-3D, developed at Delft University of Technology.

(53,54) CAPA-3D was developed 

as an FE-based platform to serve the computational needs of the research group Mechanics of 

Structural Systems at Delft University of Technology and international teams that cooperate with 

the group. Over the years, the software has evolved into a full-fledged FE system for static or 

dynamic analysis of large-scale three-dimensional pavement and soil engineering models. The 

software consists of a sophisticated user interface; a powerful, highly parallelized band-optimizing 

mesh generator; a high-quality, user-controlled graphical output; several material and element 

types; and a variety of specialized algorithms for more efficient analysis of pavement construction. 

These algorithms include a moving load-simulation algorithm and a contact algorithm, which are 

essential for modeling the compaction process. 

VERIFICATION OF FE IMPLEMENTATION

Researchers solved the constrained compression problems using the FE implementation of the 

model in CAPA-3D. In finite elements, these problems were solved by employing a 20-noded 

three-dimensional solid element. The representative material element was constrained to move 

only in a normal direction (z-direction). A comparison was drawn between the solutions obtained 

from the analysis (with aid from calculations in MATLAB®) and the solutions obtained using the 

method of finite elements.(54) The compressive strain applied was 0.05, or 5 percent. The material 

properties used are presented in table 1. An explanation of the parameters listed in table 1 is 

presented in chapter 5 of this report. These parameters were also used for FE calculations.

Table 1. Parameters employed for the study of the material response when subject to 

one-dimensional constrained compression. 

μμ� (MPa) n1 λ 1 q1 ηη� (MPa·s) n2 λ 2 q2

810 4.0 0.25 −15 1,400 2.5 0.25 −25

1 MPa = 0.145 ksi

Inputs to the model are given by linearly increasing the input to a constant value and then keeping

the input at this value. Removing the input is then performed by linearly decreasing from the 

maximum value to no input. The comparison for constant applied stress is presented in figure 24

and figure 25, and the comparison for constant applied strain is presented in figure 26 and figure 27.

The calculated values (using MATLAB®) of the stretch and stress in response to applied stress 

and strain, respectively, and the corresponding FE solutions agree well, as can be observed from 

figure 24 through figure 27. This agreement serves as validation for the implementation of the 

model in CAPA-3D using the method of finite elements.
22

Figure 24. Chart. FE solution to the model response to applied constant stress—creep load.

Figure 25. Chart. FE solution to the model response to applied constant stress—total stretch.

Figure 26. Chart. FE solution to the model response to applied constant strain 

(compressive 0.05)—applied strain.
23

Figure 27. Chart. FE solution to the model response to applied constant strain 

(compressive 0.05)—normal stress response.

CONSTANT SHEAR

In the developed FE implementation, a constant-shear problem is solved. The schematic 

representation for such a model is given in figure 28. The unit cube is subject to a rate-type shear

deformation, and the boundary conditions to be applied to the unit cube are such that the cube is

constrained to move only in a lateral direction, indicating a shear deformation. 

Figure 28. Illustration. Schematic of the constant shear loading applied to a unit cube.

In figure 29 through figure 31, T12 represents the shear in the xy-plane, T11–T22 represents the 

first normal stress, and T22–T33 represents the second normal stress. As shown in the figures, the

material response corresponds with that of a nonlinear material due to the exhibition of normal 

stress differences.
24

Figure 29. Chart. Shear stress (T12) observed in response to constant shear loading.

Figure 30. Chart. Comparison of the first normal stress (T11–T22) response 

to constant shear loading.
25

Figure 31. Chart. Comparison of the second normal stress (T22–T33) response 

to constant shear loading.

CONSTANT SHEAR RATE

Using the developed FE implementation, a constant-shear-rate problem is solved. The schematic 

representation for such a model is given in figure 32. The unit cube is subject to a rate-type shear 

deformation, and the boundary conditions applied to the unit cube are such that the cube is 

constrained to move only in a lateral direction, indicating a shear deformation.

Figure 32. Illustration. Schematic of the constant shear rate applied to a unit cube.

In figure 33 through figure 36, T12 represents the shear in the xy-plane, T11–T22 represents the 

first normal stress, and T22–T33 represents the second normal stress.
26

Figure 33. Chart. Shear stress (T12) observed in response to constant shear rate loading.

Figure 34. Chart. Shear stress (T12) as a function of the shear rate 

(using the model parameters in table 1).
27

Figure 35. Chart. Comparison of the first normal stress (T11–T22) response 

to constant shear rate loading.

Figure 36. Chart. Comparison of the second normal stress (T22–T33) response

to constant shear rate loading.
28

As shown in figure 33, figure 35, and figure 36, the material response corresponds with that of a 

nonlinear material due to the exhibition of normal stress differences (though small in comparison 

to the shear stress T12). For the chosen set of parameters, the model seems to exhibit near-linear 

shear-stress (T12) versus shear-rate characteristics. Such characteristics typically represent the 

response of a fluid with a rate of dissipation similar to Newtonian fluid. However, the model does 

exhibit normal stress differences, a typical characteristic of non-Newtonian fluids, as shown in 

figure 35 and figure 36. Hence, the model behaves like a non-Newtonian fluid, and the choice of 

the parameters enhances the shear-thickening nature of the model (unlike in figure 34). This feature

of the model is very useful for the purposes of modeling compaction in the laboratory or in the 

field as the material in those situations undergoes shear-thickening behavior. The choice of model 

functions for shear modulus and viscosity ensure such a transition in the material’s property.
29

CHAPTER 5. FE SIMULATION OF SGC COMPACTION

RELATIONSHIPS OF MODEL-INPUT PARAMETERS TO MATERIAL PROPERTIES

In the laboratory, the compaction process takes place at high temperatures (see figure 37 for an

example of the equipment used). Some phenomena and parameters necessary for the model only 

occur or are only relevant at higher temperatures. This makes it impossible to conduct conventional

laboratory testing of specimens (e.g., triaxial testing) to obtain these parameters. This limitation was

overcome by implementing the mathematical model in finite elements, conducting simulations of

the SGC compaction process, and then fitting the FE results to the available data in the form of

SGC curves to determine the model parameters. These parameters were then considered to include 

the information necessary to characterize material behavior under various conditions. Assuming that

the material in the laboratory is representative of that in the field, the model parameters obtained

from the SGC calibration can be used to simulate field compaction under loading and boundary 

conditions in the field (see figure 38 for an example of field rolling compaction). This does not 

imply that the SGC process can be assumed to be the same as field compaction.

Figure 37. Photo. Superpave

® gyratory compactor.
30

Figure 38. Photo. Static steel-wheel roller.

Once the model input parameters were identified, correlations to material properties could be 

found. It was expected that the model input parameters would be influenced by the following 

material properties:

Aggregate gradation. 

Aggregate shape characteristics. 

Binder viscosity at compaction temperatures.

Mixture temperature. 

Model Verification and Calibration

Several studies have been conducted into the experimental characterization of the influence of 

different laboratory and field compaction procedures on the internal structure of asphalt mixtures 

and their mechanical properties. From these studies, the following information and materials can 

be used in verification and calibration:

Pavement structure and layers. 

Equipment used in the project (screed, dump truck, roller, etc.) and its size, weight, 

pattern, sequence, number of passes, etc.

Number of passes of each compactor across the pavement width. 

Compaction data from the SGC. 

Mixture design and volumetrics. 
31

Available data were used to determine the relationship between mixture properties, temperature, 

and model parameters. The model parameters were determined by calibrating the FE model to 

match the SGC curves. This was followed by verification of the FE implementation’s prediction 

using a set of data that was not used in the model calibration. Field compaction measurements were

used for further calibration and validation. FE models that simulate the loading and boundary 

conditions in the field were developed. The model parameters that were determined using the 

laboratory calibration were used initially in the field FE models. Calibration was introduced as 

needed to better simulate the field compaction. Such calibration was necessary to account for 

possible differences in materials between the laboratory and the field.

The FE implementation of the model developed was used to simulate the compaction processes 

and the simulated compactor response, and changes in mixture stiffness at different stages of 

compaction were compared to the measured response and mixture stiffness.

The constitutive model developed was incorporated into two FE programs, Abaqus (a commercial

FE package) and CAPA-3D. The constitutive model was implemented in Abaqus through the user

material subroutine, UMAT.

(2) The CAPA-3D general structure has been presented in detail by 

Scarpas.

(53,54) For the implementation of the present model, CAPA-3D required the input of two 

interface subroutines that provide the system with the necessary model parameters and the algorithm

for the required stress updates at each material point. The required equations were obtained 

through a numerical discretization of the material constitutive relations of the model to obtain a 

system of nonlinear algebraic equations at each material-integration point. The computational 

sequences were then provided in the interface subroutines in the form of Fortran code. Further 

details on the usage of the computational system are provided in the CAPA-3D user manuals.(53,54)

The SGC was simulated through a lateral constraint along the boundary and a constant vertical 

stress of 87 psi (600 kPa). In addition, sinusoidal vertical displacements that were 120 degrees

out of phase were applied at equidistant locations around the periphery of a plate supporting the 

mix from the bottom. This was done to simulate the applied angle of gyration in the SGC,

which is responsible for inducing shear in the mix. The frequency of the applied gyrations was

30 gyrations per minute. The bottom support plate was modeled as a nearly rigid material. Also, 

a constant friction factor of 0.4 was introduced between the mixtures and the plates and mold.

The compaction curves obtained from the simulations were plotted as the normalized height of 

the specimen versus time of compaction. The normalized height is defined as the height of the 

specimen at any time divided by the initial height of the mixture in the SGC. Figure 39 shows an 

example of the FE mesh of the SGC.
32

Figure 39. Illustration. FE mesh used in modeling the SGC.

PARAMETRIC ANALYSIS OF MODEL PARAMETERS

A comprehensive parametric analysis was conducted in order to determine the sensitivity of the 

mixture compaction to the various models’ parameters. The parameters that were used in this 

analysis are shown in table 2. All the SGC simulations were conducted using an angle of

1.25 degrees.
33

Table 2. Model parameters used in the parametric study.

Set No. μμ� (MPa) λ1 n1 q1 ηη� (MPa·s) λ2 n2 q2

1 1,700 0.22 4 −24 1,400 0.25 3 −30

2 1,900 0.22 4 −24 1,400 0.25 3 −30

3 1,500 0.22 4 −24 1,400 0.25 3 −30

4 1,700 0.25 4 −24 1,400 0.25 3 −30

5 1,700 0.20 4 −24 1,400 0.25 3 −30

6 1,700 0.22 5 −24 1,400 0.25 3 −30

7 1,700 0.22 3 −24 1,400 0.25 3 −30

8 1,700 0.22 4 −26 1,400 0.25 3 −30

9 1,700 0.22 4 −21 1,400 0.25 3 −30

10 1,700 0.22 4 −24 1,600 0.25 3 −30

11 1,700 0.22 4 −24 1,300 0.25 3 −30

12 1,700 0.22 4 −24 1,400 0.26 3 −30

13 1,700 0.22 4 −24 1,400 0.22 3 −30

14 1,700 0.22 4 −24 1,400 0.25 4 −30

15 1,700 0.22 4 −24 1,400 0.25 2.5 −30

16 1,700 0.22 4 −24 1,400 0.25 3 −26

17 1,700 0.22 4 −24 1,400 0.25 3 −34

1 MPa = 0.145 ksi

An example of the influence of μμ̂, keeping all other parameters the same, is shown in figure 40. 

As expected, more compaction was achieved as μμ̂decreased. The effect of this parameter on the 

compaction process becomes more important at longer compaction times, when the mixture

becomes more or less a highly viscous fluid. The analysis of parameter n1 (shown in figure 41) 

showed that μ is significantly affected by the n1 value (see figure 8), which controls the maximum

compaction (change in height) of a mixture.

Figure 40. Chart. Analysis of the sensitivity of compaction to μμ�.
34

Figure 41. Chart. Analysis of the sensitivity of compaction to n1.

Analysis of the λ 1 parameter revealed that mixture compaction is affected very little, if any, by 

this parameter within the range of 0.2 to 0.3. This specific range for λ 1 is significant to obtain the 

initial viscous dissipation during the deformation. Therefore, λ 1 was treated as a constant chosen 

to lie within this range in the model in order to simplify obtaining the remaining parameters.

Another parameter that has a very small effect on the compaction process is q1. Therefore, this 

parameter was also treated as a constant in the model.

Figure 42 shows that the parameter ηη̂controls the point at which the material behavior starts to 

change from a very low-viscosity fluidlike behavior (rapid compaction) to a highly viscous fluid 

behavior (slow compaction). In figure 42, the material with ηη̂= 232.06 ksi·s (1,600 MPa·s)

experiences change in the compaction rate prior to the other two materials. The parameter λ 2, 

presented in figure 43, is directly related to the initial slope of the compaction curve (initial viscosity

of the mixture). As shown in figure 44, the model parameter q2 contributes to the nonlinear change

of viscosity from the start of the compaction process, when the mixture exhibits low-viscosity

fluidlike behavior. The overall compaction of the mixture is higher at lower (more negative) values

of q2. The parametric analysis showed that n1 has almost no effect on the compaction curve when 

it lies between 2 and 4. Consequently, n1 is assumed to be a constant within this range for the SGC

simulations presented in this chapter.
35

Figure 42. Chart. Analysis of the sensitivity of compaction to ηη�.

Figure 43. Chart. Analysis of the sensitivity of compaction to λ 2.
36

Figure 44. Chart. Analysis of the sensitivity of compaction to q2 .

The equations for shear-modulus-like function μ and viscosity-like function η are represented in 

figure 45 and figure 46, respectively, after assigning constant values to λ 1, q1, and n2.

μμ = μμ̂(1 + 0.25(IIIIIIGG)2nn1 )−25

Figure 45. Equation. Modified shear-modulus function.

ηη = ηη̂(1 + λ2(IIIIIIGG)5)qq2

Figure 46. Equation. Modified viscosity function.

The viscosity function primarily controls the initial compaction when the material behaves like a 

fluid with very low viscosity. The change in viscosity during compaction is associated with changes

in the reduction in asphalt film thickness, leading to an increase in the overall viscosity of the

mixture. On the other hand, the shear-modulus function primarily controls the compaction process

after some compaction when the material starts to behave like a highly viscous fluid. The increase

in the shear modulus is associated with changes in the aggregate structure due to the development

of more contacts and interlocking.

Based on the parametric analysis and the physical significance of the model’s parameters, the 

parameters are summarized as follows:

μμ̂is a shear-modulus parameter and marks the point at which the mixture starts to behave

as a highly viscous fluid and the shear modulus dominates the compaction process.

n1 is a shear-modulus parameter and affects the rate of compaction after the mixture starts 

to behave as a high-viscosity fluid.
37

ηη̂is a viscosity parameter and controls the point at which the material behavior starts to 

change from a low-viscosity fluid behavior (rapid compaction) to a highly viscous fluid 

behavior (slow compaction).

λ 2 is a viscosity parameter and controls the initial slope of the compaction curve (initial 

viscosity of the mixture).

q2 is a viscosity parameter and contributes to the nonlinear change of viscosity. 

The relationship between the model’s parameters and the compaction process is illustrated in 

figure 47.

Figure 47. Chart. Illustration of the relationship of the model’s parameters 

to the compaction process.

SENSITIVITY OF COMPACTION TO ANGLE OF GYRATION

Figure 48 presents the results of FE simulations of the SGC at different angles of gyration. The 

model captures the influence of the angle on the change in compaction. Also, the compaction 

behavior does not change in a linear manner with the change in angle.
38

Figure 48. Chart. Influence of angle of gyration on the compaction curve.

The FE model was used to determine the maximum shear stresses at the top of the specimens for 

angles of gyration of 1.25 and 2.0 degrees (shown in figure 49 and figure 50, respectively). The 

analysis was conducted using the model’s parameters determined for the SH-36 mixture presented

later in this report. Initially, the shear stress decreases rapidly when the material behaves as a

compressible fluid and then starts to increase gradually as the mixture starts to behave like a highly

viscous fluid. The shear stresses are higher for higher angles of compaction because an increase in

the angle of gyration is associated with an increase in the applied shear stresses. Also, the point at

which shear stress starts to increase occurs earlier at a 2.0-degree angle than at a 1.25-degree angle.

This is because the mixture compacts and gains strength faster at a 2.0-degree angle of gyration.
39

Figure 49. Chart. Maximum shear stress at the top of the specimen for a gyration angle 

of 1.25 degrees.

Figure 50. Chart. Maximum shear stress at the top of the specimen for a gyration angle 

of 2.0 degrees.

MODEL CALIBRATION AND VERIFICATION

Mixtures from various projects were compacted using two angles of gyration, 1.25 and 2.0 degrees,

to correspond to field compaction. The compaction curves at 1.25 degrees were used to determine
40

the model’s parameters, while the curves at an angle of 2.0 degrees were used to verify the ability

of the model to fit experimental measurements different from those used to determine the parameters.

Each of the compaction curves was divided into three parts: (1) time of compaction less than 15 s,

(2) time of compaction between 15 and 30 s, and (3) time of compaction more than 30 s. The first

part was used to estimate the parameters λ 2 and ηη̂of the viscosity function. The second part was 

used to estimate q2, and the third part was used to determine the parameters μμ̂and n1 of the modulus

function. The initial estimates were then slightly adjusted to minimize the difference between the 

model’s results and the data of the whole compaction curve. The final simulation results using the

parameters in table 3 are shown in figure 51 through figure 55 for an angle of 1.25 degrees and in 

figure 56 through figure 60 for an angle of 2.0 degrees. The results show that the model has a 

reasonable representation of the compaction curves at both angles. The nonsmooth response 

characteristics of the simulation, as can be observed from the plots, are a natural consequence of 

the boundary conditions and the fact that the calculations are made at an integration point in the 

material mesh that is at a slightly eccentric position from the axis of symmetry of the model.

Table 3. Model parameters obtained from compaction data.

Mix Projects μμ� (MPa) n1 ηη� (MPa·s) λ 2 q2

IH-35 2,000 5 1,400 0.25 −27 

SH-36 2,400 4 1,800 0.22 −28 

US-87 2,600 5 2,100 0.26 −27 

US-259 2,400 5 1,700 0.21 −28 

SH-21 2,500 4 2,000 0.22 −30 

1 MPa = 0.145 ksi

Figure 51. Chart. Fitting of the compaction data at 1.25 degrees for project IH-35.
41

Figure 52. Chart. Fitting of the compaction data at 1.25 degrees for project US-259.

Figure 53. Chart. Fitting of the compaction data at 1.25 degrees for project SH-36.
42

Figure 54. Chart. Fitting of the compaction data at 1.25 degrees for project SH-21.

Figure 55. Chart. Fitting of the compaction data at 1.25 degrees for project US-87.
43

Figure 56. Chart. Prediction of the compaction data at 2.0 degrees for project IH-35.

Figure 57. Chart. Prediction of the compaction data at 2.0 degrees for project US-259.
44

Figure 58. Chart. Prediction of the compaction data at 2.0 degrees for project SH-36.

Figure 59. Chart. Prediction of the compaction data at 2.0 degrees for project SH-21.
45

Figure 60. Chart. Prediction of the compaction data at 2.0 degrees for project US-87.

The comparison between the parameters in table 3 and the mixture characteristics did not reveal clear

and simple relationships. This could be attributed to the fact that each of the parameters, whether 

part of the modulus or viscosity function, is related to combinations of the mixture characteristics.

The determination of these relationships requires the analysis of parameters of more mixtures in 

order to develop regression equations between the model’s parameters and mixture characteristics.

Nevertheless, there are some trends supporting the idea that relationships exist between material 

properties and model parameters. For example, the lowest values of ηη̂and μμ̂are for the stone matrix 

asphalt mixture in IH-35. These low values suggest better compactability, which could be attributed

to the high asphalt content of this mixture. Also, the model results for very low angles of compaction

do not agree well with what can be observed in the laboratory at such angles. This is because of the

high compressibility built into the model to simulate the large deformations generally observed 

during compaction, which is not exhibited by actual mixes at very low angles of gyration.

The parametric analysis conducted in this study showed that the compaction process is insensitive

to changes in some of the model’s parameters. Therefore, constant values were used for these

parameters during all the FE simulations in order to simplify the process of finding the values of the

remaining parameters. The parametric analysis demonstrated that some of the model’s parameters

are related mostly to the initial stage of compaction when the material exhibits low viscosity

fluidlike behavior, while some parameters are related to the behavior of the mixture after it starts 

to exhibit a highly viscous fluid behavior and compaction rate decreases. A systematic method was

developed to determine the model’s parameters from the SGC curves. It was necessary to develop

this method due to the difficulty of conducting conventional tests (e.g., triaxial or shear tests) to 

determine the model’s parameters given the high compaction temperatures. The FE simulations 

demonstrated the capability of the constitutive model to simulate the SGC process at different 

compaction angles.

47

CHAPTER 6. FE SIMULATION OF FIELD COMPACTION

INTRODUCTION

Development of a useful simulation environment to study field compaction of HMA provides a

unique challenge because of the complex physical interactions involved. A framework composed of

a nonlinear material model coupled with an advanced numerical simulation environment (CAPA-3D)

was chosen for this study. The theoretical model employed for the FE model development was

adapted from the HMA material model presented in Koneru et al.(48)

This model was subsequently 

utilized for the simulation of the compaction mechanism of the SGC. The results of the study

have been documented in Masad et al.

(49)

Here, the chosen numerical framework is used to extend

the work presented in Koneru et al. and Masad et al. to specifically study the material model’s 

behavior when subject to the conditions pertinent to field compaction of HMA.(48,49)

The simulation of the rolling compaction process includes many challenges in terms of representing

the complex physical interactions in a numerical model. Some of the chief factors that influence 

the research related to the field compaction simulation are as follows:

Development of a contact-area-determining algorithm was needed to determine the effective

area of a given FE mesh in contact with a roller used for applying a compacting load on 

the material. Since the material model chosen from Koneru et al. was a compressible 

viscoelastic material, devising an algorithm for contact had to be considered carefully.

(48)

Because of the large variation in the material properties of the different layers of a pavement

structure, it was anticipated that numerically simulating the interactions of the layers would

require specialized techniques. For this purpose, researchers used specialized finite elements

for resolving the interfacial contact issues.

Because of the application of a load that is translated over the surface in field compaction, 

special care was taken to avoid wave reflections from boundaries and interfaces. Such 

interference waves produced within the material can affect the solution in a numerical 

experiment. A standard technique to deal with such detrimental effects is the use of 

“impedance layers” in FE meshes. (See Sluys for further details on the utility and 

development of such elements.)

(55)

The effect of the stiffness of substructure layers needed to be studied because the structural

response of any material is, in general, subject to the material it interacts with in response 

to external influences.

A study of the sensitivity of compaction to model parameters provided an opportunity to 

understand the behavior of the model in a structural setting. The change in behavior of 

the FE model with changes in parameter values aids in the estimation of the model 

behavior based on given input conditions to the system. Because of the physical interactions

involved and a lack of detailed experimental data, the researchers needed such a study to 

draw correlations between laboratory and field compaction activities. 
48

FE MODEL

A typical pavement structure used for simulating field compaction is shown in figure 61. It consists

of the following three layers: a base layer consisting of aggregate and/or soil approximately

18 inches (457.2 mm) thick, an overlay on top of the base layer consisting of old HMA approximately

2 inches (50.8 mm) thick, and a 3-inch (76.2-mm)-thick top layer of freshly laid HMA that will be

compacted by rolling.

Figure 61. Illustration. Pavement structure typically employed for studying

field compaction.

The general structure of CAPA-3D facilitates quick, flexible testing of material models and has 

been presented in detail by Scarpas.

(53,54) For the implementation of the chosen material model,

two interface subroutines are input into CAPA-3D. The necessary model parameters, boundary 

conditions, and algorithmic logic necessary for the required stress updates at each material point 

are provided to CAPA-3D through user-interface files for each type of necessary data. The required

stress update subroutine was devised from the equations obtained through a numerical discretization

of the constitutive relations of the material model chosen. The stress-update subroutine consists of

a system of nonlinear algebraic equations that are iteratively integrated at each material-integration

point and each step in the time domain. The computational sequences devised are then provided 

as Fortran instructions that can be embedded directly into the CAPA-3D interface. Details on the 

usage of the computational system are provided by Scarpas.

(53)

FE Mesh Used for Simulation

Some of the important characteristics of the FE mesh used for simulating compaction (shown in 

figure 62) are as follows:

The representative mesh has material in three different layers along with the material already 

compacted by the roller. The mesh is a half-sectional view. The symmetry of boundary 

conditions sometimes allows the half mesh to be employed instead of the full mesh.

The typical geometry of the mesh is given the typical thicknesses shown in figure 61, a 

lane width of 12 ft (3.66 m), and a total top-layer length in the direction of rolling of 

44 inches (1,117.6 mm). The top layer consists of a newly laid mix and a previously

compacted mix, each 22 inches (558.8 mm) in length in the direction of rolling.
49

The typical material properties used for the three layers are given in table 4.

All the finite elements used for the structural layers are 20-noded, solid, three-dimensional

continuum elements. Details of other choices and options for the elements available for 

use in CAPA-3D are available from Scarpas.

(53)

The interfacial contact between any pair of layers among the uncompacted HMA, the old 

HMA layer underneath, and the previously compacted layer on top is resolved through use

of a special interface element available in CAPA-3D.(53) Use of this layer is strictly to 

preserve the contact between the layers, and hence, these special elements do not have any

structural utility.

Impedance layers surrounding the base layer and other such elements surrounding the top 

HMA and old HMA layers (as shown in figure 62) minimize the dynamic loading effects.

The analyses are performed by fixing the time increments at 0.05 s. This value is based 

on the actual velocity of the roller and takes into account a need for consistent numerical 

stability of the algorithm. The time increments are calculated by first determining the 

total time per pass by dividing the total distance of travel (22 inches (558.8 mm)) by the 

speed of travel of the roller (initially set to 4.4 inches/s (111.76 mm/s)) and then

multiplying by the total number of time steps desired per pass (assumed 100 here). 

Figure 62. Illustration. Sectional view of the FE mesh used for setting up 

field compaction simulations.
50

Table 4. Typical material properties of pavement structural layers.

Layer Model Thickness Properties

New compacted 

asphalt layer

Compressible 

nonlinear viscoelastic

3 inches initially, changes 

during compaction

Presented as parameter 

sets later in this chapter

Old asphalt layer Elastic 2 inches E = 300 ksi, n = 0.35

Base layer Elastic 18 inches E = 145 ksi, n = 0.20

1 inch = 25.4 mm

1 ksi = 6.89 MPa

Boundary Conditions for Simulating Field Compaction

The boundary conditions most frequently utilized are presented in figure 63. The schematic

illustrates that the fixed edge corresponds with the inner edge of the pavement lane while the free 

edge corresponds with the outer edge of the lane. The displacement on the fixed edge is thereby 

constrained in the lateral direction (z-coordinate direction for the given mesh orientation), the free 

edge is allowed to move unconstrained, and a traction-free boundary condition is employed.

When the middle lane of a highway is to be simulated, simplifying further requires assuming 

symmetry in boundary conditions and employing a half mesh along the direction of rolling. In 

this case, the half-sectional mesh represented in figure 62 is sufficient for the simulation. 

Figure 63. Illustration. Schematic diagram illustrating the edges of the lane that 

correspond to fixed and free edges of the mesh in figure 62.
51

Response to Moving Load

This section discusses the typical response of the material to a moving load in field compaction. 

The simulations were performed with material parameters that correspond to those in table 5. 

The same set of parameters in SGC and field compactions do not provide the same level of

volumetric reduction because of significantly different boundary conditions such as confinement, 

temperature changes, uniformity and nonuniformity of applied loads, etc.

Table 5. Parameters used to study the typical material response during field compaction.

Highway Project μμ� (MPa) n1 λ 1 q1 ηη� (MPa·s) n2 λ 2 q2

SH-21 620 4 0.25 −15 1,400 2.5 0.25 −20 

1 MPa = 0.145 ksi

The typical response curve obtained from the simulations employing the previously described 

mesh subject to the specified boundary conditions is shown in figure 64. The response curve 

represents the compaction experienced at a point in the roller path as the roller passes over it and 

returns for a second pass from the opposite direction. The annotation arrows signify the forward 

motion in the first pass and the reverse-direction motion in the second pass. Finally, the load is 

removed and the material allowed to relax as much as possible. Thereafter, the material finally 

exhibits its state of permanent deformation, as seen in figure 64.

Figure 64. Chart. Typical displacement curve for a node under the cylindrical load for a 

cycle with forward and return passes followed by a forward pass with the load removed.

Figure 65 and figure 66 illustrate the relationship between the roller position with a contact area of

3.9 inches (100 mm) and the deflection experienced by a point P near the surface of the pavement.

The illustrations are for a mixture that does not have instantaneous elastic response. Due to the 

viscoelastic nature of the material, the deflection at point P starts when the center of the roller 

has passed the point. The point continues to deform vertically as the roller moves forward.

Time Steps

Vertical Displacement, mm
52

When the roller is at point 4, the mixture at P rebounds some of its vertical deflection because it 

is out of the region of influence of the roller. Also, point 4 marks the time at which the roller starts

to move back toward P. Further deflection occurs as the roller again approaches P (point 5). The 

final deflection is marked by point 6.

Figure 65. Chart. Vertical displacement observed at point P on the path of the roller.

Figure 66. Chart. Roller location from observation point P. 
53

Each pass over the material results in permanent deformation upon unloading, as shown in figure 67.

The deformation therefore builds up incrementally with each pass. However, as shown in figure 68,

the increments decrease in magnitude due to material densification as the number of passes increases.

Figure 67. Chart. Deflection at the node of interest when a load is applied for a short 

duration and then removed.

Figure 68. Chart. Permanent deformation prediction of the model in multiple-pass loading. 
54

Implementing a Contact-Area Algorithm

The objective was to implement an algorithm for load application on flexible pavement. Development

of the new algorithm furthers the ability to simulate the compaction loading. The loading algorithm

used in this project was based on translating a distributed load across the pavement. 

The load was applied in the form of a grid of loading patches independent of the FE mesh. The 

pressure in each loading patch was uniform. By combining a suitable number of loading patches, 

any desired distribution of surface loads could be simulated. The loading patches were translated 

along the surface of the pavement to simulate the moving roller.

A process was developed for determining the pressure of each of the loading patches. The process

involves the determination of the contact area between the roller and the pavement (see figure 69).

A displacement that mimics the cylindrical profile of the roller is applied on the surface of the 

pavement. The number of nodes in the contact area is determined such that the summation of 

their nodal forces is equal to the weight of the roller (see figure 70). Different numbers of 

elements are maintained in contact to balance the dead weight of the roller.

Figure 69. Illustration. Roller contact geometry for static indentation and during motion.

Figure 70. Chart. Indentation of a cylindrical roller into pavement. 
55

The plots presented in figure 71 represent the nodal reaction forces in the elements in response to 

the applied dead load of the roller as the contact area (represented by the number of elements in 

contact) is changed. A1 and A2 represent the edge and middle nodal forces of the same orientation

in an element when considering six elements in contact between the roller and the material. 

Similarly, B1 and B2 represent the respective nodal forces at the edge and middle nodes of an 

element when four elements are chosen in contact, and C1 and C2 are nodal forces when two 

elements are chosen in contact. On the basis of the nodal forces in the contact area, the pressure 

in each loading patch was determined according to FE procedures.

(56)

Figure 71. Chart. Change in nodal reaction forces as the load is applied over a smaller area.

Once the compaction load is determined, the load is translated over the pavement, providing a 

loading effect similar to sliding a roller over pavement. The exact mechanics to establish a 

rolling-type contact between the roller and the nonlinear pavement material require indepth

exploration of the field of nonlinear contact mechanics and has not yet been implemented in 

the FE model developed for this project.

The loading pattern was used to simulate the compaction of asphalt pavements with a roller travel

speed of 3.0 mi/h (4.8 km/h). This speed is in the range of values typically used in the field and 

reported in table 6. 

Table 6. Typical roller-speed ranges.

(1)

Type of Roller Breakdown Intermediate Finish

Static steel wheel 

2.0–3.5 mi/h

(3.2–5.6 km/h)

2.5–4.0 mi/h

(4.0–6.5 km/h)

3.0–5.0 mi/h

(4.8–8.1 km/h)

Pneumatic 

2.0–3.5 mi/h

(3.2–5.6 km/h)

2.5–4.0 mi/h

(4.0–6.5 km/h)

4.0-7.0 mi/h

(6.4–11.3 km/h)

Vibratory steel wheel 

2.0–3.0 mi/h

(3.2–4.8 km/h)

2.5–3.5 mi/h

(4.0–5.6 km/h) Not used 
56

The plots in figure 72 though figure 74 are indicative of the differences observed in the

displacement and stress response in an element in the direct path of the traversing load using

the currently developed contact-area algorithm and the previously available loading algorithm

in CAPA-3D, where the load is applied uniformly over all the finite elements in contact with the

roller. The plots are representative of the response of the material at a node directly under the 

load (experiencing maximum normal distress) at the beginning of the simulation. The stress 

response of the pavement while the new cylinder indentation algorithm (nonuniform pressure) 

is employed is characteristic of the nonuniform stress reactions obtained in a 20-noded element, 

whereas the algorithm previously in use strictly enforced uniform pressure at all nodes under a 

loaded area at a given instant of time. This, therefore, validates the choice of the use of the cylinder

indentation algorithm to better represent the physical contact situation between the roller and the 

pavement. The simulation results presented in this section utilize the model parameters for the 

SH-21 project simulations for the SGC (see Masad et al.).(49) The material parameters used for 

these comparisons are presented in table 7.

Figure 72. Chart. Comparison of vertical displacement for the two different patterns of 

loading.

Vertical Displacement, mm
57

Figure 73. Chart. Comparison of the normal stress distribution at the pavement’s 

top surface due to different loading patterns.

Figure 74. Chart. Comparison of the shear-stress (XY) distribution of the pavement’s 

top surface due to different loading patterns (on the plane of symmetry).

Table 7. Material parameters used to investigate the contact algorithm, interface effects, 

impedance element effects, and the effect of base stiffness.

Highway 

Project μμ� (MPa) n1 λ 1 q1 ηη� (MPa·s) n2 λ 2 q2

SH-21 2,500 4 0.25 −25 2,000 5 0.22 −30

1 MPa = 0.145 ksi
58

Normal Stresses Predicted by Model

The solution of a cylinder indented into an elastic medium (semi-infinite half space) was compared 

to the solution from the material model by balancing the self-weight of the roller with the normal 

reaction from the material. Then, the mean contact pressure used for the elastic solution was

calculated using the equation in figure 75.

ppmm = PP

ππaa2

Figure 75. Equation. Mean contact pressure.

In this equation, P is the load per unit thickness in the z-direction, and a is the radius of contact.

Here, the self-weight per unit of pavement width is used for P. The normal stress in the elastic 

medium is given by the equation in figure 76.

σσ2

ppmm

= − 4

ππ �1 − xx2

aa2�

1/2

Figure 76. Equation. Normal stress in elastic medium.

This representation of the solution indicates that the normal stress is expressed solely as a function

of the position x and the mean contact pressure Pm. This relationship is illustrated in figure 77 in 

comparison to the instantaneous normal stress produced in the HMA material used in the compaction

simulations with the same mean normal pressure applied. The comparatively higher stresses in the

elastic medium at a given mean contact pressure are due to the lack of any dissipative mechanism in

the elastic material. The model employed for representing the HMA develops dissipative responses 

in a nonlinear manner, and hence, the instantaneous response shows less stress. Also, the low level

of stress in the HMA shows the low-viscosity fluidlike precompaction behavior of the model. 
59

Figure 77. Chart. Comparison of the normal stress predicted by the model to the normal 

stress predicted for an elastic medium (for the same mean contact pressure, Pm).

Interface Elements

The FE analysis is dependent on the use of interface element layers between the top asphalt layer 

and the previously compacted asphalt layer as well as between the top newly laid asphalt and the 

old asphalt sublayer. This is due to the differences in properties of the loose mix and the stiff old 

asphalt layer. Hence, the specialized interface elements are employed to resolve the contact issues

that arise due to the large variation in properties of material on either side of the interface. With 

the establishment and implementation of a useful contact-area algorithm in CAPA-3D, the analysis

is performed to obtain insight into the effect of the interface elements. The use of interface elements

provides a numerically stable and converging simulation. These interface elements do not contribute

to the structural properties of the whole structure, and their applicability has been documented.

(54)

The interface elements to be used require the specification of two shear-stiffness and one normalstiffness moduli. These stiffness values are used in the numerical algorithm to ensure contact is

maintained between the layers by artificially attaching the corresponding nodes of the layers 

through a mathematical expression that governs their relative motion. The shear stiffness is used 

to govern the relative sliding of the layers, whereas the normal stiffness is used to prevent mesh 

penetration of the finite elements of a stiffer material into those of a soft material.

The results of FE analysis using the interface layer are presented in figure 78 through figure 81. 

The two shear-stiffness and one normal-stiffness properties are varied, and the resulting structural

behavior of the HMA material is observed. The shear moduli are individually made stiffer, and a 

comparison is drawn to a reference set of interface-element stiffnesses. The same comparison is 

made after normal stress is made less stiff. It can be observed in figure 78 through figure 81 that the

use of interface elements does not adversely affect the normal and in-plane displacement response

of the model. The slight deviation observed when z-direction shear is stiffened is attributed to the 

structure of the mesh. The mesh used is longer in the z-direction, and hence, the effects of free 
60

boundaries are more pronounced. The situation is rectified by using a higher stiffness value

because, beyond such a value, the interface exhibits similar behavior in the x-direction. 

Figure 78. Chart. Comparison of pavement vertical response as the normal stiffness of the 

interface elements is lowered by an order of magnitude from 1,451 ksi (10,000 MPa)

to 145 ksi (1,000 MPa).

Figure 79. Chart. Comparison of pavement response in the rolling direction as the 

normal stiffness of the interface elements is lowered by an order of magnitude 

from 1,450 ksi (10,000 MPa) to 145 ksi (1,000 MPa).

X Displacement, mm
61

Figure 80. Chart. Comparison of pavement response to increasing the shear stiffness 

in the x-direction.

Figure 81. Chart. Comparison of pavement response to increasing the shear stiffness 

in the z-direction.

Use of Impedance Layers

The impedance layers built into the FE mesh at this stage were not exploited for the current 

analyses since essentially nondynamic loading was studied. However, their impact on the FE

model is recorded, for the sake of completeness, through some simple examples. These elements

were included as part of the plan to study the effects of truly dynamic loads on the asphalt

pavement structure. The use of impedance layers is necessitated by the need to dampen the wave
62

reflections caused at the interface between layers when the pavement is subjected to dynamic 

loads. The techniques used to overcome the ill-posed nature of the initial value problem, which 

causes wave reflections at interfaces in highly viscous and elastic media, involve progressively 

dampening the displacement waves propagating through the medium through use of numerically 

cheap impedance layers. These layers are primarily composed of elastic material, with each 

enclosing layer proportionally less stiff in comparison to the inner one it is in contact with. This 

proportion, the “acoustic impedance,” is maintained to keep an interface between different materials

from reflecting waves that can interfere with the solution. Further details of the uses and applications

of such special FE layers have been documented by Sluys.

(55)

In the mesh employed, three impedance layers are used to encase the pavement structure, each 

1 inch (25.4 mm). The maximum value of the innermost layer is finally set at 725 ksi (5,000 MPa),

and the densities are calculated such that the ratio of moduli between consecutive layers is kept at 

10 (the cross-sectional area is fixed by the thickness of the layers). These values eliminate the wave 

reflections in the asphalt mix and hence are used for further analysis for all impedance layers.

The response plots that have impedance layers around the top fresh mix, old asphalt, and base (as 

shown in figure 82 and figure 83) are useful in reducing the vertical and horizontal reflections to 

the point that they are completely eliminated, as indicated in figure 84 and figure 85. Because of

this observation, a high enough impedance is sought when use of such layers is warranted, keeping 

the thickness of layers and the density of the material of the layers fixed.

Figure 82. Chart. Comparison of the dampening effect provided by impedance layers 

surrounding the structure laterally in the x-direction.
63

Figure 83. Chart. Comparison of the dampening effect provided by impedance layers 

surrounding the structure laterally in the y-direction.

Figure 84. Chart. Comparison of the horizontal dampening effect provided by impedance 

layers surrounding the top layer (pavement and old asphalt) laterally in the x-direction.
64

Figure 85. Chart. Comparison of the vertical dampening effect provided by impedance 

layers surrounding the top layer laterally in the x-direction.

Effect of Base Stiffness

FE analyses were performed by assigning properties for the base material in order to study the

effect of having a single base layer underneath the pavement. The results indicate that the 

material response does not seem affected much when the material properties of the single base 

layer are varied over a wide range with a minimum base stiffness value of 72.5 ksi (500 MPa)

and a maximum stiffness value of 290 ksi (2,000 MPa).

The plots of the volumetric response shown in figure 86 and figure 87 indicate that the material 

response is not very sensitive to changes in stiffness of the base when the response is viewed in 

an averaged sense. The volumetric response is examined through the evolution of the volumetric 

component of the viscous-deformation gradient (see Koneru et al. and Masad et al.).(48,49) However,

the FE model is responsive to a drastic increase in base modulus (as shown in figure 88 through 

figure 90) due to increasing rigidity in the base structure leading to interference waves (as a result

of structural rigidity and not due to any other dynamic characteristic) reflected internally. This is 

in contrast to what is observed in the field, where mix designs with stiff base materials produce 

more compaction in the pavement layer because of the compressive reaction acting on the HMA 

layer from below. As a result, for further FE analysis, a value between the two extremes (i.e.,

between the minimum base stiffness value of 72.5 ksi (500 MPa) and the maximum stiffness 

value of 290 ksi (2,000 MPa), for example 145 ksi (1,000 MPa)), is chosen to conduct the 

parametric and sensitivity analysis presented later. Vertical Displacement, mm
65

Figure 86. Chart. Comparison of the volumetric component of the viscous-evolution 

gradient for soft (72.5 ksi (500 MPa)) versus stiff (290 ksi (2,000 MPa)) bases.

Figure 87. Chart. Comparison of the volumetric component of the viscous-evolution 

gradient at two base-stiffness moduli of interest.
66

Figure 88. Chart. Comparison of the effect on the x-displacement of a node in the roller 

path as the base stiffness is varied from 72.5 ksi (500 MPa) (soft base) 

to 290 ksi (2,000 MPa) (stiff base).

Figure 89. Chart. Comparison of the effect on the y-displacement of a node in the 

roller path as the base stiffness is varied from 72.5 ksi (500 MPa) (soft base) 

to 290 ksi (2,000 MPa) (stiff base).

Vertical Displacement, mm
67

Figure 90. Chart. Comparison of the deflection for two base-stiffness moduli of interest.

Compaction Equipment Characteristics

As a further development to the simulation model, a vibratory-loading algorithm was developed 

for use in CAPA-3D specifically for the compaction project. The vibratory loads take into account

the effect of the displacement amplitude and frequency specifications of the machine used for 

rolling in a field project by calculating the equivalent load to be used to achieve the same level of 

compaction. Analytical expressions for calculating the contact forces for cylinder indentation

into an elastic pavement material are adapted from Johnson.(57) Thereafter, the equivalent loads 

are calculated for different amplitudes of indentation as necessary. Then, the force necessary for 

the vibrations is calculated using these analytical expressions at an instant of time, and the pressure

is redistributed over the contact area.

To summarize, the algorithm has the following features:

The loading area is the same as the one determined from the static indentation.

The amplitude of the half-sine, vibratory load is determined based on vibratory-system 

specifications.

The time taken for each load cycle is adjusted to match the loading frequency. 

Simulations were set up on a mesh representative of a 12-ft (3.66-m)-wide pavement lane with a 

roller weight of 27,000 lb (12.258 kg), a diameter of 51 inches (1,295.4 mm), and a width of 7 ft

(2.135 m). The results of the simulations (shown in figure 91 through figure 94) show that the 

material responds according to the mechanics of the loading algorithm implemented. The results 

are presented at a point in the interior of the material along the width (approximately 5 ft (1.525 m)

from the edge). In the results that follow, the term compaction is used as a measure to represent Vertical Displacement, mm
68

the viscous evolution, det(G)-1 (see Masad et al. for clarification of the term).(49) The effect of a 

cyclic load vibrating at a certain angular frequency on the deformation induced on the material 

model of choice, with a certain relaxation time associated with it, is also presented (see figure 94).

The following observations can be made regarding the response of the material:

An increase in the amplitude of the vibratory load leads to an increase in the amount of 

compaction achieved (see figure 91).

The increase in compaction as a result of the increase in amplitude is due to the roller 

indenting deeper into the material, thereby increasing the pressure on the surface.

As shown in figure 92 and figure 94, increasing the frequency of the vibration or decreasing 

the material’s viscous nature results in an increase in the percent of compaction of the 

material, as is observed in the field.

The response of the material to small steps of increase in load (keeping the roller area of

contact constant) results in a uniform increase in compaction level, as indicated in figure 93. 

Figure 91. Chart. Compaction over a sequence of passes as the amplitude of vibration 

increases.
69

Figure 92. Chart. Compaction over a sequence of passes at different frequencies.

Figure 93. Chart. Material response to change in dead load carried by each roller.
70

Figure 94. Chart. Field compaction response at a constant frequency over multiple passes 

on a point. 
71

CHAPTER 7. SIMULATION OF FIELD COMPACTION

SENSITIVITY ANALYSIS

The purpose of a sensitivity analysis is to correlate various features of the output of the

mathematical model to the different input factors and parameters of the model. The sensitivity 

analysis is constrained by a lack of experimental evidence for the mechanical response of HMA 

while compaction is being performed in the field. In order to establish a reference material response

and aid in the sensitivity study, a reference parameter set was chosen. Since compaction of HMA 

manifests primarily in the volumetric changes that occur during the process, this study is based 

on a volumetric measure that quantifies material response. A comparison of the general volumetric

viscous evolution response of the material (measured through the mathematical invariant det(G)-1

that represents the volumetric change in the viscous evolution tensor G) and its response when 

using a reference set of material parameters is presented.

The sensitivity analysis provides a general qualitative understanding of the trends of the material 

behavior in relation to the variation in parameters. The viscous evolution of the model as the 

parameters are varied is presented in figure 95 and figure 96. The model parameters used to

study sensitivity are presented in table 8.

Figure 95. Chart. Evolution of the volumetric viscous gradient with a change 

in values of individual parameters μμ̂, ηη̂, λ 1, and λ 2.
72

Figure 96. Chart. Evolution of the volumetric viscous gradient with a change 

in values of individual parameters n1, n2, q1, and q2.

Table 8. Parameters employed for the sensitivity study of the material.

Mix Parameter Sets μμ� (MPa) n1 λ 1 q1 ηη� (MPa·s) n2 λ 2 q2

Reference set 620 4.0 0.25 −15 1,400 2.5 0.25 −20 

Set 1 (change ηη̂) 620 4.0 0.25 −15 1,700 2.5 0.25 −20 

Set 2 (change μμ̂) 810 4.0 0.25 −15 1,400 2.5 0.25 −20 

Set 3 (change λ1) 620 4.0 0.20 −15 1,400 2.5 0.25 −20 

Set 4 (change λ2) 620 4.0 0.25 −15 1,400 2.5 0.2 −20 

Set 5 (change n1) 620 3.0 0.25 −15 1,400 2.5 0.25 −20 

Set 6 (change n2) 620 4.0 0.25 −15 1,400 1.5 0.25 −20 

Set 7 (change q1) 620 4.0 0.25 −10 1,400 2.5 0.25 −20 

Set 8 (change q2) 620 4.0 0.25 −15 1,400 2.5 0.25 −15 

1 MPa = 0.145 ksi

Observations from the Sensitivity Analysis

The response of the material is stiffer (see figure 95) as the viscosity parameter ηη̂is increased. 

This correlates to the response observed in the SGC simulations and is understood to be due to 

the material becoming a more viscous fluid as compaction progresses.

The shear-modulus parameter μμ̂causes a significant stiffening of the material in the initial stages

of compaction. This can be observed in the response of the material in the first 50 time steps of 

the simulation. The comparative increase in the elastic rebound experienced by the mix in the 

initial stages, as μμ̂increases, corresponds to the constitutive assumption that the stored energy 
73

in the material is dependent on the shear-modulus function. Thus, as in the SGC simulations, the 

shear modulus affects the material response, but the viscous nature of the material dominates the 

response as the load traverses the pavement surface. A decrease in the λ 1 and λ 2 values also serves

to make the material stiffer. λ 1 has an effect on the model response only at the very beginning of 

compaction. Even at the initial stages, the influence of λ 2 is just as significant, thereby drawing 

attention to the response through all compaction stages. After the initial stages, the response to 

changes in λ 2 is dominant, as is the case in SGC compaction. As shown in figure 97, the regions 

of the parameters influence the gyratory compaction. This observation is useful for drawing

generalized correlations between the material response expected in SGC compaction modeling and

that expected in field compaction using the present model. As shown in figure 95, the increase in q2

results in the pavement responding significantly differently than when q1 is increased individually.

This points toward q2 being the control parameter in a nonlinear transition of the material from the

initially loose mix in the forward pass to a stiffer mix in the return pass. The lesser the magnitude

of q2, the less permanent compaction is achieved. The pavement exhibits a relaxation mechanism 

as shown in figure 95 and figure 96. Note that the material is again compacted on the return pass,

midway through the total number of time steps. 

Figure 97. Chart. Regions of influence of model parameters in gyratory compaction.

(51)

Parametric Analysis

The parametric analysis focused on the influence the parameters found not to affect the SGC curves 

(λ 1, q1, and n2) had on field compaction. These results are presented in figure 98 through figure 103

by plotting the viscous deformation parameter (det(G)-1) for different model parameters.
74

Figure 98. Chart. Evolution of the volumetric viscous gradient with a change in λ1.

Figure 99. Chart. Evolution of the volumetric viscous gradient with a change in q1.
75

Figure 100. Chart. Evolution of the volumetric viscous gradient with a change in n2.

Figure 101. Chart. Evolution of the volumetric viscous gradient with a change in λ2.
76

Figure 102. Chart. Evolution of the volumetric viscous gradient with a change in q2.

Figure 103. Chart. Evolution of the volumetric viscous gradient with a change in n1.
77

Correlation with Laboratory Compaction

Researchers sought to draw correlations between the behavior of model parameters during the 

SGC simulation and their behavior during the simulation of field compaction. The results were 

as follows:

The profiles of the deflection and volumetric part of viscous evolution indicate that the 

model response is not sensitive to λ 1 (see figure 98), and hence, the parameter can be a 

model constant. Therefore, λ 1 can be the same constant as in the SGC simulations.

The response observed from varying q1 indicates that the model is sensitive to this 

parameter, as shown in figure 96 and figure 99. The volumetric part of the viscous 

evolution exhibits no sensitivity to q1.

Because of the response exhibited by the model to q1 in terms of the thermodynamic 

quantity det(G), q1 can be considered a model constant. Therefore, the same constant 

used in the SGC simulations can be used.

Figure 96 and figure 100 show that the model exhibits considerable sensitivity to changes 

in n2. This response from the model shows deviation from the response exhibited by the 

model during the simulation of the SGC.

Parameter n2 is understood to be related to the model’s initial rate of dissipation (at the 

low-viscosity, loose-mix stage of the compaction). The significance of n2 for field 

compaction is attributed to the unbound nature of the pavement being simulated. This 

causes the material in the loose-mix state to be more sensitive to applied stresses as 

compared to material in a confining mold in the SGC.

Parameter n2 is calibrated from available field data.

Figure 101 through figure 103 indicate that the significant parameters in the SGC 

simulation are also significant for field compaction.

The behavior of parameter n1 should be examined. While n1 ≤ 4.0, the model does not 

exhibit any sensitivity to this parameter, whereas at a higher value, the model exhibits 

significantly higher stiffness, as shown in figure 96 and figure 103.

A correlation exists in parametric behavior for both compaction processes.

ANALYZING LONGITUDINAL JOINTS

Asphalt pavements close to longitudinal joints tend to be compacted less than the center of the 

pavements. This is due to the tendency to apply fewer passes at the joints. The low confinement 

at some types of joints (unrestricted or unconfined joints) and the higher rate at which the mixture

at the joint loses heat reduce the efficiency of compaction at the joint compared with the pavement

center. Here, the compaction of longitudinal joints is studied by simulating the process using the 

model developed in finite elements. FE simulations were conducted to study the effects of
78

longitudinal joints. Figure 63 shows the boundary conditions that are to be applied on the edges 

that have been fixed and on the edges that have been unrestrained, along with the longitudinal 

joints of the pavement. The model predicts a higher level of compaction close to the fixed edge 

of a lane (0.5 ft (0.1525 m)). The compaction predicted at the outer part of the lane, close to the 

unconfined/free edge (0.5 ft (0.1525 m)), is significantly lower because of the lack of a confining 

pressure close to the edge. The overlap zone between two roller passes contains the longitudinal 

joint where material shoving occurs to accommodate the compaction of the material. This kind 

of behavior is also observed in the field and is depicted in figure 104.

Figure 104. Chart. Plot representing the final compacted state of the material 

along the width of the pavement. 

VERIFICATION OF SIMULATION MODEL

Highway SH-21

A section of highway SH-21 was compacted using a vibratory roller compactor to provide data 

for preliminary evaluation of the model developed in this study. The compactor had two steel 

drums with a width of 7 ft (2.135 m) and a distance of 12 ft (3.66 m) between the centers of the 

drums. The total weight of the roller was 27,783 lb (12,613 kg); the front drum had a weight of 

13,980 lb (6,347 kg), and the back drum had a weight of 13,803 lb (6,266 kg). The roller was

operated in vibration mode with about 3,000 vibrations per minute at a speed of 2 to 3 mi/h (3.2 to 

4.8 km/h). Aggregate characteristics (angularity, texture, and sphericity) were measured using the 

Aggregate Imaging System. Higher numbers for angularity and texture indices mean higher 

aggregate angularity and texture. Higher sphericity values indicate that particles are less flat and 

elongated (a sphere has a value of 1). The properties of the various layers are shown in table 9.
79

Table 9. Material properties used for the SH-21 project.

Layer Modulus, E psi Poisson Ratio, n

2.0-inch Type C HMA Refer to table 10

2.0-inch of Type D HMA 450,000 0.30

18.0-inch flexible (granular) base 30,000 0.35

1 inch = 25.4 mm

The FE model was used to simulate field rolling compaction on material 12-ft (3.66-m) wide with

the roller centered on the material, as shown in figure 105. The initial simulations were conducted

using the material model parameters from SGC simulations, shown in table 10. However, the

percent of compaction was much lower than the experimental measurements in terms of percent 

air void (%AV). Computations of percent compaction close to field air voids were obtained only 

by reducing the viscosity coefficient ˆ to a value less than 72.51 ksi·s (500 MPa·s). The results 

shown in figure 106 are for ˆ of 29 ksi·s (200 MPa·s). The laboratory compaction curves were 

described well using a ˆ value of 290.07 ksi·s (2000 MPa·s), as shown in table 10.

The measured %AV is shown in figure 106, and the change in %AV is shown in figure 107. A

comparison of the typical simulated responses using the FE model shows that the model developed

predicts a trend of compaction over multiple passes similar to trends measured in the field (see 

figure 108 through figure 110).

Figure 105. Illustration. Schematic of a roller on a material with three locations 

for density measurements.

(49)

Table 10. Model parameters used for projects SH-21, US-87, and US-259.

Highway 

Projects

Parameter 

Sets μμ� (MPa) n1 λ1 q1 ηη� (MPa·s) n2 λ2 q2

SH-21

SGC 2,500 4 0.25 −25 2,000 2.5 0.22 −30

Field 2,500 4 0.25 −25 2,00 2.5 0.22 −30

US-87

SGC 2,600 5 0.25 −25 2,100 2.5 0.26 −27 

Field 2,300 5 0.25 −25 2,100 2.5 0.26 −27 

US-259 

SGC 2,400 5 0.25 −25 1,700 2.5 0.21 −28 

Field 2,150 5 0.25 −25 1,700 2.5 0.21 −28 

1 MPa = 0.145 ksi
80

Figure 106. Chart. Measurements of the %AV in the asphalt mix.(49)

Figure 107. Chart. Measurements of change in %AV in the asphalt mix.(49)

0.00

2.00

4.00

6.00

8.00

10.00

12.00

14.00

16.00

18.00

0 1 2 3 4

Pass Number

Percent Air Voids %

Measurement

Model

Figure 108. Chart. Measurements and modeling results of %AV at point A of the pavement 

locations shown in figure 105.
81

0.00

2.00

4.00

6.00

8.00

10.00

12.00

14.00

16.00

18.00

0 1 2 3 4

Pass Number

Percent Air Voids %

Measurement

Model

Figure 109. Chart. Measurements and modeling results of %AV at point B of the pavement 

locations shown in figure 105.

0.00

2.00

4.00

6.00

8.00

10.00

12.00

14.00

16.00

18.00

0 1 2 3 4

Pass Number

Percent Air Voids %

Measurement

Model

Figure 110. Chart. Measurements and modeling results of %AV at point C of the pavement 

locations shown in figure 105.

The difference between the simulation parameters obtained from the laboratory and field compaction

processes raises important points regarding the modeling aspects that need to be addressed to 

accurately describe both processes. One aspect is the accuracy of representing the boundary and 

loading conditions of the laboratory and the field. For example, the load used in simulating field 

compaction is quasi-static and does not accurately represent the actual dynamic forces applied by 

a vibratory roller. There is also the mathematically and computationally complex contact problem

between the roller drum and the asphalt pavement surface that is not addressed in the FE simulations.

Another important cause of the differences in the model’s parameters is the material model’s limited

ability to represent the various mechanisms involved in the laboratory and field compaction 

processes. This cause is less obvious and harder to address compared with the differences in

boundary conditions. For example, the asphalt mix behavior is rate dependent, but the model might

be limited in capturing the rate dependency. This rate dependency would not significantly affect 

the SGC simulations because the compaction is conducted under continuous application of one 

loading rate of 30 revolutions per minute. Such a limitation, however, would be clearly manifested

in field simulations because the rate of loading could vary considerably. Consequently, calibration 

with field compaction data would be required as the field conditions could vary beyond the 

laboratory experiments used to determine the model’s parameters.
82

Highway US-87

In order to further examine the utility of the model in simulating field compaction trends, field 

rolling compaction was simulated to match the rolling patterns in two highway projects, a US-87 

paving job near Port Lavaca in Calhoun County, TX, in October 2006 and a section of US-259 

located in Rusk County, TX, in February 2007.

(7) Specimens from both projects were obtained 

from the wheel path, between the wheel path, and on the longitudinal joint (restrained and

unrestrained). The asphalt mixtures and compaction data were part of a study funded by the 

Texas Department of Transportation to evaluate the influence of various field compaction 

methods on asphalt-mixture properties. The properties of the top asphalt mixture are reported 

in table 11 and table 12.

Table 11. Summary of mixture designs.

Highway 

Project

Mixture 

Type

Aggregate 

Type

Binder 

Grade

AC 

(Percent) Gmm VMA

Design 

Air Voids 

(Percent)

SH-21 Type C Limestone PG 70-22 4.7 2.467 14.3 3.0

US-87 Type C 

Siliceous 

River Gravel PG 76-22S 4.3 2.439 13.8 4.0 

US-259 Type C 

Sandstone and 

Limestone PG 70-22S 4.3 2.478 13.1 3.0 

AC = Asphalt content

Gmm = Maximum specific gravity

PG = Performance grading

VMA = Voids in the mineral aggregate

Table 12. Summary of properties of mixture constituents.

Highway 

Project

Binder Viscosity, 

Pa·s (135 °C)

Aggregate 

Texture 

Index

Aggregate 

Angularity 

Index

Aggregate 

Sphericity

SH-21 0.883 106 2,811 0.708

US-87 2.258 112 3,062 0.651 

US-259 0.818 189 2,791 0.618 

1 Pa = 0.00014 psi

US-87 is a four-lane divided highway. Test sections were located on the northbound outside lane. 

The mixture was Type C (Texas Department of Transportation 1993 specification) designed with 

Fordyce Gravel and Colorado Materials limestone screening with 4.3 percent PG 76-22S binder. 

The Type C mixture is similar to the coarse-graded Superpave® mixture. The Type C mixture was

laid on a Type B material, which primarily included crushed river gravel. The thickness of the 

Type C layer was 2 inches (50 mm). The mixture was laid with a 16-ft (4.8-m) material width,

with 1.5 ft (4.5 m) tapered on one side. The structure of the pavement for the US-87 pavement 

project is presented in figure 111, the layer properties are given in table 13, and the width of the 

top layer of material is shown in figure 112.
83

Figure 111. Illustration. Pavement structure for the US-87 project.

Table 13. Material properties used for the US-87 project.

Layer Modulus, E psi Poisson Ratio, n

2.0-inch Type C HMA Refer to table 10

3.5-inch Type B HMA 450,000 0.30 

6.0-inch flexible (granular) base 30,000 0.35 

6.0-inch lime-treated subgrade 

(5 percent lime) 12,000 0.45 

1 inch = 25.4 mm

Figure 112. Chart. Schematic for the rolling patterns for the US-87 project.
84

The general rolling pattern can be described as breakdown by a steel-wheel vibratory roller and 

pneumatic wheel roller for both intermediate and finish rolling. In addition, the rolling pattern 

consisted of progressively moving the vibratory steel-wheel roller in transverse directions, with 

the pneumatic wheel roller again acting as both the intermediate and finish roller. The sequence 

and locations of the roller are simulated approximately, along with the boundary conditions

representative of the restrained and unrestrained edges of the pavement, and are indicated by the 

bars in the schematic in figure 112. In the schematic, the inner edge of a lane is indicated by the 

vertical dotted line passing through 0. The roller locations are represented with respect to this 

line as the datum to measure distance. The line segments represent rollers with their rolling 

directions. An upward arrow indicates forward rolling, and a downward arrow indicates the reverse.

The rolling pattern and the measured %AV are presented in table 14. There were 11 total passes 

in the compaction process involving the vibratory roller, with 9 passes in vibratory mode and 

2 passes in static mode.

Table 14. Rolling pattern and %AV measured in the field for the US-87 project.

Core Group

Distance from 

Edge (ft)

Vibratory 

Mode Static Mode %AV

1 1 3 1 9.65

2 4 5 1 6.77

3 7 3 1 7.33

4 10 6 2 5.01

1 ft = 0.305 m

Highway US-259

The overlay in the US-259 project used a Type C surface mixture compacted in a 2-inch (5-cm)

lift thickness. The coarse part of the aggregate was sandstone, and the intermediate- and fine-size 

particles were limestone. The mix had 11 percent field sand and 4.3 percent PG 70-22S binder. 

The test sections were in the southbound outside lane. Type C mix was laid on top of a recently 

compacted Type D level-up course where the Type D mixture was similar to a Superpave®

fine-graded mixture. The paving width was approximately 15 ft (4.57 m) (including shoulder),

and a vertical longitudinal joint was maintained. The structure of the pavement for the US-259

pavement project is presented in figure 113, the layer properties are given in table 15, and the 

width of the top layer of material is shown in figure 114.

Figure 113. Illustration. Pavement structure for the US-259 project.
85

Table 15. Material properties used for the US-259 project.

Layer Modulus, E psi Poisson Ratio, n

2.0-inch Type C HMA Refer to table 10

1.25-inch Type D HMA 450,000 0.30

7- to 8-inch asphalt concrete pavement 400,000 0.30

10-inch flexible (granular) base 34,000 0.35

1 inch = 25.4 mm

Figure 114. Chart. Schematic for the rolling patterns for the US-259 project.

The general rolling pattern can be described as breakdown by a steel-wheel vibratory roller and 

pneumatic wheel roller for both intermediate and finish rolling. In addition, the rolling pattern 

consisted of progressively moving the vibratory steel-wheel roller in transverse directions, with the

pneumatic wheel roller acting as both the intermediate and finish roller. The sequence and locations

of the roller are simulated approximately, along with the boundary conditions representative of 

the restrained and unrestrained edges of the pavement, and are indicated by the bars in the schematic

in figure 114. In the schematic, the inner edge of a lane is indicated by the vertical dotted line 

passing through 0. The roller locations are represented with respect to this line as the datum to 

measure distance. The line segments represent rollers with their rolling directions. The rolling 

patterns and the %AV are presented in table 16. There were nine total passes in the compaction 

process involving the vibratory roller, with five passes in vibratory mode and four in static mode.
86

Table 16. Rolling pattern and %AV measured in the field for the US-259 project.

Core Group

Distance from 

Edge (ft)

Vibratory 

Mode Static Mode %AV

1 1 1 2 11.26

2 5 3 3 7.70

3 8 3 2 8.12

4 11 2 2 10.78

5 14 2 1 9.27

1 ft = 0.305 m

Observations from Field Compaction Simulations for US-87 and US-259

As in actual field compaction, the initial air voids generally vary between 15 and 20 percent of the

volume of the mix, and the measured %AV is adjusted by using these two initial air-void values 

to calculate the change in measured %AV, which is more indicative of the compaction the material

has experienced. The change in measured %AV is then compared directly to the percent-compaction

values predicted by the FE model (see table 17 and table 18). A compaction zone is enclosed by 

the use of two initial %AV values, as indicated in figure 115 and figure 116, where the measured 

and calculated percent compaction at different locations across the material are compared at the 

end of rolling pattern cycles. In the figures, the comparisons are made per core group, which

represents different locations across the material relative to the edge (see table 9 through table 16).

Table 17. Change in measured %AV in the field and calculated percent compaction

for the US-87 project.

Core Group

15 Percent

Initial 

Air Voids

20 Percent

Initial 

Air Voids

Compaction Using 

SGC Parameters 

(Percent)

Compaction Using 

Field Parameters 

(Percent)

1 5.35 10.35 3.52 7.50

2 8.23 13.23 5.17 10.39

3 7.67 12.67 4.42 9.49

4 9.99 14.99 6.04 12.83

Table 18. Change in measured %AV in the field and calculated percent compaction 

for the US-259 project.

Core Group

15 Percent

Initial 

Air Voids

20 Percent

Initial 

Air Voids

Compaction Using 

SGC Parameters 

(Percent)

Compaction Using 

Field Parameters 

(Percent)

1 3.74 8.74 4.02 7.48

2 7.30 12.30 5.92 10.95

3 6.88 11.88 3.36 7.25

4 4.22 9.22 2.64 5.28

5 5.73 10.73 3.26 7.22
87

Figure 115. Chart. Comparison of the total percent compaction from simulations with the 

general trend of the %AV measured at the end of the field compaction process for US-87.

Figure 116. Chart. Total percent compaction from simulations compared to the general 

trend of the %AV measured at the end of the field compaction process for US-259.
88

The percent-compaction values were calculated using parameters from SGC compaction of these 

two mixes as well as field-adjusted parameters. The compaction predicted by the SGC parameters 

was outside the range of the change in %AV measured in the field (see figure 115 and figure 116).

Therefore, the parameters were adjusted so that the compaction obtained in simulations was

contained within the range of measured %AV values. The parameters used for simulation purposes

are shown in table 10. As shown in figure 115 and figure 116, within the compaction zone the 

behavior of the mix in the simulations correlates well with the trends observed in the field. The

type of compaction in each project (based on the combination of the type of roller and the rolling 

pattern) correlates directly to the measured change in %AV for the corresponding core groups (see

table 17 and table 18). The greater the measured change in %AV, the more compaction occurs at 

that location. This behavior is reflected in figure 115 and figure 116 for simulations of both projects

and, in that sense, agrees with reality.

The simulation model is capable of providing further insight into field compaction trends when 

using the measured %AV to gauge material response. As can be seen in figure 117 (comparing 

the response of two materials at the same location relative to the edge of the material), the

simulations predict that the US-87 material will undergo more compaction by the end of the 

whole process than the US-259 material. This is in direct correlation to the higher change in 

measured %AV for the US-87 pavement (table 17, core group 1) as compared to the change in 

measured %AV for the US-259 pavement (table 18, core group 1).

Figure 117. Chart. Comparison of prediction of percent compaction per roller pass

for US-87 and US-259. 
89

Some interesting observations can be drawn from the response of the material as it is subject to 

different roller passes. From figure 118 (with the roller patterns from table 14), the following 

observations of the material behavior of US-87 pavement undergoing compaction can be made:

Core group 4 received more passes and accordingly underwent the most compaction. 

Core groups 1 and 3 received the same types of load and the same number of passes; 

however, group 3 underwent more compaction. This is attributed to the rolling pattern 

chosen for this project, which subjected the material in this region to more overlapping 

loads by virtue of its centralized location along the material width. 

Core group 2 underwent compaction in between group 4 and groups 1 and 3.

Figure 118. Chart. Prediction of percent compaction per roller pass across the material 

for US-87 (cores taken at four locations).

From figure 119, the following observations for the compaction response of the pavement 

material for US-259 can be made:

Core group 2 was subjected to the most vibratory and static-load passes and was centrally 

located; hence, it underwent the most compaction. 

Core group 1 was subjected to rapid compaction when static rollers were brought on after 

it was subject to only one vibratory loading. 

Core groups 4 and 5 were influenced by the same number of vibratory loads, with group 5

subjected to fewer static loads. However, group 5 achieved more compaction. This is in 

agreement with the field compaction as evidenced by the %AV for groups 4 and 5 in table 16.
90

Group 3 underwent uniform loading (in terms of the percent compaction per pass) as a 

result of the vibratory loads and less static loading (compared to group 2) due to its 

central location in the material. 

Figure 119. Chart. Prediction of percent compaction per roller pass across the material 

for US-259 (cores taken at four locations).

In conclusion, researchers observed that compaction simulated across the pavement over the

different roller passes correlates with the compaction measured (through a change in %AV) in 

the field. The compaction simulation using model parameters from the SGC simulations is similar

to the measured values. These parameters, with a simple shift in one of the parameter values (μμ̂),

can therefore be used to predict field compaction, allowing a correlation between laboratory and

field compaction simulations. A compaction-zone characteristic of each mix type can be predicted.

Applying vibratory-mode loads initially causes a more uniform compaction (in terms of the percent 

compaction per pass). Rolling patterns and load locations are of utmost importance in determining 

the amount and type of compaction (more uniform or less uniform). The simulation model provides

an opportunity to experiment with different rolling patterns to achieve a specific percent compaction

or a certain trend in the uniformity of compaction from pass to pass. Considerable variability exists

in the compaction trends. This has been shown through the study of two projects. The variability in 

trends due to loading patterns is a natural reflection of the nonlinear nature of the materials 

employed. Only vibratory and static wheel loads have been considered. Simulation of finish 

rolling through use of pneumatic rollers was not considered. 
91

CHAPTER 8. SUMMARY AND CONCLUSIONS

The main reasons for the selection of the thermodynamics-based nonlinear viscoelastic model are

as follows:

The compaction process involves large strains, which can be modeled using the adopted 

thermodynamic framework. 

The asphalt mix experiences significant microstructure changes during the compaction 

process. These changes include aggregate reorientation, reduction in air voids, and increase

in aggregate contacts. The continuum model was formulated to account for these phenomena

through the specification of appropriate evolution functions of the Helmholtz potential 

and rate of dissipation. 

A few parameters are used in calibrating the model. These parameters are the ones used 

in determining the evolutions of Helmholtz potential and rate of dissipation. 

Continuum models similar to the one that was developed in this study have been used 

successfully in modeling asphalt-material behavior under complex loading and

environmental conditions. The temperature of asphalt mix changes during compaction,

and this phenomenon needs to be accounted for in the model in future research. 

FE IMPLEMENTATION

The FE implementation of the material model enables numerical investigations into the behavior 

of the material when subjected to compaction. The following points summarize the main findings 

from the FE simulations:

The model was implemented in CAPA-3D using the FE method. The implementation 

was then validated by comparing the FE solution for one-dimensional constant stress and 

constant strain loading with the corresponding solutions from analytical techniques (aided

by calculations using MATLAB®). 

The model exhibits a nonlinear response in shear through the exhibition of normal stress 

differences. 

The continuum model lends itself to implementation in the FE method to simulate the 

compaction process under various conditions. Therefore, it can be used to relate the 

laboratory compaction process to the field compaction process.

The material model was implemented successfully in the FE package CAPA-3D and was 

used to simulate compaction in the SGC and in the field under various loading and 

boundary conditions.
92

SIMULATION OF GYRATORY COMPACTION

The FE method was used to simulate the laboratory compaction of asphalt mixtures. The main 

findings from the simulation of gyratory compaction are as follows:

The parametric analysis conducted in this study showed that the compaction process is 

insensitive to changes in some of the model’s parameters. Therefore, researchers used

constant values for these parameters during all the FE simulations in order to simplify the 

process of finding the values of the remaining parameters. However, the model needs to 

be further examined for other processes in which these parameters might be significant.

The parametric analysis demonstrated that some of the model’s parameters are related 

mostly to the initial stage of compaction, when the material exhibits low-viscosity 

fluidlike behavior, while other parameters are related to the behavior of the mixture after 

it starts to exhibit a high-viscosity fluidlike behavior and the compaction rate decreases.

A systematic method was developed to determine the model’s parameters from the SGC

curves. It was necessary to develop this method because of the difficulty of conducting 

conventional tests (e.g., triaxial or shear tests) to determine the model’s parameters given 

the high compaction temperatures.

The FE simulations demonstrated the capability of the constitutive model to simulate the 

Superpave® gyratory compaction process at different compaction angles.

The developed model and FE implementation will allow the prediction of mixture 

compactability in the field based on laboratory measurements. It will also allow the study 

of the influence of changes in material properties and mixture designs on mixture

compactability under various laboratory and field conditions. 

SIMULATION OF FIELD COMPACTION

The details of the FE model with regards to the mesh, geometry, loading methodology, and 

boundary conditions were presented earlier in this report. The primary findings from the 

simulation of field compaction are as follows:

A loading algorithm was developed for use in CAPA-3D for pavement compaction 

simulations. The loading algorithm utilizes nonuniform pressure loads applied over the 

surface of the finite elements in contact with a roller.

Moving loads were applied in a quasi-static manner by incrementally translating the 

necessary element surface loads during the simulation.

The vibratory loads were also applied in a quasi-static manner. The inertial effects of the 

asphalt mix layer were not considered in these simulations.

Interface-layer elements were utilized to resolve contact issues pertaining to the layer 

interactions.
93

The FE model developed in this study is capable of capturing the influence of the boundary

conditions (fixed/free) near the longitudinal joint on the compaction process. Consistent with

the field experience, fixed longitudinal joints exhibited higher compaction than free edges.

The FE model was validated against data from three field projects. The FE results agree 

well with the general trends of compaction observed in the field.

The FE model is capable of providing insight into which mix design and rolling pattern 

combination provides the most compaction at different locations across the pavement.

RECOMMENDATIONS FOR FUTURE INVESTIGATIONS

This study can be extended in several ways. The model developed can be used to develop a database

of parameters for different mix designs. This database of model parameters can be used for further

numerical investigations to better understand the IC process.

The material model presented in this report has the capability to include the effect of temperature 

on the material behavior. However, due to the lack of experimental data on the influence of 

temperature on material properties within the range of compaction temperatures, this effect was 

not accounted for in the simulations conducted in this study. Conducting experimental measures 

is recommended to monitor temperature changes during compaction and to measure the effects 

of these changes on material properties. In addition, the model needs to be expanded to simulate 

nonisothermal conditions during the compaction process.

The material’s model parameters were related qualitatively to the characteristics of the SGC 

curves. This approach was followed due to the lack of experimental methods to measure the 

asphalt mixture properties at the range of temperatures involved in the compaction process. 

Future research should focus on the development of a consistent experimental program to 

quantitatively assess the model’s parameters.

The continuum model implemented in the FE method is a powerful approach for simulating the 

material and structure response under various loading conditions. The model’s parameters are 

obtained from experiments on the asphalt mixtures. The relationships between the model’s 

parameters and the properties of the mixture constituents can be estimated by observing changes

in the model’s parameters and material response due to changes in mixture constituents. However,

it cannot be used to directly determine the influence of changes in material properties on compaction

characteristics. Researchers at the Turner-Fairbank Highway Research Center have worked on the

development of a micromechanical model for asphalt mixture compaction. Such micromechanical

models are computationally intensive and do not lend themselves to simulations of field

compaction at various compaction operations and field conditions. However, the advantage of

the micromechanical approach is its ability to directly account for the characteristics and 

properties of mixture constituents in modeling compaction.

A combination of the continuum and micromechanical models could be used to directly estimate 

the continuum model’s parameters using the micromechanical model’s results. This can be 

achieved by reformulating the constitutive model developed in the framework of this project 
94

along the lines of multiplicative split “multislip” models that simulate, on a continuum level, the 

micromechanical interactions that occur at the constituent materials level. 

As shown in figure 120, the new model will utilize features of the micromechanical response of 

an asphalt mix (e.g., dilation, aggregate rotation, etc.) that are controlled by the compositional 

characteristics of the mix (e.g., gradation, angularity, etc.) to determine the in-time development 

of material compaction as represented by the inelastic deformation gradient G of the continuum 

constitutive model.

G sn nn =g Ä +d Ä

δ

g θ

1

3

s

n

X,x 1 1

X ,x 2 2

X ,x 3 3

0

dx

F

Fe

G

r

G sn nn =g Ä +d Ä

δ

g θ

1

3

s

n

X,x 1 1

X ,x 2 2

X ,x 3 3

0

dx

F

Fe

G

r

X,x 1 1

X ,x 2 2

X ,x 3 3

0

dx

F

Fe

G

r

Figure 120. Illustration. Micromechanical response of an asphalt mix.

The micromechanical model will be used to directly examine the influence of changes in 

materials (binder viscosity, aggregate shape characteristics, aggregate gradation) on compaction. 

Then, it will be possible to develop relationships that relate the continuum model’s parameters to 

the compositional characteristics of the materials. This approach is illustrated in figure 121.
95

Figure 121. Chart. Representation of tasks involved in modeling asphalt mix compaction.

Micromechanical

model 

simulations

Finite element 

simulations using 

continuum 

model

Relationships between continuum model parameters 

and mixture compositional properties

Simulate field and laboratory compaction using continuum model

Relationships between compositional properties and field 

compaction processes

Simulate laboratory compaction

97

REFERENCES

1. The Asphalt Institute. (1989). The Asphalt Handbook, Lexington, KY. 

2. SIMULIATM. (2007). Abaqus User’s Manual, Dassault Systèmes, Providence, RI.

3. Brown, S.F., and Snaith, M.S. (1978). “The Permanent Deformation Characteristics of a 

Dense Bituminous Macadam Subjected to Repeated Loading,” Proceedings of the 4th 

Conference on Structural Design of Asphalt Pavements, 225–248, University of Michigan, 

Ann Arbor, MI.

4. Krishnan, J.M. and Rajagopal, K. (2003). “Review of the Uses and Modeling of Bitumen 

from Ancient to Modern Times,” Applied Mechanics Reviews, Vol. 56, No. 2, 149–214.

5. U.S. Army Corps of Engineers. (2000). Hot-Mix Asphalt Paving Handbook, Washington, DC. 

6. Roberts, F., Kandhal, P., Brown, E.R., Lee, D., and Kennedy, T. (1991). Hot Mix Asphalt 

Materials, Mixture Design and Construction, National Asphalt Pavement Association, 

Research and Education Foundation, Lanham, MD.

7. Kassem, E.A.R. (2008). “Compaction Effects on Uniformity, Moisture Diffusion, and 

Mechanical Properties of Asphalt Pavements,” Dissertation, Texas A&M University, College 

Station, TX. 

8. Consuegra, A., Little, D.N., Von Quintus, H., and Burati Jr, J.L. (1989.) “Comparative 

Evaluation of Laboratory Compaction Devices Based on Their Ability to Produce Mixtures 

with Engineering Properties Similar to Those Produced in the Field,” Transportation 

Research Record 1228, 80–87.

9. Harvey, J. and Monismith, C.L. (1993). “Effects of Laboratory Asphalt Concrete Specimen 

Preparation Variables on Fatigue and Permanent Deformation Test Results Using Strategic 

Highway Research Program A-003: A Proposed Testing Equipment,” Transportation 

Research Record 1417, 38–48.

10. Peterson, B., Mahboub, K., Anderson, M., Masad, E., and Tashman, L. (2004). “Comparing

Superpave Gyratory Compactor Data to Field Cores,” Journal of Materials in Civil 

Engineering, Vol. 16, No. 1, 78–83. 

11. Masad, E., Muhunthan, B., Shashidhar, N., and Harman, T. (1999). “Internal Structure 

Characterization of Asphalt Concrete Using Image Analysis,” Journal of Computing in Civil 

Engineering, Vol. 13, 88. 

12. Tashman, L., Masad, E., Peterson, B., and Saleh, H. (2001). “Internal Structure Analysis of 

Asphalt Mixes to Improve the Simulation of Superpave Gyratory Compaction to Field 

Conditions,” Journal of the Association of Asphalt Paving Technologists, Vol. 70, 605–645. 

13. Huerne, H.L. (2004). “Compaction of Asphalt Road Pavements,” Dissertation, University of 

Twente, Netherlands. 
98

14. Lee, A.R. and Markwick, A.H.D. (1937). “The Mechanical Properties of Bituminous 

Surfacing Materials Under Constant Stress,” J. Soc. Chem. Ind., Vol. 56, 146–154. 

15. Saal, R. and Labout, J. (1958). “Rheological Properties of Asphalts,” Rheology: Theory and 

Applications, Vol. 2, 363–400. 

16. Saal, R. (1950) “Rheological Properties,” in The Physical Properties of Asphaltic Bitumen:

with Reference to its Technical Applications, J.P. Pfeiffer, ed., Elsevier, Netherlands.

17. Van der Poel, C. (1954). Road Asphalt, Building Materials: Their Elasticity and Inelasticity, 

M. Reiner, ed., North-Holland Publishing, Netherlands. 

18. Reiner, M. (1960). Deformation, Strain and Flow: An Elementary Introduction to Rheology, 

H.K. Lewis & Co., England.

19. Monismith, C.L. and Secor, K.E. (1962). “Viscoelastic Behavior of Asphalt Concrete 

Pavements,” Proceedings of the International Conference on the Structural Design of 

Asphalt Concrete Pavements, 476–498, University of Michigan, Ann Arbor, MI. 

20. Secor, K. and Monismith, C. (1964). “Analysis and Interrelation of Stress-strain-time Data

for Asphalt Concrete,” Journal of Rheology, Vol. 8, 19.

21. Pagen, C.A. (1965). “Rheological Response of Bituminous Concrete,” Highway Research

Record, Vol. 67, 1–26.

22. Monismith, C.L., Alexander, R.L., and Secor, K. E. (1966). “Rheologic Behavior of Asphalt 

Concrete,” Proceedings of the Association of Asphalt Paving Technologists, Vol. 35, 400–450.

23. Huschek, S. (1985). “The Deformation Behavior of Asphalt Concrete Under Triaxial 

Compression,” Proceedings of the Association of Asphalt Paving Technologists, Vol. 54,

407–431.

24. Davis, E.F., Krokosky, E.M., and Tons, E. (1965) “Stress Relaxation of Bituminous Concrete 

in Tension,” Highway Research Record, Vol. 67, 38–58. 

25. Huang, Y.H. (1968). Deformation and Volume Change Characteristics of a Sand-asphalt

Mixture Under Constant Direct and Triaxial Compressive Stresses. Virginia Highway 

Research Council, Charlottesville, VA. 

26. Moavenzadeh, F. and Soussou, J. (1968). “Viscoelastic Constitutive Equation for Sandasphalt Mixture,” Highway Research Record, Vol. 256, 36–52.

27. Perl, M., Uzan, J., and Sides, A. (1983). “Visco-elasto-plastic Constitutive Law for a 

Bituminous Mixture Under Repeated Loading,” Transportation Research Record 911, 

20–27, Transportation Research Board, Washington, DC.

28. Kim, Y.R. (1990). “One-Dimensional Constitutive Modeling of Asphalt Concrete,” Journal 

of Engineering Mechanics, Vol. 116, 751.
99

29. Collop, A.C., Cebon, D., and Hardy, M.S.A. (1995). “Viscoelastic Approach to Rutting in 

Flexible Pavements,” Journal of Transportation Engineering, Vol. 121, No. 1, 82–93.

30. Eisenmann, J. and Hilmer, A. (1987). “Influence of Wheel Load and Inflation Pressure on the 

Rutting Effect at Asphalt Pavements—Experiments and Theoretical Investigations,” 

Proceedings of the Sixth International Conference on the Structural Design of Asphalt 

Pavements, 392–403, University of Michigan, Ann Arbor, MI.

31. Nijboer, L.W. (1948). Plasticity as a Factor in the Design of Dense Bituminous Road 

Carpets, Elsevier Science Ltd., New York, NY. 

32. Van der Poel, C. (1958). “On the Rheology of Concentrated Dispersions,” Rheologica Acta, 

Vol. 1, No. 2, 198–205.

33. Frohlich, H. and Sack, R. (1946). “Theory of the Rheological Properties of Dispersions,”

Proceedings of the Royal Society of London. Series A, Mathematical and Physical Sciences, 

Vol. 185, No. 1003, 415–430. 

34. Hills, J.F. (1973). “The Creep of Asphalt Mixes,” Journal of the Institute of Petroleum, 

Vol. 59, No. 570, 247.

35. Cheung, C.Y., Cocks, A.C.F., and Cebon, D. (1999). “Isolated Contact Model of an Idealized 

Asphalt Mix,” International Journal of Mechanical Sciences, Vol. 41, No. 7, 767–792. 

36. Deshpande, V.S. and Cebon, D. (1999). “Steady-state Constitutive Relationship for Idealized 

Asphalt Mixes,” Mechanics of Materials, Vol. 31, No. 4, 271–287.

37. Deshpande, V.S. and Cebon, D. (2000). “Uniaxial Experiments on Idealized Asphalt Mixes,” 

Journal of Materials in Civil Engineering, Vol. 12, No. 3, 262–271. 

38. Boutin, C. and Auriault, J.L. (1990). “Dynamic Behavior of Porous Media Saturated by a 

Viscoelastic Fluid. Application to Bituminous Concretes,” International Journal of 

Engineering Science, Vol. 28, No. 11, 1157–1181. 

39. Florea, D. (1994). “Associated Elastic/Viscoplastic Model for Bituminous Concrete,” 

International Journal of Engineering Science, Vol. 32, No. 1, 79–86. 

40. Florea, D. (1994). “Nonassociated Elastic/viscoplastic Model for Bituminous Concrete,” 

International Journal of Engineering Science, Vol. 32, No. 1, 87–93.

41. Krishnan, J.M. and Rao, C.L. (2000). “Mechanics of Air Voids Reduction of Asphalt 

Concrete Using Mixture Theory,” International Journal of Engineering Science, Vol. 38, 

No. 12, 1331–1354. 

42. Krishnan, J.M. and Rao, C.L. (2001). “Permeability and Bleeding of Asphalt Concrete Using 

Mixture Theory,” International Journal of Engineering Science, Vol. 39, No. 6, 611–627.

43. Rajagopal, K. (1995). “Multiple Configurations in Continuum Mechanics,” Reports of the

Institute for Computational and Applied Mechanics, University of Pittsburgh, Pittsburgh, PA.
100

44. Rajagopal, K. and Srinivasa, A. (2004). “On the Thermomechanics of Materials That Have 

Multiple Natural Configurations, Part I: Viscoelasticity and Classical Plasticity,” Zeitschrift 

fur Angewandte Mathematik und Physik (ZAMP), Vol. 55, No. 5, 861–893.

45. Rajagopal, K. and Srinivasa, A. (2004). “On the Thermomechanics of Materials That Have 

Multiple Natural Configurations, Part II: Twinning and Solid to Solid Phase Transformation,” 

Zeitschrift fur Angewandte Mathematik und Physik (ZAMP), Vol. 55, No. 6, 1074–1093. 

46. Krishnan, J. and Rajagopal, K. (2005). “On the Mechanical Behavior of Asphalt,” Mechanics 

of Materials, Vol. 37, No. 11, 1085–1100. 

47. Guler, M., Bosscher, P.J., and Plesha, M.E. (2002). “A Porous Elasto-plastic Compaction 

Model for Asphalt Mixtures with Parameter Estimation Algorithm,” Geotechnical Special 

Publication, Vol. 123, 126–143. 

48. Koneru, S., Masad, E., and Rajagopal, K.R. (2008). “A Thermomechanical Framework for 

Modeling the Compaction of Asphalt Mixes,” Mechanics of Materials, Vol. 40, No. 10,

846–864. 

49. Masad, E., Koneru, S., Rajagopal, K.R., Scarpas, A., and Kasbergen, C. (2009). “Modeling 

of Asphalt Mixture Laboratory and Field Compaction Using a Thermodynamics Framework,”

Journal of the Association of Asphalt Paving Technologists, Vol.78, 639–678. 

50. Moore, W. (2006). “Intelligent Compaction: Outsmarting Soil and Asphalt,” Construction 

Equipment, April issue, 38–48.

51. Rahman, F., Hossain, M., and Hunt, M. (2007). “Intelligent Compaction Control of Highway 

Embankment Soil,” Presented at the 86th Annual Transportation Research Board Meeting,

Washington, DC.

52. The MathWorks. (2007). MATLAB® User Guide, Natick, MA. 

53. Scarpas, A. (2000). CAPA-3D Finite Element System Users Manual I, II, and III, Delft 

University of Technology, Netherlands.

54. Scarpas, A. (2005). “A Mechanics Based Computational Platform for Pavement 

Engineering,” Dissertation, Delft University of Technology, Netherlands.

55. Sluys, L. (1992). “Wave Propagation, Localisation and Dispersion in Softening Solids,” 

Dissertation, Delft University of Technology, Netherlands.

56. Bathe, K. (1995). Finite Element Procedures, Prentice Hall, Upper Saddle River, NJ.

57. Johnson, K. (1987). Contact Mechanics, Cambridge University Press, United Kingdom.

HRDi-20/12-10(100)E Recycled

Recyclable