Permanent deformation modelling of large-size unbound pavement materials tested in a heavy vehicle simulator under different moisture conditions

ABSTRACT Climate changes alter the environmental conditions which pavement design is based on, invalidating empirical design methods. Transition to mechanistic design requires the ability to model the behaviour of pavement materials under relevant environmental conditions. An accelerated pavement test (APT) is designed to test two instrumented pavement structures under moisture conditions which are altered by raising the groundwater table (GWT). Open-graded and well-graded subbase materials are used to investigate the effect of gradation on moisture dependency. Pavement response behaviour is modelled using a non-linear elastic (NLE) approach. Accumulation of permanent deformation under different moisture conditions is calculated by two models and compared to measured surface rutting. Moisture transport through the structures differs due to the subbase gradation. Increased GWT accelerates the accumulation of permanent deformations in both structures, identified by both models. One model provides a significantly better fit to the subgrade deformations and the width of the rutting profile.


Introduction
Environmental conditions are important inputs to pavement design, combined with traffic load, material choice and layer thicknesses (Yoder & Witczak, 1975;Zapata et al., 2007). Climate changes are in many areas expected to lead to increased precipitation and more intense rainfall events (Seneviratne et al., 2012). Such changes will lead to increased moisture within road structures and possible overloading of road drainage systems. As the moisture content increases, the friction between aggregate particles becomes lower, and the resistance to differential particle deformation is reduced, leading to a reduced resilient modulus of unbound aggregates (ARA Inc, 2004;Erlingsson, 2010;Lekarp & Dawson, 1998). Pavement saturation during flooding is one of the key deterioration processes that result in degradation of pavement materials (Lu et al., 2020). When unbound materials are saturated, pavement structures observe a significant loss of structural capacity (Elshaer et al., 2019). Pavement design has traditionally been done using empirical methods, based on long-term experience with similar materials and conditions (ARA Inc, 2004). When climate changes cause a shift towards increased moisture levels, empirical methods can no longer be used to predict pavement performance. A transition from empirical to mechanistic design is needed, and to achieve this, we need to be able to model the behaviour of the pavement structures under all conditions (Erlingsson, 2007;Mamlouk, 2006).
In the Nordic countries, most roads are constructed using flexible pavements, with relatively thin hot mix asphalt (HMA) layers above thicker unbound base and subbase providing a substantial part of the bearing capacity. High-quality aggregate resources are readily available, making unbound aggregates an affordable solution. In Norway, local aggregate resources originating from tunnels and road cuts are often utilised in construction projects (Aarstad et al., 2019). Subbase materials with particle sizes > 90 mm are commonly used, even though such materials are outside the scope of European standards and lack suitable quality assurance methods (Fladvad & Ulvik, 2019). Due to experience with frost heave problems caused by excess fines, a practice has developed where all fine material is sorted out from the large-size aggregates used in the subbase layer (Aksnes et al., 2013). When aggregate resources are transported in from quarries, the production methods allow for full quality assessment of the products and large-size subbase materials with a controlled fines content can be used in the subbase layer. The Nordic practice of using large-size unbound pavement materials in pavement structures is uncommon (Fladvad et al., 2017), and has not been a focus in the international development of pavement prediction models. Because of this, there is a need for local calibration of mechanistic design models to Nordic conditions.
In mechanistic design, pavement performance is calculated using numerical methods. Permanent deformation modelling is a central part of pavement performance modelling, and numerous models for unbound granular materials have been developed (e.g. ARA Inc, 2004;Korkiala-Tanttu, 2009;Luo et al., 2017;Rahman & Erlingsson, 2015). Additionally, well-known models e.g. from the mechanisticempirical pavement design guide (MEPDG) are continuously evaluated and improved (e.g. Luo et al., 2017;Schwartz et al., 2011;Zhou et al., 2010). Although open-graded unbound aggregates have been subject for previous research (e.g. Heydinger et al., 1996;Horak & Triebel, 1986;Nguyen & Ahn, 2019), there is a knowledge gap regarding the performance of large-size open-graded aggregates used in pavements.
An accelerated pavement test (APT) using a heavy vehicle simulator (HVS) was set up to enable a full-scale investigation of the differences between two types of large-size subbase materials. Two subbase materials of equal geological composition were tested, one material containing fines, while the other material had lower sieve size 22 mm. Albeit being composed of the same rock type, such materials will have different properties regarding density, void ratio, permeability and stiffness. Depending on gradation, varying groundwater levels result in varying moisture content in the pavement materials also above the groundwater table (GWT); hence, the influence of moisture content on pavement performance can be evaluated (Erlingsson, 2010;Li & Baus, 2005;Saevarsdottir & Erlingsson, 2013b).
The response behaviour of the pavement structures was modelled by Fladvad and Erlingsson (2021), and although differences were found in the stress and strain response in the two structures, both structures showed similar rut development throughout the test. Additionally, the development of permanent deformations in both structures measured as surface rut development was surprisingly linear. The development of permanent deformations is investigated further in the current research.
The research presented in this paper has three specific objectives: • Investigate the development of permanent deformations in pavement structures with opengraded and well-graded subbase materials • Investigate the influence of increased moisture levels following raised GWT on the open-graded and well-graded subbase structures • Investigate the feasibility of two models for calculating the accumulation of permanent deformations in large-size unbound aggregates The APT was conducted at the Swedish National Road and Transport Research Institute (VTI) test facility in Linköping, Sweden in 2018-2019.

Pavement structures
The APT was conducted using a HVS type Mark IV (Figure 1(a)) at the Swedish National Road and Transport Research Institute (VTI) full-scale pavement testing facility in Linköping, Sweden. The pavement structures were constructed in a concrete test pit which is 3 m deep, 5 m wide, and 15 m long (Figure 1(b)). In the test pit, GWT can be controlled and adjusted.
Two pavement structures were tested simultaneously in order to investigate the difference between an open-graded and a well-graded subbase. Both structures were constructed in the same test pit, each part with a length of 7.5 m. The pavement structures were constructed of identical materials, except for the gradation of the subbase materials: 0/90 and 22/90 mm. The surface course and bituminous base course were constructed from conventional dense-graded asphalt concrete (AC), where the surface course had upper aggregate size 16 mm and pen 70/100 binder. The bituminous base course had upper aggregate size 22 mm and pen 160/220 binder. Unbound base course (0/32 mm) and subbase were constructed from crushed rock. Figure 2 shows the particle size distribution curves for the unbound materials. The subgrade consisted of silty sand. Figure 3 shows the layer thicknesses for the full pavement structures. The difference in layer thickness between the two structures is due to practical adjustments in the construction process, such as differences in compaction. Due to the geometry of the test pit, the unbound materials were installed partly manually as the test pit could not be accessed by full-size construction equipment normally used to install large-size materials. The unbound materials were compacted until they met requirements for static plate load tests given by Norwegian Public Roads Administration (2018), thus ensuring that   design stiffness of the materials was reached. The AC was placed in the test pit using conventional paving equipment.
The subbase materials were supplied from a tunnel construction site in the construction project E39 Svegatjørn-Rådal south of Bergen, Norway. The materials fulfilled requirements for subbase materials (Norwegian Public Roads Administration, 2018) regarding gradation, fines content and physical properties. The open-graded material was selected from the ongoing production at the construction site, while the well-graded material was created by mixing 22/90 and 0/32 mm material at a 1:1 weight ratio.

Heavy vehicle simulator
In the main test, the HVS was set up with a dual-wheel configuration. Two 295/80R22.5 tyres were used with a 34 cm centre-to-centre wheel spacing. The dual-wheel load during testing was 60 kN, corresponding to a 120 kN axle load, and the tyre pressure was 800 kPa. The lateral wander of the wheel was ±25 cm in 5 cm increments, following the normal distribution shown in Figure 4.
Before the main test started, a pre-loading phase with lower load was applied in order to postcompact the structure. For the initial 20,000 load repetitions, a single 425/65R22.5 wheel was used  with a load of 30 kN (axle load 60 kN) and tyre pressure 700 kPa. In this phase, the load was evenly distributed over the wheel path, ±35 cm in 5 cm increments. All traffic load was applied with bidirectional loading at a constant rolling speed of 12 km/h. A climate chamber was installed to keep the temperature constant at 10 • C.

Groundwater table
The APT was divided into three phases, each with a separate location of the GWT. In phase w1, GWT was located at great depth, > 3 m below the pavement surface. GWT was stable at this level for the first 550,000 load repetitions. For phase w2, GWT was raised to 30 cm below the formation level (Figure 3), corresponding to the depth of the drainage level of a pavement structure in operation. GWT was stable at this level for 368,000 load repetitions. For phase w3, GWT was raised further to a level of about 5 cm above the formation level. This level simulates a situation where the drainage system is overloaded and unable to keep the GWT at the designed drainage level. GWT was stable at this level for 286,000 load repetitions. No traffic load was applied while GWT was raised. Opening a road before the floodwater is drained out to an appropriate extent before saturation determines the following pavement deterioration trend after a flood hazard (Lu et al., 2020). Phase w3 is an approximation of this kind of situation.

Instrumentation
Both pavement structures were instrumented using asphalt strain gauges (ASG), soil pressure cells (SPC), strain measuring units (εMU) and moisture and temperature sensors ( Figure 3). ASGs placed at the bottom of the AC measured longitudinal and transversal strain, while SPCs at three levels in the subbase layers measured vertical stress. εMUs measured vertical strain at seven depth intervals; over the bound base layer, unbound base layer, subbase divided between three sensors and two sensors in the subgrade, the lowest reaching to approximately 30 cm below the formation level. Each structure was instrumented with two sets of εMUs and SPCs, and three sets of ASGs. Moisture content was measured at four different levels; in the subgrade approximately 15 cm below the formation level, at two levels in the subbase layer, and in the middle of the base layer. Temperature sensors were installed in the bituminous base 4 and 8 cm below the pavement surface.The pavement response was measured regularly during the test, and permanent deformations were measured 13 times, as shown in Table 1. Moisture content was registered continuously.

Surface measurements
The surface rut profile was measured using a laser which registers the surface cross-section by 250 measurements over 2.5 m width. Three profiles were registered for each structure. Laser profiles were measured on average once per 24,000 load repetitions. In the presentation of rut depth values and profiles, the measured rut depth is corrected for the surface variations measured before any traffic load was applied. Average rut depth is calculated as the average over the middle 200 mm of the rut profile (average of 20 individual measurements).

Models
The modelling of pavement performance is conducted in two steps. First, a response model is established, describing the material parameters of each layer. The response model is calibrated according to measured vertical strain and stress throughout the structure. Next, the vertical strain as a function of depth calculated from the chosen response model is implemented in two permanent deformation models. Both phases of modelling take the changing moisture levels due to the raised GWT into account. The theoretical background for the modelling is described in the following.

Stiffness of unbound materials
The resilient modulus M r is a common characteristic of the stiffness of unbound materials, as defined in Equation (1), where ε r is the resilient (recoverable) strain in the material under the deviatoric stress σ d (Equation (2)). A stiffer material shows less strain under a certain load and thus has a higher M r .
where θ bulk stress; θ = σ 1 +σ 2 +σ 3 3 ; p a reference pressure (100 kPa) k 1 ; k 2 experimentally determined constants In the research presented in this paper, the pavement structures' response behaviour is modelled using multi-layer elastic theory (MLET) in an axisymmetrical system using the ERAPave software (Ahmed & Erlingsson, 2012;Erlingsson & Ahmed, 2013) to find the appropriate material properties for the unbound pavement materials. ERAPave is used to calculate both linear elastic and non-linear elastic material behaviour.

Moisture dependency
Increasing moisture content results in a decreasing M r (Lekarp et al., 2000). The variation in M r depending on moisture can be expressed by Equation (4) given in the mechanistic-empirical pavement design guide (MEPDG) (ARA Inc, 2004;Cary & Zapata, 2011). Degree of saturation S can be calculated from the water content and porosity of the material.
where S opt degree of saturation at a reference condition M r opt resilient modulus at a reference condition a minimum of log 10 ( M r The calibrated MEPDG M r -moisture model shows that coarse-grained materials are less affected by a change in S than fine-grained materials. Due to the large upper size of the subbase materials used in the current research, optimum water content or saturation could not be evaluated, and established relations regarding variations in M r as a function of saturation (Equation (4)) could not be used directly. Due to these constraints, M r was back-calculated directly from measurements.

Permanent deformation models
Two models are employed to calculate the development of permanent deformations throughout the APT. Both models calculate the permanent deformation from the number of load repetitions (N) and resilient strain (ε r ). In the calculations, the pavement layers can be divided into sublayers. Both models calculate accumulated permanent strainε p , which can be used to calculate the total accumulated permanent deformationδ p of the structure using Equation (5).
In Equation (5),ε p jk is the average permanent strain of sublayer j of layer k, and h jk is the thickness of sublayer j of layer k. Furthermore, n is the total number of layers, and m is the total number of sublayers within each layer.
In HVS testing similar to the current, Saevarsdottir and Erlingsson (2013b) found the deformation of AC layers to be less than 1 mm, and deemed it negligible. Ahmed and Erlingsson (2015) tested AC materials with similar binder properties as the current test at different temperatures in an extra-large wheel-tracking test and a heavy vehicle simulator, and found that at 10 • C, hardly any permanent deformation is accumulated. As the temperature was held constant at 10 • C throughout the current test, the total permanent deformation of the AC layers is expected to be less than 1 mm. From this, the calculation of permanent modelling in the following is limited to unbound materials. Tseng and Lytton (1989) presented a method characterising the accumulated permanent strainε p of pavement materials in terms of the material parameters ε 0 ε r , β and ρ. The Tseng & Lytton model was further developed and presented in MEPDG (ARA Inc, 2004):

MEPDG model
The MEPDG model is calibrated for unbound pavement materials using repeated load triaxial (RLT) tests. In Equation (6), ε 0 ε r , ρ and β can be calculated from the water content using Equation (7)-(10) (ARA Inc, 2004). β 1 is a calibration factor for each layer. The average vertical resilient strain in the midpoint of the layer (ε r ) is denoted ε v to avoid confusion with the ε 0 ε r factor.
Gravimetric water content w c is calculated from volumetric water content w vol using Equation (8), where γ water and γ is the densities of water and soil. Figure 5. The time-hardening approach for calculating accumulated permanent strain from a sequence of stress levels (Erlingsson & Rahman, 2013).
In Equation (9) and (10), a 1 = 0.15 and a 9 = 20.0. A time-hardening approach can be used to enable the calculation of permanent deformation from a sequence of stress levels such as in the multi-stage repeated load triaxial (MS RLT) tests (Erlingsson & Rahman, 2013;Lytton et al., 1993;Zhou et al., 2010). The time-hardening concept takes stress history from previous loading into account when calculating the permanent deformation from the present loading ( Figure 5). An equivalent number of load cycles N eq is added to the number of load cycles at each load level, corresponding to the already accumulated permanent strain at the start of each level. The illustration in Figure 5 assumes an increasing stress level from one stress path to the next. If the stress level in one cycle is significantly lower than the previous cycle, N eq approaches infinity, and no further permanent strain development is calculated from the current cycle.
For the MEPDG model, N eq is calculated from Equation (11), where i refers to the current stress path, and i−1 refers to the previous stress path. Rahman and Erlingsson (2015) introduced a stress-based empirical model for predicting the permanent deformation behaviour of unbound materials (Equation (12)) based on the power law presented by Sweere (1990). The model is based on MS RLT testing.

RE model
The term S f in Equation (12) is a descriptor of the stress conditions, defined in Equation (13).
The factors a, b and α are regression-based material parameters. The term p a is a reference stress assumed equal to the atmospheric pressure (p a = 100 kPa), applied to make the S f expression dimensionless. Rahman and Erlingsson (2016) proposed that the factor a is moisture dependent and increases linearly with increasing water content, as described by Equation (14).
The S f model is well suited to RLT tests, but gives wrongful results when applied for cases where negative values for the principal stresses occur from the response model. Hence, the model in Equation (12) has been further developed to the model in Equation (15) (Rahman & Erlingsson, 2021), where the stress factor S f is exchanged for resilient strain ε r .
The coefficients a and b are material parameters depending on moisture content, degree of compaction and particle size distribution. The time hardening approach in Figure 5 can also be used for this model, employing Equation (16).
The model described by Equations (15) and (16) will in the following be referred to as the RE model.

Incorporation of lateral wander and cross-section width
For both models, the permanent deformations are calculated in a width of 1 m from the centre line in 5 cm increments. The resilient strain is calculated for all wheel positions relative to each cross-section position. The calculations encompass the lateral wander of the loading wheel by using the time hardening procedure (Equations (11) and (16)) and the resilient strain calculated for each wheel position. The wander is incorporated stepwise using the normal distribution ( Figure 4) in steps of 500 load repetitions where all positions from −25 to +25 cm is covered, with 120 load repetitions in the centre position reducing down to two load repetitions in position ±25 cm. The time hardening procedure is also used to calculate the transition between phases, where both resilient strain and calculation coefficients (β, ρ, ε 0 ε r , β 1 , a) changes. Only resilient strain changes between side positions, the full cross-section is calculated using the same coefficients as the centre position.

Water content
The raised GWT resulted in an increased volumetric water content w vol in the unbound materials, as shown in Figure 6. In the unbound base, w vol is constant at about 7.3% for both structures in phase w1. In phase w2, w vol remains constant in the 22/90 mm structure, while it increases to 8.6% in the 0/90 mm structure. This increase finishes about one month after GWT is stable at the new level. In phase w3, both structures reach a final w vol of 10.0%. Again, this increase takes about one month in the 0/90 mm structure. In the 22/90 mm structure, on the other hand, the increase from about 7.3% to 9.9% takes place over about 100 days. The 22/90 mm crushed rock subbase shows a limited capacity of storing moisture, as w vol starts at 1.5-2.0% in phase w1 and only increases to 2.0-2.5% at the end of phase w3. For this material, the highest w vol is registered by the sensor located at the higher level in the structure. For the 0/90 mm crushed rock material, the tendency is opposite, the lower sensor registers w vol from 8.8% in phase w1 via 9.1% (w2) to 10.3% in phase w3, while the upper sensor spans from 6.1 to 6.8%.
In the subgrade, the registrations show a decreasing tendency in phase w1, likely caused by drainage of water added during compaction of the structures. From the end of phase w1 to the beginning of phase w2 the increase in w vol is very rapid from 11.4 to 26.0%. The difference between phase w2 and w3 is about 1%, showing that the subgrade is almost fully saturated already in phase w2, even though GWT is located about 15 cm below the moisture sensor at this time.
The average w vol during traffic loading is displayed in Table 2, where each phase corresponds to the following time periods in Figure 6: w1 -days 26-89; w2 -days 105-131; w3 -days 209-303. Overall, the moisture sensor registrations show that the water content is affected by the raised GWT far above the actual GWT. The open-graded subbase material does transport moisture vertically, although it takes much longer time than in the well-graded subbase. In the end, w vol in unbound base stabilises at the same level for both structures, indicating that w vol has reached equilibrium.

Measured surface rutting
The development of the rut profile throughout the APT is shown in Figure 7, where the laser-measured surface profile is displayed for approximately every 100,000 load repetitions. The displayed profiles represent the average of three laser profiles for each structure. The applied lateral wander distribution   of the loading is given in Figure 4. There are no signs of uplift in the outer parts of the profile as the rut depth increases. The rutting mainly takes place within 1 m from the centre line. The development of average rut depth for both structures throughout the test is shown in Figure 8. The rut depth is calculated as the average depth of the mid 200 mm of three laser profiles for each structure, with error bars showing the variation between the three profiles. The 22/90 mm structure show a higher rut development in phase w1, but is less affected by the increased GWT between phase w1 and w2. At the end of phase w2, both structures have a similar rut depth. Both structures show a similar rut depth in phase w3. As the error bars are overlapping, no real distinguishment can be made between the structures in phase w3.
A linear trend for the rut development is established after 50,000 load repetitions in phase w1, after some initial post-compaction of the structure. The rut development continues a linear trend  throughout phase w2 and w3, but with increased slope after each phase transition. In phase w2, the rut development is clearly higher for the 0/90 mm structure. Figure 9 shows that a linear trend for the rut development is established after 50,000 load repetitions in phase w1, after some initial post-compaction of the structure. The rut development continues a linear trend throughout phase w2 and w3, but with increased slope after each phase transition. In phase w2, the rut development is clearly higher for the 0/90 mm structure. Table 3 shows the rutting rate calculated as the average rut increase per 100,000 load repetitions for the final part of each phase, when the development is linear. The rutting rate corresponds to the slope of the lines in Figure 8, and shows how the increased GWT causes a substantial acceleration of rut development between each phase. Both structures show a similar rutting rate in phase w1. From w1 to w2, the increase is significantly larger for the 0/90 mm structure compared to the 22/90 mm structure. This may indicate that the well-graded subbase structure is more moisture susceptible, but in the transition to phase w3, the trend is opposite, and the 22/90 mm structure ends up at a slightly higher rutting rate in the final phase. Figure 10 shows the permanent deformations calculated from static εMU registrations. The registrations are from one sensor in the unbound base, three sensors in the subbase and one sensor representing the top 15 cm of the subgrade. The pavement structures were instrumented with εMU sensors for the top 30 cm of the subgrade, but the static registrations from the lowest span proved invalid after the GWT was raised in phase w2. Hence, the Sg values presented in Figure 10 represent only the top 15 cm of the 2.4 m thick subgrade.

Permanent deformation in unbound layers
For both structures, the unbound base shows no reaction to the increased GWT between phase w1 and w2, while a clear acceleration of the rut development is seen in phase w3. The 22/90 mm subbase show a large increase in permanent deformation throughout all phases, whereas the 0/90 mm subbase show a clear acceleration between phase w2 and w3. These measurements contrast the tendency suggested by the water content measurements, where the change in moisture content for subbases are small, while larger changes are seen for the unbound base. However, as Equation (4) suggests, the dependency of moisture is not calculated from moisture content directly, but depends on the difference between observed and optimum moisture content. The effect of changed moisture on material stiffness is larger close to optimum moisture content.

Modelling results
The MEPDG and RE models both calculate permanent deformations based on induced resilient strain throughout the structures, called ε v in the MEPDG model and ε r in the RE model. The response model must, therefore, provide a reliable estimation of vertical strain as a function of depth.

Response modelling
The pavement response in the APT has previously been modelled by Fladvad and Erlingsson (2021), where two different approaches were applied to model the responses; a linear elastic (LE) modelling approach using material stiffnesses back-calculated from falling weight deflectometer (FWD) measurements, and a non-linear elastic (NLE) modelling approach based on Equation (3) fitted to the measured strain levels. In the NLE model, AC and subgrade are treated as linear elastic materials (Equation (1)), while the unbound base and subbase materials are treated as non-linear elastic (Equation (3)).
As the proposed models for permanent deformations both use resilient strain from the response model, the most important result is the response models' fit to measured induced vertical strain, shown as a function of depth in Figure 11.
The results showed that the LE approach based on FWD was not able to replicate the induced vertical strain measured by εMUs in the pavement structures. The NLE model with parameters shown in Table 4 fitted the measured vertical strain very well, and was selected as the primary response model in the further modelling of permanent deformations in the structures. In the NLE calculations, k 2 and Poisson's ratio is assumed constant for all materials at 0.6 and 0.35, respectively. The AC stiffness value in Table 4 is reduced due to the low speed of the HVS (12 km/h), and corresponds to AC stiffness of 6500 MPa under normal traffic speed at 10 • C (Kim, 2011). Figure 11 shows the resilient strain down to 1 m below the surface, but the model was used to calculate strain through the uninstrumented subgrade down to the bottom of the test pit, at 3 m depth.  For both structures, the response model supports the findings from Figure 6, that the subgrade is almost fully saturated already in phase w2, as the increase in strain between phase w2 and w3 is very small. For the 0/90 mm structure, the response model shows a gradual increase in strain in all unbound layers as GWT is raised. For the 22/90 mm structure, strain in the subgrade and two lower layers of the subbase increases as GWT is raised. In the unbound base and upper part of the subbase, on the other hand, strain decreases between phase w2 and w3, to a level lower than in phase w1. In the unbound base, there is no difference in strain between phase w1 and w2, corresponding well to the moisture content, which does not change between these phases.

MEPDG model
The material parameters used to calculate permanent deformations using the MEPDG model is shown in Table 5. The structures are divided into sublayers following Table 4, except for the uninstrumented subgrade which is divided into four sublayers of 53 cm.
(2) Calculate ε v for each GWT phase from response model at the midpoint of each sublayer.
(9) Verify model by comparing sum ofδ p,calc from all layers (Equation (5)) to measured rut depth from surface laser measurements. If poor fit, adjust calculation ofδ p,calc from uninstrumented subgrade by adjusting β 1 .
Initially, the calculation procedure aimed to calculate one average β 1 for each material covering all GWT phases, as the calculation in step 1 and 2 should account for the moisture variations. However, these calculations provided a poor fit to the measuredε p for all materials, and showed that a separate calculation of β 1 for each GWT phase was necessary.
The subgrade is constructed of the same material for both structures, hence, the same material parameters are used for both structures. In phase w2, subgrade 3-6 is below the GWT and is assumed to have the same water content as the subgrade moisture sensor registers in phase w3 when the sensor is below GWT. Following step 9, a gradual reduction of β 1 with depth was calculated for the uninstrumented subgrade.

RE model
The calculations of permanent deformations using the RE model following Equation (15) are defined by the a and b values given in Table 6. The structures are divided into sublayers following Table 4, except for the uninstrumented subgrade which is divided into four sublayers of 53 cm.
Factor a is calculated using Equation (14), using linear regression to the factors found by fitting a for each sublayer to the permanent strain found from static εMU measurements. The modelling aimed to establish a a(w) relationship per material. Calculation procedure: (1) Calculate ε r for each GWT phase from response model at the midpoint of each sublayer.
(2) Calculate accumulated permanent strainε p,calc from elastic strain using Equations (15) and (16) using a set b value. (3) Fitε p,calc toε p,meas from εMU for each phase by adjusting a.
(7) Verify model by comparingδ p,calc toδ p,meas from εMU for each instrumented layer. If poor fit, adjust b and repeat procedure from step 2. (8) Verify model by comparing sum ofδ p,calc from all layers (Equation (5)) to measured rut depth from surface laser measurements. If poor fit, adjust calculation ofδ p,calc from uninstrumented subgrade.
The starting point in step 2 was b = 250, based on recommendations by Rahman and Erlingsson (2021). The subgrade is constructed of the same material for both structures, and is modelled with one set of a(w) and b parameters. For the 0/90 mm subbase, b had to be reduced in order to achieve correlation between a and water content. For the 22/90 mm subbase, no correlation between a and water content was found, regardless of b. In the analysis, b was varied from 50 to 400. Hence, b for the 22/90 mm subbase was set to the same value as for the 0/90 mm subbase, and the a values fitted to the permanent strain in each sublayer was used in the modelling of permanent deformations.
In phase w2, subgrade 3-6 is below the GWT and is assumed to have the same water content as the subgrade moisture sensor registers in phase w3 when the sensor is below GWT. The model provided a reasonable fit to the rut depth after step 7, and no adjustments for the uninstrumented subgrade was necessary in step 8. Figure 12 shows the calculated permanent deformations from both models compared to measured permanent deformations for each instrumented layer in each structure. The black line in each figure represents the sum of the UB, Sb and Sg lines. For the subbase, the presented values represent the sum of three sublayers.

Permanent deformations in individual layers
Both models provide a reasonable fit to the measured permanent deformations. Note that the good fit of the 22/90 mm subbase from the RE model is somewhat misleading, as no a(w) relationship was found for this material, and the deformation development is fitted directly to the deformation measured by εMUs. Both models calculate a deformation development with a fast increase at the start of each phase which slows down as more load repetitions are added.
The difference between the total modelled deformation in Figure 12 and the total modelled rut depth in Figure 13, is the calculated permanent deformations for the uninstrumented subgrade (0.73-3.0 m). It is clearly visible that the MEPDG calculates more deformation in the uninstrumented subgrade than the RE model, especially in phase w1 and w3. The results are similar for both structures, although the values are slightly higher for the 0/90 mm structure in phase w2 and w3, due to higher resilient strain calculated from the response model. The transition between w2 and w3 should not make a large impact on the uninstrumented layers, as most of the subgrade was below GWT already in phase w2. The resilient strain is however affected, as the load distribution into the subgrade is affected by the moisture and stiffness of the layers above, which change significantly between w2 and w3 (Figures 6 and 11). The results show that the MEPDG model calculates a bigger difference between phase w2 and w3 than the RE model, which cannot be explained by moisture dependency, as the moisture level is constant. The results thus show that the two models judge the dependency between resilient and permanent strain differently. Figure 13 compares the modelled permanent deformations from both models and structures to the measured rut depth from the pavement surface. The modelled data in Figure 13 is the sum of the permanent deformations in the instrumented and uninstrumented parts of the structures. Neither model is able to replicate the linear rut development, especially for phase w2. The RE model is close to a linear development in phase w3. Both models overestimate the rut development in phase w2.    Figure 14 shows the cross-sections of the modelled permanent deformations to the measured laser rut profiles at the end of each GW phase. The overestimation of rut depth, especially from the MEPDG model in Figure 13 is clearly visible also here. The RE model provides the better fit to the overall shape of the rut profile, calculating a narrower rut profile than the MEPDG model. Figure 13 shows that the best fit is found close to the phase transitions, so in that sense, Figure 14 shows the profiles when they are closest to the measured profile.

Calibration of the MEPDG model
In the definitions of the MEPDG model, ARA Inc (2004) presents global calibration factors for unbound granular base (β GB = 1.673) and subgrade soils (β SG = 1.35). Several authors have found that the MEPDG model shows low sensitivity to inputs from unbound layers (e.g. Luo et al., 2017;Schwartz et al., 2011). Gu et al. (2016) used the global calibration, and found that the MEPDG model generally underestimated the permanent deformation behaviour of unbound materials. In the present research as well as similar previous research (Saevarsdottir & Erlingsson, 2013a, these global factors have been substituted for β 1 factors fitted to measurements. In the present research, β 1 values for unbound base and subbase are generally lower than the suggested calibration factor β GB , varying from 0.25 to 1.37. Using the original factors would heavily overestimate the permanent deformations in unbound layers, contrary to Gu et al. (2016). In the subgrade, β 1 values vary from 0.5 to 5.0, also showing that using a single β SG would not provide a good fit to the measured permanent deformation. The current research does not encompass a sensitivity analysis, and can not be used to evaluate the models' sensitivity to parameter changes. However, we find that with sufficient measurement data to verify the models against, a reasonable fit to full-scale traffic deterioration can be found ( Figure 12).
A linear calibration factor such as β 1 should also be implementable for the RE model, scaling the accumulation of permanent deformation in addition to the a(w) relationship. However, such implementation would require laboratory measurements and more material data than was available in this research. For the present research, the direct relationship between a and water content was sufficient to provide a reasonable fit to the accelerated accumulation of permanent deformation as GWT was raised.

Fit to linear rut development
Neither model is able to replicate the linear rut development seen from the surface laser measurements, although the RE model is close in phase w3. During construction, special emphasis was laid on fulfilling requirements for compaction measured by the static plate load test even though the pavement structures were not constructed using equipment corresponding to a regular road construction project. The stiffness requirements were met by adding extra compaction effort. As a result, the unbound materials in the APT are likely more well-compacted than the triaxial test samples the models are based on.
In the calculations, the average moisture content for each phase is used (Table 2), assuming constant moisture conditions in each phase. We know from Figure 6 that the moisture content varied (increased) during the periods of accelerated traffic. The assumption of constant water content may have affected the calculations for the unbound base, where the change in moisture was greatest. For the subbase and subgrade, the change in moisture content during traffic load was very small, and would likely not have impact on the calculations. The 22/90 mm subbase had a very low and nearly constant moisture content throughout the test.
In order to change the slope of the MEPDG prediction, factors β and ρ must be changed. In the current research, these factors were calculated directly from moisture content using Equations (7) and (9). The slope differences between the measured and modelled rut development suggest that these calculations should be calibrated as well. Zapata et al. (2007) found that the linear β(w)-relationship suggested in MEPDG (Equation (7)) was insufficient, and should be exchanged for a exponential relationship. Such changes would, in turn, change the β 1 calibration discussed previously. To change the slope of the RE model, factor b must change. In this research, we chose to keep b as constant as possible to isolate the impact of moisture on factor a. In further analysis of the data, the impact of varying b should be investigated. Rahman and Erlingsson (2016) identified a linear relationship between a and water content for the model described by Equation (12), and a similar relation is expected for the RE model. As strain increases when moisture content increases, a should have a positive correlation to water content. However, for the unbound base, the current data provide a positive correlation for the 22/90 mm structure and negative correlation for the 0/90 mm structure. As the unbound base is constructed from the same material for both structures, the relation should be of the same type for both structures even though the water content is different. The response model for the 22/90 mm structure introduces a complication for the unbound base calculation, as the induced vertical strain is lower in phase w3 than in phases w1 and w2, even though the water content has increased. The conducted analyses have not uncovered the reason for why the unbound base material cannot be calculated as a single material in the RE model, while it works well for the MEPDG model.

Models' fit to open-graded large-size material
For the 22/90 mm subbase, the RE model a values fitted to permanent strain measurements varied by factor 50 from 0.2 to 11.2, while the water content varied between 1.5 and 2.5%. The open-graded subbase material is by design a draining material and thus not moisture susceptible, which may explain why no correlation is found between a and water content. For sublayer 1 and 3 for the subbase, a increases with increasing water content. In sublayer 2, a increases from 6.8 to 11.2 from phase w1 to w2, before it decreases to 5.7 for phase w3 (Table 6). Both the relation between phases and the values are far from all other layers, which may indicate some kind of measurement error in the permanent deformations which are used to verify the model. However, when the permanent deformations are calculated with these high values, the resulting deformation of the full subbase layer is lower than calculated with the MEPDG method.

Contributions from uninstrumented subgrade
The calculations have shown that the response model and water content is not sufficient basis for calculating permanent deformations using the MEPDG models, as the variation in β when the model is verified against measured deformation is considerable. These findings suggest that the ε 0 ε r relationship is more affected by the GWT changes than the model predicts.
Both models provide a reasonable fit to the measured permanent deformations for both structures in Figure 12, but end up overestimating the rut depth in Figure 13. As the difference between the total rut depth in Figures 12 and 13 is the contribution from the uninstrumented subgrade, the subgrade rutting in the final phase must be overestimated. The RE model shows less overestimation, suggesting a more accurate subgrade calculation. This is also supported by the fact that no adjustments to the modelling coefficients were necessary when calculating the deformation using this model. Without reduction of the β 1 values, the MEPDG model calculated deformation values of the same magnitude as the surface rut from the uninstrumented subgrade alone, which undoubtedly is wrong.

Conclusions
The development of permanent deformation during accelerated traffic loading was investigated for two pavement structures with different gradations of the subbase layers. During the APT, GWT was raised twice, causing moisture variations in the pavement structures. The permanent deformation behaviour was modelled using two strain-based models. The following conclusions can be drawn from the analyses: • The raised GWT resulted in a clear acceleration of the permanent deformations in both structures.
• The NLE response model captured the resilient strain behaviour of the structures very well, and formed a reliable basis for the permanent deformation models. • The time-hardening approach succeeded in implementing both lateral wander and changing moisture conditions. • Both permanent deformation models identified the acceleration of permanent deformations as GWT was raised. • Neither model replicated the linear rut development measured on the surface. • The RE model provided a better fit to the width of the rutting profile, showing a narrower profile than the MEPDG model. • The RE model gave a better estimation of the subgrade deformations, where the MEPDG model required significant adjustments in order to reduce the modelled deformation. • The RE model generally provided a simpler moisture dependency calculation compared to the MEPDG model. For the MEPDG model, the suggested calculation of moisture dependency was not sufficient to cover the difference between GWT phases. However, no moisture dependency was found for the open-graded subbase material in the RE model. • The subbase gradation did not make a great impact on the total rut depth, but these results were likely governed by the high degree of compaction of the unbound materials achieved during construction. • The subgrade deformation constituted a significant part of the total rut depth, showing that knowledge of the subgrade soil is essential to achieve a reliable modelling result.
The current research clearly showed how increased moisture in pavement structures impact the pavement performance. The new RE model was very well suited for cross-sectional calculation of rut development from full-scale testing, and easily adaptable to changing moisture conditions. For wellcompacted large-size pavement materials, gradation was of less importance for the overall pavement performance, as both structures showed similar rut development throughout the test.