**Suggested Citation:**"CHAPTER 4. FINDINGS." National Academies of Sciences, Engineering, and Medicine. 2019.

*Proposed Enhancements to Pavement ME Design: Improved Consideration of the Influence of Subgrade and Unbound Layers on Pavement Performance*. Washington, DC: The National Academies Press. doi: 10.17226/25583.

**Suggested Citation:**"CHAPTER 4. FINDINGS." National Academies of Sciences, Engineering, and Medicine. 2019.

*Proposed Enhancements to Pavement ME Design: Improved Consideration of the Influence of Subgrade and Unbound Layers on Pavement Performance*. Washington, DC: The National Academies Press. doi: 10.17226/25583.

**Suggested Citation:**"CHAPTER 4. FINDINGS." National Academies of Sciences, Engineering, and Medicine. 2019.

*Proposed Enhancements to Pavement ME Design: Improved Consideration of the Influence of Subgrade and Unbound Layers on Pavement Performance*. Washington, DC: The National Academies Press. doi: 10.17226/25583.

**Suggested Citation:**"CHAPTER 4. FINDINGS." National Academies of Sciences, Engineering, and Medicine. 2019.

**Suggested Citation:**"CHAPTER 4. FINDINGS." National Academies of Sciences, Engineering, and Medicine. 2019.

**Suggested Citation:**"CHAPTER 4. FINDINGS." National Academies of Sciences, Engineering, and Medicine. 2019.

**Suggested Citation:**"CHAPTER 4. FINDINGS." National Academies of Sciences, Engineering, and Medicine. 2019.

**Suggested Citation:**"CHAPTER 4. FINDINGS." National Academies of Sciences, Engineering, and Medicine. 2019.

**Suggested Citation:**"CHAPTER 4. FINDINGS." National Academies of Sciences, Engineering, and Medicine. 2019.

**Suggested Citation:**"CHAPTER 4. FINDINGS." National Academies of Sciences, Engineering, and Medicine. 2019.

**Suggested Citation:**"CHAPTER 4. FINDINGS." National Academies of Sciences, Engineering, and Medicine. 2019.

**Suggested Citation:**"CHAPTER 4. FINDINGS." National Academies of Sciences, Engineering, and Medicine. 2019.

**Suggested Citation:**"CHAPTER 4. FINDINGS." National Academies of Sciences, Engineering, and Medicine. 2019.

**Suggested Citation:**"CHAPTER 4. FINDINGS." National Academies of Sciences, Engineering, and Medicine. 2019.

**Suggested Citation:**"CHAPTER 4. FINDINGS." National Academies of Sciences, Engineering, and Medicine. 2019.

**Suggested Citation:**"CHAPTER 4. FINDINGS." National Academies of Sciences, Engineering, and Medicine. 2019.

**Suggested Citation:**"CHAPTER 4. FINDINGS." National Academies of Sciences, Engineering, and Medicine. 2019.

**Suggested Citation:**"CHAPTER 4. FINDINGS." National Academies of Sciences, Engineering, and Medicine. 2019.

**Suggested Citation:**"CHAPTER 4. FINDINGS." National Academies of Sciences, Engineering, and Medicine. 2019.

**Suggested Citation:**"CHAPTER 4. FINDINGS." National Academies of Sciences, Engineering, and Medicine. 2019.

**Suggested Citation:**"CHAPTER 4. FINDINGS." National Academies of Sciences, Engineering, and Medicine. 2019.

**Suggested Citation:**"CHAPTER 4. FINDINGS." National Academies of Sciences, Engineering, and Medicine. 2019.

**Suggested Citation:**"CHAPTER 4. FINDINGS." National Academies of Sciences, Engineering, and Medicine. 2019.

**Suggested Citation:**"CHAPTER 4. FINDINGS." National Academies of Sciences, Engineering, and Medicine. 2019.

**Suggested Citation:**"CHAPTER 4. FINDINGS." National Academies of Sciences, Engineering, and Medicine. 2019.

**Suggested Citation:**"CHAPTER 4. FINDINGS." National Academies of Sciences, Engineering, and Medicine. 2019.

**Suggested Citation:**"CHAPTER 4. FINDINGS." National Academies of Sciences, Engineering, and Medicine. 2019.

**Suggested Citation:**"CHAPTER 4. FINDINGS." National Academies of Sciences, Engineering, and Medicine. 2019.

**Suggested Citation:**"CHAPTER 4. FINDINGS." National Academies of Sciences, Engineering, and Medicine. 2019.

**Suggested Citation:**"CHAPTER 4. FINDINGS." National Academies of Sciences, Engineering, and Medicine. 2019.

**Suggested Citation:**"CHAPTER 4. FINDINGS." National Academies of Sciences, Engineering, and Medicine. 2019.

**Suggested Citation:**"CHAPTER 4. FINDINGS." National Academies of Sciences, Engineering, and Medicine. 2019.

**Suggested Citation:**"CHAPTER 4. FINDINGS." National Academies of Sciences, Engineering, and Medicine. 2019.

**Suggested Citation:**"CHAPTER 4. FINDINGS." National Academies of Sciences, Engineering, and Medicine. 2019.

**Suggested Citation:**"CHAPTER 4. FINDINGS." National Academies of Sciences, Engineering, and Medicine. 2019.

**Suggested Citation:**"CHAPTER 4. FINDINGS." National Academies of Sciences, Engineering, and Medicine. 2019.

**Suggested Citation:**"CHAPTER 4. FINDINGS." National Academies of Sciences, Engineering, and Medicine. 2019.

**Suggested Citation:**"CHAPTER 4. FINDINGS." National Academies of Sciences, Engineering, and Medicine. 2019.

**Suggested Citation:**"CHAPTER 4. FINDINGS." National Academies of Sciences, Engineering, and Medicine. 2019.

**Suggested Citation:**"CHAPTER 4. FINDINGS." National Academies of Sciences, Engineering, and Medicine. 2019.

**Suggested Citation:**"CHAPTER 4. FINDINGS." National Academies of Sciences, Engineering, and Medicine. 2019.

**Suggested Citation:**"CHAPTER 4. FINDINGS." National Academies of Sciences, Engineering, and Medicine. 2019.

**Suggested Citation:**"CHAPTER 4. FINDINGS." National Academies of Sciences, Engineering, and Medicine. 2019.

**Suggested Citation:**"CHAPTER 4. FINDINGS." National Academies of Sciences, Engineering, and Medicine. 2019.

**Suggested Citation:**"CHAPTER 4. FINDINGS." National Academies of Sciences, Engineering, and Medicine. 2019.

**Suggested Citation:**"CHAPTER 4. FINDINGS." National Academies of Sciences, Engineering, and Medicine. 2019.

**Suggested Citation:**"CHAPTER 4. FINDINGS." National Academies of Sciences, Engineering, and Medicine. 2019.

**Suggested Citation:**"CHAPTER 4. FINDINGS." National Academies of Sciences, Engineering, and Medicine. 2019.

**Suggested Citation:**"CHAPTER 4. FINDINGS." National Academies of Sciences, Engineering, and Medicine. 2019.

**Suggested Citation:**"CHAPTER 4. FINDINGS." National Academies of Sciences, Engineering, and Medicine. 2019.

**Suggested Citation:**"CHAPTER 4. FINDINGS." National Academies of Sciences, Engineering, and Medicine. 2019.

**Suggested Citation:**"CHAPTER 4. FINDINGS." National Academies of Sciences, Engineering, and Medicine. 2019.

**Suggested Citation:**"CHAPTER 4. FINDINGS." National Academies of Sciences, Engineering, and Medicine. 2019.

**Suggested Citation:**"CHAPTER 4. FINDINGS." National Academies of Sciences, Engineering, and Medicine. 2019.

**Suggested Citation:**"CHAPTER 4. FINDINGS." National Academies of Sciences, Engineering, and Medicine. 2019.

**Suggested Citation:**"CHAPTER 4. FINDINGS." National Academies of Sciences, Engineering, and Medicine. 2019.

**Suggested Citation:**"CHAPTER 4. FINDINGS." National Academies of Sciences, Engineering, and Medicine. 2019.

**Suggested Citation:**"CHAPTER 4. FINDINGS." National Academies of Sciences, Engineering, and Medicine. 2019.

**Suggested Citation:**"CHAPTER 4. FINDINGS." National Academies of Sciences, Engineering, and Medicine. 2019.

**Suggested Citation:**"CHAPTER 4. FINDINGS." National Academies of Sciences, Engineering, and Medicine. 2019.

**Suggested Citation:**"CHAPTER 4. FINDINGS." National Academies of Sciences, Engineering, and Medicine. 2019.

**Suggested Citation:**"CHAPTER 4. FINDINGS." National Academies of Sciences, Engineering, and Medicine. 2019.

**Suggested Citation:**"CHAPTER 4. FINDINGS." National Academies of Sciences, Engineering, and Medicine. 2019.

**Suggested Citation:**"CHAPTER 4. FINDINGS." National Academies of Sciences, Engineering, and Medicine. 2019.

**Suggested Citation:**"CHAPTER 4. FINDINGS." National Academies of Sciences, Engineering, and Medicine. 2019.

**Suggested Citation:**"CHAPTER 4. FINDINGS." National Academies of Sciences, Engineering, and Medicine. 2019.

**Suggested Citation:**"CHAPTER 4. FINDINGS." National Academies of Sciences, Engineering, and Medicine. 2019.

**Suggested Citation:**"CHAPTER 4. FINDINGS." National Academies of Sciences, Engineering, and Medicine. 2019.

**Suggested Citation:**"CHAPTER 4. FINDINGS." National Academies of Sciences, Engineering, and Medicine. 2019.

**Suggested Citation:**"CHAPTER 4. FINDINGS." National Academies of Sciences, Engineering, and Medicine. 2019.

**Suggested Citation:**"CHAPTER 4. FINDINGS." National Academies of Sciences, Engineering, and Medicine. 2019.

**Suggested Citation:**"CHAPTER 4. FINDINGS." National Academies of Sciences, Engineering, and Medicine. 2019.

**Suggested Citation:**"CHAPTER 4. FINDINGS." National Academies of Sciences, Engineering, and Medicine. 2019.

**Suggested Citation:**"CHAPTER 4. FINDINGS." National Academies of Sciences, Engineering, and Medicine. 2019.

Below is the uncorrected machine-read text of this chapter, intended to provide our own search engines and external engines with highly rich, chapter-representative searchable text of each book. Because it is UNCORRECTED material, please consider the following text as a useful but insufficient proxy for the authoritative book pages.

23 CHAPTER 4. FINDINGS INTRODUCTION This project addressed the issues associated with the sensitivity of pavement performance to base layers and subgrade for flexible and rigid pavements. Therefore, necessary enhancements to the Pavement ME Design procedures are needed to better reflect the influence of subgrade and unbound layers (properties and thicknesses) on the performance. These enhancements include modifications of the models contained in Pavement ME Design and the development of new models. Following the research plan presented in Chapter 3, this chapter gives detailed explanations for each aspect in the proposed plan. SWCC OF BASE AND SUBGRADE FOR FLEXIBLE AND RIGID PAVEMENTS The SWCC represents a relationship between soil suction and water content (or degree of saturation). The soil suction due to the presence of water between soil particles is an important indicator, which is directly related to the strength, volume change, and fluid flow characteristics of unsaturated soil material (77). The unbound material is usually in an unsaturated condition when it is exposed in the field. Currently, there are two common ways to determine the SWCC of unbound material. One way is by conducting laboratory tests (e.g., filter paper test and pressure plate test) to measure the matric suction at different water contents. This experimental method is usually time consuming and requires special test equipment and test expertise. An alternative way is to correlate the suction-water relationship to other unbound material indicators, such as grain size distribution (GSD), particle diameter, porosity, liquid limit (LL) and PI, mean annual air temperature (MAAT), etc. There are many prediction methods available, which can be classified into three categories (94): ï· Category 1: Direct prediction of SWCC, which is a curve-fitting method to correlate the soil properties (e.g., GSD and porosity) to the matric suction at the different tested water contents (95â98). ï· Category 2: Prediction of SWCC model parameters, which is a statistical method to correlate the soil properties (e.g., GSD, LL, PI, and MAAT) to the parameters of existing SWCC models (e.g., Fredlund-Xing equation). These estimated model parameters are further used to generate the SWCC of unbound material (94, 99â102). ï· Category 3: Development of a micromechanics-based SWCC model, which uses the microstructure of unbound material (e.g., porosity distribution) and the water-particle contact (or capillary rise) model to determine the matric suction (103â105). Among these approaches, Category 1 usually has a low prediction accuracy, since it ignores the existing relationship between matric suction and water content. Category 2 recognizes the existing SWCC models and uses an indirect method to predict the SWCC via soil properties. The prediction accuracy of Category 2 relies on the selected statistical method and the quantity of the database used in the analysis. However, most of the aforementioned models were developed based on either a simple regression analysis or a limited number of test results (101). Compared to Categories 1 and 2, Category 3 is still in a state-of-development status, which is limited to the soils with large pores (94). Therefore, approaches in Category 2 are more accurate and feasible to predict the SWCC of unbound material, as long as the issues of statistical method and quantity of database are solved.

24 Currently, ANN approach becomes a more popular tool for the development of prediction models. Compared to regression models, the main advantage of the ANN approach is that it can capture nonlinear and complex scattered relationships between input and output parameters. In the past, several studies attempted to use the ANN approach to predict the SWCC of soil, which mainly focused on the Category 1âbased method (i.e., direct prediction of SWCC) (106â109). Currently, there is no Category 2âbased ANN model available to predict the SWCC of soil. The main reason is that the development of Category 2âbased ANN model requires a much larger database of SWCC results. In addition, most of these existing ANN models dealt with unsaturated soil, but few ANN models can predict the SWCC of unbound granular material that is commonly used in pavement engineering. To overcome the aforementioned problems, a Category 2âbased ANN approach is recommended for estimating the SWCC of unbound material. In this project, a large database collected from the NCHRP 9-23A project is used to develop the ANN models. The Fredlund and Xing equation is used to estimate the relationship between matric suction and water content, as shown in Eqs. 4.1 and 4.2 (110): f fb c 1(h) x[ ] {ln[e ( ) ]} f S C h a ï½ ï« (4.1) where C(h) is a correction factor defined as: 5 ln(1 ) C( ) 1 1.45 10ln[1 ( )] r r h hh x h ï« ï½ ï ï« (4.2) where S is the degree of saturation (unit: %); h is the soil matric suction (unit: psi); and af, bf, cf, and hr are soil fitting parameters. The four fitting parameters are the outputs of ANN models and are predicted by the input variables. Development of ANN Models The ANN approach is a computer-based information processing technique, which allows the correlations to be established between the input variables Xi and the output variables Yj through the inter-connected neurons (i.e., weight factor, wji). Note that the input variables Xi and the output variables Yj are usually normalized to xi and yj, respectively. The correlations developed by the ANN models between the normalized input parameters xi and the normalized output variables yj are shown in Eq. 4.3: 1 n j ji i i y f w x ï½ ï¦ ï¶ ï§ ï· ï¨ ï¸ ï½ ï¥ (4.3) where f is a transfer function, which normally uses a sigmoidal, Gaussian, or threshold functional form, and wji are the unknown weight factors. Developing a neural network model specifically refers to the determination of the weight factors wji in Eq. 4.3. In this project, the output variables Yj represent the four fitting parameters in the Fredlund-Xing equation (i.e., af, bf, cf, and hr). The

25 input variables Xi are selected from the SWCC-related material indicators, which will be elaborated in the following section. In general, the development of ANN models includes two critical steps: 1) data collection, and 2) construction of ANN architecture. Data Collection The database used in this project is collected from the NCHRP 9-23A project entitled âA National Catalog for Subgrade SWCC Default Inputs and Selected Soil Properties for Use with the ME-PDGâ carried out at Arizona State University in 2010 (111). In addition to the measured SWCC fitting parameters, the database includes material size distribution, soil LL, PI, and saturated volumetric water content. To better evaluate the influence of soil plasticity, the collected database is separated into two groups: plastic soil data set (PI>0), and non-plastic soil data set (PI=0). Construction of ANN Architecture As shown in Figure 1a and Figure 1b, two three-layered neural network architectures consisting of one input layer, one hidden layer, and one output layer were constructed for plastic soil and non-plastic soil, respectively. The input variables in the ANN models were selected based on previous studies and the mechanics of unsaturated soil. The SWCC of a soil is greatly dependent on its pore structure (110, 112, 113). The pore size is related to the height of each capillary tube using the YoungâLaplace equation, which is equivalent to a suction value (77). The capillary height of the corresponding pore can be expressed in terms of the radius of that pore. Thus, soil GSD (i.e., percent passing No. 4 sieve and percent passing No. 200 sieve) and Ï´, and Ñ°, and saturated volumetric content were selected as input parameters. Saturated volumetric water content is a measure of the total porosity in soil structure. The equilibrium soil suction of plastic soil is proportional to its specific surface area (98). The PI is a simple indicator that reflects the specific surface area of soil. The relationship between soil surface area and PI changes, particularly when the LL is above 30 percent to 40 percent (114). Hence, LL was also selected as an input parameter in the ANN models. Previous studies investigated the effects of temperature on soil suction-water retention relationship (115â117). Soil water is composed of continuous water and isolated pockets of water. When the temperature increases, water flows from isolated pockets to the continuous phase, which results in a shift in the SWCC. In this study, the MAAT is identified as one of the inputs for ANN models and is collected from the NRCS soil data mart. In Figure 1a, the input variables for plastic soil include percent passing No. 4 sieve, percent passing No. 200 sieve, LL, PI, saturated volumetric water content (Sat. vol. wc.) and local MAAT (Unit: Â°C). Figure 1b shows that the input variables for non-plastic soil are D30, D60, D90, Ï´, and Ñ°, saturated volumetric water content, and local MAAT (unit: Â°C). Herein, Ó¨ and Ñ° are estimated by fitting a power law model (Ó¨xÑ°) to the curve of the cumulative percent passing versus sieve size in mm. The hidden layer assigns 20 neurons to establish the connection between the output layer and the input layer. Both the non-plastic and plastic soil models use the sigmoidal transfer function, as shown in Eq. 4.4 (118): ï¨ ï© ï¨ ï© 1 1 expi i f I Iïª ï½ ï« ï (4.4)

26 where Ii is the input quantity. Ï is a positive scaling constant, which controls the steepness between the two asymptotic values 0 and 1. The ANN model determines these weight factors wji through the two major functions: training and validating. The training data set is used to determine the trial weight factors, wji, and the validating data set is employed to examine the accuracy of the model prediction. In this study, 80 percent of the data set is used for training and 20 percent of the data set for validation. The training algorithm uses the Levenberg-Marquardt back propagation method to minimize the mean squared error (MSE). The gradient descent weight function is employed as a learning algorithm to adjust the weight factors wji. (a) Plastic soil (b) Non-plastic soil Figure 1. Illustration of Three-Layered Neural Network Architecture. Prediction of SWCC Fitting Parameters Using ANN Model The development of ANN models is to predict the SWCC parameters (i.e., af, bf, cf, and hr) in the Fredlund-Xing equation. These parameters were correlated to the soil physical properties. In this study, the ANN models were programmed using the Matlab software. The

27 database used for training of ANN models consists of 3600 samples of plastic soil and 250 samples of non-plastic soil. A statistical analysis was performed to determine the root mean squared error (RMSE) and R2 associated with the predicted SWCC fitting parameters. The RMSE and R2 are computed by using Eqs. 4.5 and 4.8: Root mean squared error: 2( )x y RMSE n ï ï½ ï¥ (4.5) Residual sum of squares: 2( )resSS x yï½ ïï¥ (4.6) Total sum of squares: 2( )totSS x xï½ ïï¥ (4.7) Coefficient of determination: 2 1 res tot SSR SS ï½ ï (4.8) where x is the measured SWCC fitting parameter. y is the predicted SWCC fitting parameter. n is the number of data points. Figure 2 and Figure 3 compare the predicted SWCC parameters of plastic and non-plastic soils, respectively, against the measured ones using the training data sets. As shown in Figure 2 and Figure 3, the predicted af, bf, cf, and hr parameters are closely coincident with the measured results. This indicates that the developed ANN models have the desirable accuracy to predict the SWCC for plastic and non-plastic soils. The next section compares these ANN models against the existing regression models for model prediction accuracy.

28 (a) (b) (c) (d) Figure 2. Comparison of Measured versus Predicted SWCC Fitting Parameters Using ANN Model for Plastic Soils.

29 (a) (b) (c) (d) Figure 3. Comparison of Measured versus Predicted SWCC Fitting Parameters Using ANN Model for Non-plastic Soils. Comparison of ANN Model with Other Regression Models Several regression models have been developed to predict the SWCC fitting parameters. In this study, a comparison analysis is performed between the developed ANN models and those of the existing prediction models, including the Zapata model (94) and Perera model (100). Table 12 lists the equations of these regression models. The predictability of the SWCC fitting parameters using the regression models above are analyzed in this study. Table 13a and Table 13b list the computed RMSE and R2 values for each model. The values of R2 and RMSE indicate the model prediction accuracy. The higher R2 value and the smaller RMSE value represent a higher prediction accuracy for the developed model (73). As presented in Table 13a and Table 13b, the RMSE values associated with all of the existing SWCC parameter prediction models are very high, and the corresponding R2 values are extremely low. This indicates that the Zapata and Perera models have low prediction accuracy of the SWCC fitting parameters for both plastic and non-plastic soils. Among these two models, the computed RMSE and R2 values are comparable with each other. This demonstrates that these models have the same level of prediction accuracy. Compared to these existing models, the developed ANN models have much higher R2 values and smaller RMSE values. This indicates

30 that the developed ANN models outperform these existing regression modelsâ prediction accuracies. Table 12. List of Existing SWCC Fitting Parameter Prediction Models. Plastic Soil Non-plastic Soil Zapata Model 3.350.00364( ) 4( ) 11fa wPI wPIï½ ï« ï« 0.142.313( ) 5f f b wPI c ï½ ï ï« 0.4650.0514( ) 0.5fc wPIï½ ï« 0.0186( )32.44 wPIr f h e a ï½ 0.751 600.8627( )fa D ïï½ 7.5fb ï½ 600.1772ln( ) 0.7734fc Dï½ ï« 4 60 1 9.7 r f h a D eï ï½ ï« where: w = percent of material passing No. 200 sieve PI = Plasticity Index D60 = material diameter corresponding to 60% passing by weight of material Perera Model 3 2 . 8 3 5 { l n ( ) } 3 2 . 4 3 8fa w P Iï½ ï« 0.31851.421( )fb wPI ïï½ - 0 . 2 1 5 4 { l n ( ) } 0 . 7 1 4 5fc w P Iï½ ï« 500rh ï½ 1 . 1 4 0 . 5fa aï½ ï 6 4.34 20 200 30 100 2.79 14.11log( ) 1.9 10 7log( ) 0.055 a D x P D D ïï½ ï ï ï ï« ï« 0 . 9 3 6 3 . 8fb bï½ ï 0.57 1.19 0.190 200 0 200 1 10 {5.39 0.29ln[P ( )] 3D 0.021P }mDb D ï½ ï ï« ï« 1 90 60 30 [log( ) - log( )] m D D ï½ , 2 30 10 20 [log( ) log( )] m D D ï½ ï 0.758 100.26 1.4 c fc e Dï½ ï« 1.15 2 1log( ) (1 ) f c m b ï½ ï ï 100rh ï½ where: D10, D20, D30, D60, D90, D100 = material diameter corresponding to 10%, 20%, 30%, 60%, 90%, and 100% passing by weight of material, respectively

31 Table 13. Prediction Accuracy of SWCC Fitting Parameter Models. (a) Plastic Soil Model RMSE R 2 af bf cf hr af bf cf hr Zapata 4.86 0.46 0.21 2798 0.0078 0.059 0.29 0.003 Perera 7.86 0.51 0.33 2544 0.0033 0.04 0.21 0.008 ANN 1.51 0.12 0.068 0.395 0.68 0.87 0.78 0.45 (b) Non-plastic Soil Although the developed ANN models can accurately predict the SWCC fitting parameters, it is still necessary to examine the prediction accuracy of the SWCC. In the laboratory, these fitting parameters (i.e., af, bf, cf, and hr) in the Fredlund-Xing equation are calculated by fitting the SWCC with experimentally available suction-saturation data points. A nonlinear least squared regression analysis is performed to determine the best fit parameters. In this project, the predicted SWCC fitting parameters are input into the Fredlund-Xing equation to estimate the soil suction at various water contents. The estimated soil suction values are then compared against the measured ones to investigate the model prediction accuracy. To evaluate the suction predictability over the full range of the saturation level, the matric suction values are calculated at three different degrees of saturation (i.e., 10 percent, 40 percent, and 80 percent). The predicted matric suction values using different models were compared with experimental data at the same saturation level. Figure 4 shows the plots of measured versus predicted matric suction for plastic soils. Figure 4a, Figure 4b, and Figure 4c correspond to comparisons made by experimental data with the Zapata model, Perera model, and ANN model, respectively. Compared to the two existing models, the developed ANN model has the least RMSE value and the highest R2 value. This indicates that the developed ANN model outperforms other models to accurately estimate the matric suction of plastic soils at various degrees of saturation. The obtained R2 value of 0.95 demonstrates that the developed ANN model can accurately predict the SWCC of plastic soil. Figure 5 presents the measured versus predicted suction for non-plastic soils. As seen from Figure 5a and Figure 5b, the predicted matric suction values using the existing models deviate significantly from the measured values. This indicates that the existing models yield inaccurate matric suction of non-plastic soils. As shown in Figure 5c, the developed ANN model has an R2 value of 0.91, which significantly improves the prediction accuracy of the matric suction for the non-plastic soils. Researchers concluded from Figure 4 and Figure 5 that the developed ANN models are capable of accurately predicting the SWCC of both plastic and non- plastic soils. Model RMSE R 2 af bf cf hr af bf cf hr Zapata 5.76 5.34 0.33 3000 0.009 0.00 0.018 0.01 Perera 7.95 2.19 0.68 2899 0.14 0.12 2.42E-04 0.00 ANN 0.033 0.028 0.063 0.398 0.92 0.76 0.83 0.48

32 (a) Zapata Model (b) Perera Model (c) ANN Model Figure 4. Comparison of Measured versus Predicted Suction at Various Saturation Levels for Plastic Soils.

33 (a) Zapata Model (b) Perera Model (c) ANN Model Figure 5. Comparison of Measured versus Predicted Suction at Various Saturation Levels for Non-plastic Soils. Validation of the Developed ANN Models The validation of the prediction accuracy of the developed ANN models involves two steps: 1) the validation through the collected data from NCHRP 9-23A database, and 2) the validation via the independent data from other literature sources. At first, a new data set of 500 plastic soils and 33 non-plastic soils were selected from the NCHRP 9-23A database. The selected soil physical and climatic properties were input into the developed ANN models. Figure 6 and Figure 7 compare the measured SWCC fitting parameters to those predicted by the ANN models for the selected plastic and non-plastic soils, respectively. The ANN model predictions have relatively small RMSE and high R2 values for both cases. The determined R2 values for plastic soils are in the range from 0.44 to 0.87, non-plastic soils show R2 ranging from 0.49 to 0.87. This indicates that the model predicted SWCC fitting parameters match well with the measured results. The model predicted fitting parameters are then input into the Fredlund-Xing equation to estimate the matric suction values at the selected degrees of saturation.

34 (a) (b) (c) (d) Figure 6. Validation of Measured versus Predicted SWCC Fitting Parameters Using ANN Model for Plastic Soils.

35 (a) (b) (c) (d) Figure 7. Validation of Measured versus Predicted SWCC Fitting Parameters Using ANN Model for Non-plastic Soils. Figure 8a and Figure 8b plot the measured versus predicted matric suction values at 10 percent, 40 percent, and 80 percent saturation level for plastic and non-plastic soils, respectively. Results show that the predicted suction values fit well with the measured ones. Both plastic and non-plastic soil data sets show a R2 value greater than 0.90. This validates that the developed ANN models provide a desirable prediction accuracy of the SWCC.

36 (a) Plastic soil (b) Non-plastic soil Figure 8. Validation of Measured versus ANN Predicted Suction at Various Saturation Levels for Unbound Materials. To further validate the prediction accuracy of the developed ANN models, the independent data sources (i.e., two plastic and two non-plastic soil data) were collected from the literature (18, 119, 120). All these soil physical and climatic properties are shown in Table 14a and Table 14b, which are used as inputs of the developed ANN models. Table 14. Input Parameters Collected from Literature for Model Validation. (a) Plastic soil Reference Soil source #4 sieve #200 sieve LL PI Sat. vol. wc* MAAT(Â°C) (119) Louisiana 100 45 23 12 5.4 18.5 (18) Red Lake Falls 100 93.9 29 10 8.17 3.9 Note: *Sat.vol.wc = Saturated volumetric water content (b) Non-plastic soil Reference Soil source D30 D60 D90 Sat. vol. wc MAAT(Â°C) (119) Mississippi 0.075 0.1021 1.41 18.76 17.9 (120) Torpsbruk 0.0911 0.6151 11.23 15.73 6.91 Figure 9a and Figure 9b compare the measured and predicted SWCCs for plastic and non-plastic soils, respectively. As illustrated, the SWCC curves generated by using the predicted parameters from the predicted ANN model produce good fits with the experimental curves. In Figure 9a, the model predicted matric suction values of plastic soils (i.e., Red Lake Falls and Louisiana soils) are generally coincident with the experimental results. A small deviation is observed in the residual region for the Louisiana soil. This is due primarily to the effect of the cf model, which has a lower R2 than the af and bf parameters. Similarly, Figure 9b presents that the predicted matric suction values of non-plastic soils (i.e., Torpsbruk and Mississippi soils) were in good agreement with the measured results. Considering the SWCC curve generated from experiment involves the high level of variability, the difference between predicted and experimental curves were sufficiently accurate.

37 (a) Plastic Soil (b) Non-plastic soil Figure 9.Â Comparison of Measured versus Predicted SWCC Curves for Unbound Materials.Â EQUILIBRIUM SUCTION OF BASE AND SUBGRADE FOR FLEXIBLE AND RIGID PAVEMENTS In Pavement ME Design, the control moisture level is determined from the depth of ground water table. This approach has limitations when there is no water table data available and moisture level is controlled by other factors. Researchers proposed a new approach stating that when the water table is deeper than the moisture level underneath the pavement surface is controlled by the equilibrium suction and its variations. A relationship has been developed in this study to predict the equilibrium suction from readily available parameters (i.e., TMI, PI, and dry suction value at surface ]) based on the data collected from the continental United States. The TMI (121) is a moisture index that reflects the aridity or humidity of the soil and climate. This index is calculated using actual data such as the collective effects of precipitation, evapotranspiration, soil water storage, etc. The value of TMI can be found from the TMI map of the United States (Figure 10) (121).

38 Figure 10. TMI Distribution in United States (121). However, since the development of the original TMI map in 1948, various studies have been performed to improve the water balance approach and simplify the calculation process. TMI values vary depending on the methods adopted for calculation and the time span of data used in computation (122). Olaiz et al. (123) performed a comparison study of the TMI contour maps developed by Thornthwaite et al. (121, 124â126) and concluded that the method developed by Witczak et al. (126) matched most closely with the original moisture balance approach of the Thornthwaite method. Witczakâs equations were used in this study to calculate the TMI value in the continental United States. Detailed description of the development of the TMI map in the GIS platform is presented in the next section. GIS-Based Map of TMI The spatial climatic data sets (precipitation and temperature) that were used to generate a GIS-based map of TMI were collected from Oregon State Universityâs Parameter-elevation Regressions on Independent Slopes Model (PRISM) Climate Group, named for the PRISM climate mapping system (127â129). In situ measurements of temperature and precipitation from 5852 weather stations across the continental United States are obtained from the National Oceanic and Atmospheric Administration climate database and ingested into the PRISM statistical mapping system. The PRISM products use a weighted regression scheme to account for complex climate regimes associated with orography, rain shadows, temperature inversions, slope aspect, coastal proximity, and other factors and develop climatological normal (from 1981 to 2010) at 30-arcsec (800 meters) resolution. Figure 11 shows the average yearly precipitation (P) and temperature map of the United States in GIS platform.

39 (a) (b) Figure 11. Average Annual (a) Precipitation; (b) Mean Temperature GIS Map (1981 to 2010). The potential evapotranspiration was calculated empirically using the mean monthly temperature (ti) of the site (130). A yearly heat index (Hy) is determined first, as shown in Eq. 4.9: 12 1.54 1 (0.2 )y i i H t ï½ ï½ï¥ (4.9) where, yH is the yearly heat index. it is the mean monthly temperature (i denotes the given month; from 1 to 12).

40 The monthly potential evapotranspiration is quantified using Eq. 4.10: 101.6( )aii y tpe H ï½ (4.10) where 7 3 5 26 . 7 5 1 0 7 . 7 1 1 0 0 . 0 1 7 9 2 0 . 4 9 2 3 9y y ya x H x H Hï ïï½ ï ï« ï« . The monthly ipe values determined using Eq. 4.10 represents a quantity of 30-day month with 12-hour daylight duration. To consider the variations in daylight due to latitude and varying days in each month, a day length correction factor is introduced (131) to obtain the corrected potential evapotranspiration in Eq. 4.11: ( ) 30 i i i i d nPE peï½ (4.11) where iP E is the monthly corrected potential evapotranspiration (unit: cm). di is the daylight duration for each month i (unit: hrs). in is the number of days in the given month. The yearly average evapotranspiration is calculated by taking the sum of iP E for 12 months: 12 1 y i i PE PE ï½ ï½ï¥ (4.12) Figure 12 depicts the calculated average yearly potential evapotranspiration map of the United States on GIS. Figure 12. Average Annual Potential Evapotranspiration GIS Map (1981 to 2010).

41 Finally, the TMI value is calculated following Eq. 4.13 according to Witczak et al. (126) as shown in Figure 13: 75( 1) 10 y PTMI PE ï½ ï ï« (4.13) Figure 13. GIS-Based Contour Map of TMI (1981 to 2010). In addition, researchers have developed a fundamental method of calculating equilibrium suction that involves a steady state diffusion equation proposed by Mitchell (93), as shown in Eq. 4.14: 0( ) * n z eu z u u e ï° ï¡ ï ï½ ï± (4.14) where ( )u z is the suction at depth z expressed on the pF scale. eu is the equilibrium value of suction on the pF scale. 0u is the suction profile amplitude. n is the number of suction cycles per second (1 year = 365 Ã 24 Ã 60 Ã 60 seconds). ï¡ is the unsaturated soil diffusion coefficient (unit: cm2/sec). z is the depth below the surface (unit: cm). The terms dryu and w etu represent soil suction at the surface corresponding to an air dry state and soil field capacity, respectively. dryu is controlled by vegetation (4.5 pF) and by bare soil (5.7 pF). The term w etu corresponds to 3 pF suction level (132). The soil water characteristics curve prediction model developed by Saha et al. (133) was used in this study to determine the volumetric water content at specific suction levels. Finally, the equilibrium suction value is calculated for the continental United States following Eq. 4.14 and shown in Figure 14.

42 Figure 14. GIS-Based Equilibrium Suction Map. Next, the developed TMI and equilibrium suction map of United States is divided into 72332 map units. The mapunit shapefile of the United States is collected from the U.S. Department of Agriculture (USDA) NRCS STATSGO2 database. A zonal distribution analysis is performed in ArcGIS to determine the average TMI and pF value for each mapunit. The following section describes the relationship between TMI and equilibrium suction for various AASHTO soil classes. Relationship Between TMI and Equilibrium Suction All the NRCS mapunits are divided according to the AASHTO soil classes, and Figure 15 plots the equilibrium suction versus TMI values. For each soil class, the equilibrium suction value decreases with an increase of TMI. A polynomial trend line is fitted for each case to express the relationship between the equilibrium suction and TMI. The A-3 soil class has the lowest equilibrium suction range whereas the A-7-6 class shows the highest range. In general, the equilibrium suction increases with the increase of fine material and especially with clay fines. Equilibrium Suction (cm) <3000 3,000 - 4,000 4,000 - 5,000 5,000 - 6,000 >6000

43 (a) AASHTO soil type: A-1 (b) AASHTO soil type: A-2 (c) AASHTO soil type: A-3 (d) AASHTO soil type: A-4 (e) AASHTO soil type: A-6 (f) AASHTO soil type: A-7-6 Figure 15. TMI versus Equilibrium Suction (pF) (cm) for (a) A-1; (b) A-2; (c) A-3; (d) A-4; (e) A-6; and (f) A-7-6 Soil. Ue = 3.527e-3E-04(TMI) 3 3.5 4 4.5 5 5.5 â60 â20 20 60 100 Eq ui lib riu m S uc tio n (p F) TMI St. Deviation,Ï = 0.18 Ue = 3.582e-3E-04(TMI) 3 3.5 4 4.5 5 5.5 â60 â20 20 60 100 Eq ui lib riu m S uc tio n (p F) TMI St. Deviation,Ï = 0.25 Ue = 3.370e-3E-04(TMI) 3 3.5 4 4.5 5 5.5 â40 0 40 80 Eq ui lib riu m S uc tio n (p F) TMI St. Deviation,Ï = 0.09 Ue= 3.621e-3E-04(TMI) 3 3.5 4 4.5 5 5.5 â60 â20 20 60 100 Eq ui lib riu m S uc tio n (p F) TMI St. Deviation,Ï = 0.17 Ue = 3.706e-9E-05(TMI) 3 3.5 4 4.5 5 5.5 â60 â20 20 60 100 Eq ui lib riu m S uc tio n (p F) TMI St. Deviation,Ï = 0.26 Ue = 3.786e-3E-04(TMI) 3 3.5 4 4.5 5 5.5 â60 â20 20 60 100 Eq ui lib riu m S uc tio n (p F) TMI St. Deviation,Ï = 0.34

44 A linear regression analysis was performed to predict the equilibrium suction value from the TMI and PI. Eqs. 4.15 and 4.16 express the relationship for vegetative and bare soil, respectively. Figure 16 shows the plot of the calculated versus predicted equilibrium suction using Eqs. 4.15 and 4.16. An R2 value of 0.8 indicates that the regression equations have very good prediction accuracy: ue (pF) [Vegetative Soil] = 0.0199*PI-0.0047*TMI+3.3976 (4.15) ue (pF) [Bare Soil] = 0.0214*PI-0.00607*TMI+3.5318 (4.16) Figure 16. Calculated versus Predicted Equilibrium Suction (pF). RESILIENT MODULUS OF BASE AND SUBGRADE FOR FLEXIBLE AND RIGID PAVEMENTS Unbound base and subgrade materials beneath the pavement surface are characterized by the MR in the Pavement ME Design (67), which directly affects the design and analysis of flexible pavement structures. The advantage of MR over other material properties is that it represents the response to dynamic loading coming from vehicle movement. It is defined as the ratio of repeated cyclic loading ( c y cï³ ) divided by recoverable strain ( rï¥ ) as shown in Eq. 4.16: cyc R r M ï³ ï¥ ï½ (4.16) The MR can be estimated in three different approaches: 1) conduct lab tests with compacted soil specimens (e.g., repeated load triaxial test), 2) follow a back calculation method from an in situ device (e.g., FWD), and 3) predict MR model coefficients from soil physical properties. The first and second approaches are expensive, time consuming, and require expert labor. An MR model based on the third approach was used in this study to characterize the unbound and subgrade material under different conditions. Various models have been developed to predict the MR of unbound base and subgrade materials. The most widely used generalized model developed in the NCHRP Project 1-28A, is shown in Eq. 4.17 (67):

45 32 k1 1 ( ) ( 1) k oct R a a a IM k P P P ï´ ï½ ï« (4.17) where 1I is first invariant of the stress tensor, octï´ is octahedral shear stress, aP is atmospheric pressure, 1k , 2k , and 3k are regression coefficients. However, several studies have shown that the MR of unbound granular material is not only stress dependent but also moisture dependent (134). The AASHTO employed a moisture- dependent MR model, which is later adopted by the Pavement ME Design: m opt log 1 exp[ln k (S S )] R Ropt M b aa bM a ï ï½ ï« ïï« ï« ï (4.18) where RM is resilient modulus at a given degree of saturation, R o p tM is resilient modulus at a reference condition, a and b are minimum and maximum of log(MR/MROPT), respectively, km is regression coefficient, o p t(S S )ï is variation of degree of saturation, expressed as a decimal. Other researchers (19, 61, 135) also proposed different moisture-dependent MR models for both unbound aggregates and subgrade. However, previous studies also demonstrated the importance of the matric suction in the MR model. The Pavement ME Design incorporated the suction in the MR model through NCHRP Project 9-23A (111): 32 1 * _ ( ) ( 1)kk octR a a a wc matric suction M k P P P ï´ï± ï« ï½ ï« (4.19) where ï± is bulk stress, wc is water content. There have been several other studies recommending the use of soil suction in modeling MR (17, 18, 40, 58, 136). The matric suction can be different for the same variation of degree of saturation (i.e., S-Sopt), from one material to another (66). Therefore, the degree of saturation and soil suction alone cannot accurately estimate the change of MR due to moisture. In order to consider the moisture and suction dependency of resilient modulus characteristics, a new model is proposed by Lytton (60): 321 1 3( ) ( )kkm octy a a a I fhE k P P P ï± ï´ï ï½ (4.20) where 1I is first invariant of the stress tensor, aP is atmospheric pressure, ï± is volumetric water content, mh is matric suction in the base material, f is saturation factor, 1 < f < 1/Î¸, octï´ is octahedral shear stress, 1k , 2k , and 3k are regression coefficients. In this model, 1I and octï´ vary with the stress state, and mh is related to the moisture content of unbound aggregates. Appendix D provides more details on this model. Three steps are followed to determine the coefficients 1k , 2k , and 3k using Eq. 4.20:

46 1. Matric suction hm in Eq. 4.21 is estimated first at the tested specimen water content. A separate ANN model is used to predict the coefficients of the SWCC parameters in the Fredlund-Xing equation as described in subsequent sections. 2. Volumetric water content, Î¸, and saturation factor, f, are calculated for each base and subgrade material from the collected moisture content, MDD, and specific gravity data. The saturation factor is an indicator of degree of saturation (60). It is multiplied with Î¸ and hm to express the stress on soil mineral skeleton in the transition zone. The transition zone occurs between the air entry point suction value and the unsaturation point suction value. The saturation corresponding to the air entry point and the unsaturation point were considered as 100 percent and 85 percent, respectively, in this study. Eq. 4.21 is followed to calculate the saturation factor, f in the transition zone: 85 11 ( 1) 15 Sf ï± ï ï½ ï« ï (4.21) where S is degree of saturation (unit: %), ï± is volumetric water content (unit: decimal); it is considered to be 2 pF and 3.5 pF at the air entry point and the unsaturation point, respectively. 3. Coefficients k1, k2, and k3 are calculated from collected MR data with different combinations of confining pressure and deviator stresses. A Matlab code was generated to fit the curve for each base and subgrade material. Researchers have developed many correlation models to predict k1, k2, and k3 from soil gradation, index, and strength properties (137â140). However, most of these models are either confined to a limited number of data sets or have a poor prediction accuracy. As a result, an ANN model is developed in this study to predict the coefficients from soil properties. The main advantage of a neural network model over nonlinear regression models are that it can capture complex nonlinear scattered relationships between input and output parameters and train the model based on the evaluation of error function. The application of the ANN approach for predicting MR of subgrade soil has been studied by a few researchers (141, 142). Soil index properties (i.e., material percent passing No. 200 sieve, PI, OMC, MDD, variation in compaction water content) and stress variables (i.e., confining pressure, bulk stress, deviator stress) were used as input variables to predict MR (141). A different neural network analysis was performed with nine different soil sources in Georgia (142). A similar study has been performed for ln(k1), k2, and k3 coefficients by Nazzal and Mohammad (138). Input variables for the correlation models consisted of material percent passing No. 200 sieve, LL, PI, MDD, test dry density (TDD), OMC, TMC, and clay%. However, previous ANN approaches have had some limitations. The resilient modulus values were predicted either as model output for a single stress combination or using only the locally available soil data as input (143). To overcome these limitations, an ANN model is developed in this project using a large database of unbound base and subgrade materials collected from the LTPP. Moreover, the resilient modulus equation coefficients (k1, k2, and k3) were predicted as output to cover any stress combination. Prediction of Suction at Test Specimen Water Content Matric suction of unbound granular base and subgrade material was required to accurately estimate the resilient modulus in an unsaturated condition as shown in Eq. 4.20. However, the quantity of available SWCC data of unbound base materials from the available

47 literature was not sufficient to train an ANN model. Hence the Soil Survey Geographic Database was used for data collection from USDA-NRCS. The database included soil physical and climatic properties as well as water content at specific suction pressures. Collected suction- saturation data points for each soil were fitted with the Fredlund-Xing equation. A Matlab code was generated to fit the Fredlund-Xing equation with measured data points and to calculate the SWCC fitting parameters. Out of a total 34,237 data sets on soils collected from the NRCS, soils that satisfy gravel characteristics according to Unified Soil Classification System (more than 50 percent material retained on the No. 200 sieve is defined as coarse grained material and more than 50 percent of the coarse grained material retained on the No. 4 sieve) were separated out to train and validate the ANN models for the SWCC fitting parameters. The 2598 plastic (PI>0) soil data and 311 non-plastic (PI=0) soil data were collected and used to train the neural network to predict the SWCC fitting parameters for unbound materials. The ANN models for the fitting parameters (i.e., af, bf, cf, and hr) of the SWCC in the Fredlund-Xing equation are detailed above. The hidden layer in the ANN model assigned 20 neurons to establish the connection between input and output layers. The 80 percent of the total data set was used for training and 20 percent for validation of ANN models. Water contents at 0.1, 0.33, and 15 bars (1 bar = 100 kPa) were available in the USDA-NRCS database. The predicted SWCC parameters were input into the Fredlund-Xing equation to estimate the water content at the same suction levels. The estimated water contents were compared against the measured values to investigate the prediction accuracy. Figure 17a and Figure 17b show the plots of the measured versus the predicted degrees of saturation for plastic and non-plastic soils, respectively. The developed ANN models had a very good prediction accuracy generating a R2 value of 0.91 for both plastic and non-plastic soils. Hence the same trained models were used to predict the SWCC fitting parameters for base and subgrade materials collected from the LTPP database. Predicted af, bf, cf, and hr parameters of the LTPP materials were input into the Fredlund-Xing equation to estimate matric suction of the specimen at the tested water content. (a) Plastic soils (b) Non-plastic soils Figure 17. Comparison of Measured versus Predicted Saturation (%) at 0.1, 0.33, and 15 Bars Suction Level for Unbound Granular Base Materials Using ANN Model.

48 Development of ANN Models for MR Model Coefficients The ANN approach is an adaptive information processing technique, which allows to establish the correlations between the input variables Xi and the output variables Yj through the inter-connected neurons (i.e., weight factor, wji). Note that the input variables Xi and the output variables Yj are usually normalized to xi and yj, respectively. The correlations developed by the ANN models between the normalized input parameters xi and the normalized output variables yj are shown in Eq. 4.22: 1 j n j ji i i by f w x ï½ ï¦ ï¶ ï«ï§ ï· ï¨ ï¸ ï½ ï¥ (4.22) where f is a transfer function, which normally uses a sigmoidal, Gaussian, or threshold functional form, wji and bj are the unknown weight factors and bias term, respectively. A neural network model specifically adjusts the weight factors wji and bias bj in Eq. 4.22 based on the minimum error function. In pavement engineering, the ANN approach is usually used to develop prediction models based on a large number of data collected from experiments and numerical analysis. ANN models have been successfully developed to predict the crack growth function (e.g., reflective cracking and top-down cracking) in asphalt concrete (144, 145). The ANN approach is also used to predict the geogrid-reinforced flexible pavement performance (118). In general, the development of ANN models includes two critical steps: 1) data collection; and 2) construction of ANN architecture. Data Collection A large collection of unbound base materials was available in the LTPP InfoPave database. Out of total 3010 unbound base materials, 717 materials were selected that pass gravel specifications. Researchers separated out 217 plastic (PI>0) and 500 non-plastic (PI=0) base materials to develop two different sets of ANN models for k1, k2, and k3. In this way, the subgrade ANN models for k1, k2, and k3 were developed based on non-plastic subgrade data from 180 LTPP sections and plastic subgrade data from 400 LTPP sections. MR tests were conducted on 15 different combinations of confining pressure and nominal maximum axial stress level. In addition to the MR test data, physical properties of base materials such as gradation, Atterberg limits, moisture content, density, and specific gravity were collected from the LTPP database. Construction of ANN Architecture Two three-layered ANN models, one for plastic and another for non-plastic soil were constructed as shown in Figure 18a and Figure 18b, respectively. Three layers of the ANN model consisted of a input layer, a hidden layer, and an output layer. Soil physical properties (i.e., percent of material passing 3/8 in. sieve, percent of material passing No. 200 sieve, PL, PI, OMC, MDD, and TMC) were introduced as input parameters for plastic soils. Input parameters for non- plastic soils included the percent of material passing the 3/8 in. sieve, percent of material passing the No. 200 sieve, Ï´, Ñ°, OMC, MDD, and TMC. The hidden layer assigned 10 neurons to establish the connection between the output layer and the input layer. The number of neurons was selected based on the number of data points. Too many hidden neurons allow the network to memorize instead of generalizing the training set (146). Both the non-plastic and plastic soil

49 models used the same sigmoidal transfer function as that for the SWCC ANN models, 80 percent of the data set was used for training and 20 percent of the data set for validation. (a) Plastic soil (b) Non-plastic soil Figure 18. Illustration of Three-Layered Neural Network Architecture (a) Plastic; (b) Non- plastic Soil. Prediction of MR Model Coefficients Using ANN Model The ANN model was developed to predict the MR model coefficients (i.e., k1, k2, and k3) from base physical properties. The coefficients of the MR model were correlated with base material physical properties and used to estimate the MR. Figure 19âFigure 20 show the comparison between measured and predicted MR model coefficients for plastic soils. The predicted MR coefficients show high R2 and low RMSE, which indicate a good prediction accuracy. The R2 values for MR coefficients of plastic base and subgrade materials were in the range of 0.66 to 0.92. Similarly, Figure 21âFigure 22 plotted the comparisons between the measured and predicted MR model coefficients for non-plastic soils.

50 The R2 value for k1, k2, and k3 coefficients of non-plastic base and subgrade materials were in the range of 0.51â0.97. Note that a high level of variability was involved in estimated coefficients due to variations in cyclic load triaxial test results and seasonal changes in selected physical properties (147). Hence, the predicted k1, k2, and k3 values need to be scrutinized in the validation process. If the validation data set provides comparable results, then the predicted coefficients can be sufficiently accurate to use in MR model equation. (a) k1 (b) k2 (c) k3 Figure 19. Predicted MR Model Coefficients of Plastic Base Materials from Physical Properties Using the ANN Approach.

51 (a) k1 (b) k2 (c) k3 Figure 20. Predicted MR Model Coefficients of Plastic Subgrade Materials from Physical Properties Using the ANN Approach.

52 (a) k1 (b) k2 (c) k3 Figure 21. Predicted MR Model Coefficients for Non-plastic Base Materials from Physical Properties Using the ANN Approach.

53 (a) k1 (b) k2 (c) k3 Figure 22. Predicted MR Model Coefficients for Non-plastic Subgrade Materials from Physical Properties Using ANN Approach. Comparison of ANN Model with Other Regression Models Several regression models have been developed to correlate the soil physical properties with MR model coefficients. In this study, the correlation models (137, 139, 148) were selected and compared against measured values to investigate the accuracy of the MR coefficients. They were named as Yau model (137), Malla model (139), and Soliman model (148), respectively. All these three regression models used the generalized MR constitutive equation shown in Eq. 4.17. Yau and Quintus (137) proposed seven different sets of regression equations for crushed stone materials, crushed gravel, uncrushed gravel, sand, coarse grained soil-aggregate mixture, fine grained soil-aggregate mixture, and fine grained soil. The physical properties used to correlate with the generalized MR model coefficients were material percent passing No. 3/8 in. sieve, No. 4 sieve, No. 40 sieve, No. 200 sieve, percentage of silt (% silt), percentage of clay (% clay), LL, PI, OMC, MDD, TMC, and TDD. Malla and Joshi (139) studied the correlation between the MR model coefficients and the gradation and compaction parameters of base materials. Regression analysis was performed on MR model coefficients for uncrushed gravel and crushed limestone. The variables were selected based on the literature and the logical influence of gradation, compaction, and moisture content parameters on the regression constants (k1, k2, and k3). These variables included OMC, TMC,

54 MDD, material percent passing No. 4 sieve, material percent passing No. 40 sieve, and material percent passing No. 200 sieve. Soliman and Shalaby (148) presented regression equations of the generalized MR model coefficients for six AASHTO soil types (A-1-b, A-3, A-2-4, A-4, A-6, and A-7-6). Data used in the regression analysis were collected from 19 states in New England and the nearby regions in the United States and two provinces in Canada. The soil properties that were considered in the regression analysis include: OMC, TMC, MDD, LL, PI, TDD, percentage of coarse sand (% CSAND), percentage of fine sand (% FSAND), percentage of silt (% silt), percentage of clay (% clay), material percent passing 3 in. sieve, 2 in. sieve, 1 Â½ in. sieve, 1 in. sieve, Â¾ in. sieve, Â½ in. sieve, 3/8 in. sieve, No. 4 sieve, No. 10 sieve, No. 40 sieve, No. 80 sieve, and No. 200 sieve. The accuracy of the prediction of the MR model coefficients by the existing correlation models were analyzed in this study. Table 15 listed all the predicted MR coefficients using Yau, Soliman, Malla, and ANN models. Table 15. Prediction Accuracy of SWCC Fitting Parameter Models. Model R 2 RMSE k1 k2 k3 k1 k2 k3 Yau (137) â4.51 â2.15 â3.06 0.31 0.25 0.54 Soliman (148) â2.27 â0.16 â1.42 0.43 0.59 0.89 Malla (139) â0.023 â0.03 â1.86 28.82 0.94 0.61 ANN 0.76 0.72 0.63 0.10 0.05 0.04 Compared to the ANN model, the estimated RMSE values from the three regression models were high, and the R2 values were extremely low. High RMSE and low R2 values indicated a poor prediction accuracy. The R2 values for k1, k2, and k3 from all the three regression models were less than zero, while R2 values from ANN model were in the range of 0.63 to 0.76. The negative R2-values indicated the slope of the best fit line relating the predicted to the measured coefficients was negative. Therefore, the developed ANN models had a significant improvement in predicting the MR model coefficients from base physical and compaction properties. The accuracy of MR values using the predicted coefficients were examined in the next section. Predicted coefficients from the three regression models were input into Eq. 4.17, and the ANN model coefficients were input into Eq. 4.20 to calculate MR values, which were compared to the measured values. Figure 23âFigure 24 showed the comparisons between measured and estimated MR values using the ANN predicted coefficients for base and subgrade materials.

55 (a) Plastic base materials (b) Non-plastic base materials Figure 23. Comparison of ANN Model Predicted Resilient Moduli against Measured Values for Base Materials. (a) Plastic subgrade materials (b) Non-plastic subgrade materials Figure 24. Comparison of ANN Model Predicted Resilient Moduli against Measured Values for Subgrade Materials. Clearly the MR values estimated using the ANN coefficients matched well with the test data. The R2 value estimated for both plastic and non-plastic base materials were 0.91, which demonstrate a high accuracy. Figure 25a, Figure 25b, and Figure 25c presented the results from Yau, Soliman, and Malla models, respectively. The R2 value estimated from Yau, Soliman, and Malla models were 0.32, â1.53, and â1.62e06, respectively, which indicate that the regression models failed to predict the coefficients. The negative R2 value from Soilman and Malla model indicate that the best fit line relating the predicted to the measured values had a negative slope.

56 (a) Yau model (b) Soliman model (c) Malla model Figure 25. Comparison of Measured versus Predicted Resilient Moduli Using Regression Models. Validation of the Developed ANN Models To validate the developed ANN models for the MR model coefficients, test data were collected from the Texas A&M Transportation Institute (TTI) (74, 82) and University of Illinois Urbana-Champaign (UIUC) (149). Epps et al. (74) conducted triaxial repeated load tests to evaluate the laboratory performance of unbound granular materials. In addition to the MR test results at different combinations of confining pressure and deviator stresses, aggregate physical, strength, and moisture properties were also collected. Titus-Glover and Fernando (82) performed the compressive creep and recovery tests for collected unbound base materials in the TTI laboratory. Tutumluer et al. conducted MR tests at the 15 various stress states. Since no suction values were provided by UIUC researchers, the suction is estimated using the developed ANN models from the base physical properties. Table 16 listed all the collected physical properties used as the ANN inputs for the validation data set.

57 Table 16. Input Parameters Collected from Literature for Model Validation. Data Sources Soil source #3/8 sieve #200 sieve PL PI MDD OMC TMC Epps et al. (74) E-02-2-3-2 45.9 12.36 10 5 136.2 7.2 7.2 E-08-2-1-6 49.2 6.3 9 5 140.4 6.5 6.5 Titus-Glover and Fernando (82) Limestone (â2% OMC) 64.71 17.61 12.8 7.9 127.4 9.045 7.045 Limestone (OMC) 64.71 17.61 12.8 7.9 127.4 9.045 9.045 Tutumluer et al. (149) Dolomite FC 4% (OMC) 69.9 4 N/A N/A 133 10.5 10.5 Dolomite FC 8% (OMC) 71.1 8 N/A N/A 133 8.48 8.48 The collected base properties were input into the trained ANN model to predict k1, k2, and k3. MR values were calculated using the predicted coefficients from the ANN models and compared with the test results. Figure 26 showed the comparison between measured and estimated MR values using the ANN predicted coefficients. The MR values estimated using the predicted coefficients from the ANN models produced a good fit with the test results. Figure 26 depicted an R2 value of 0.8 between measured and predicted data set. Hence, the developed ANN model can accurately estimate the MR model coefficients that are used in pavement design and analysis. Figure 26. Validation of Measured versus ANN Predicted MR at Various Stress Levels for Collected Unbound Materials. MODULUS OF SUBGRADE REACTION FOR RIGID PAVEMENTS The subgrade layer acts as a foundation in a rigid pavement structure. Proper characterization of subgrade strength is necessary for pavement design and construction. The modulus of subgrade reaction (k-value) is widely used to evaluate subgrade strength and soil structure interaction, and design the rigid pavements. It is represented as the reaction pressure (P) sustained by the soil under a rigid plate per unit settlement (â).

58 The k-value can be measured either from a field plate load test conducted on top of the subgrade (150, 151) or correlation with other load bearing capacity tests (e.g., consolidation test, triaxial test, and California Bearing Ratio test) (152). The k-value from the plate load test is calculated using the Winkler model, which assumes an elastic plate resting on a liquid foundation (88). This model considers the soil behavior as a series of linear elastic springs, as shown in Figure 27a. Many researchers have used the model to characterize soil-foundation interaction (153â157). In the Winkler model, the submodel for the elastic plate is well established, but the submodel of the foundation is particularly in need of modification. Researchers attempted to improve the Winkler model by adding other forms of interaction among the spring elements (89, 91, 158â160). The Pasternak model shown in Figure 27b allows the transverse connection in the supporting media and introduces shear interaction on slab-base interface (161, 162). (a) (b) Figure 27. Foundation Models for Rigid Pavement (a) Winkler Model; (b) Pasternak Model. The problem of beams on foundations was studied using the FE method and compared the k-values with an analytical model (163, 164). In the FE analysis, they considered the complex phenomena such as the realistic contact stresses between slab and foundation, and the nonlinear stress-dependent elastoplastic behavior of the soil foundation. Recently, the concept of backcalculating the modulus values was proposed using the deflection basin generated by a FWD (165). Hall (166) derived the method of determining of the dynamic k-value of the subgrade from the FWD deflection basin. Subsequently, the backcalculated best-fit (BBF) approach was adopted by the LTPP program to determine k-values for rigid pavements (167). For the rigid pavement analysis in the Pavement ME Design, the load bearing capacity of the subgrade foundation is characterized by an effective dynamic k-value. The actual rigid pavement structure is transformed into an equivalent structure with a k-value, which represents the compressibility of all layers beneath the concrete slab. The subbase and subgrade materials in the design guide have been characterized by the Winkler model, which is easy for implementation but neglects the shear interaction with the supporting media. In the BBF approach, the degree of bonding on the slab-base interface is assumed to be either 1 for completely bonded or 0 for non-bonded situation. However, none of these assumptions are realistic for base courses (168). Shear restraints exist between the slab and base layer because of rough interfaces. This causes an interlocking of the slab and base, which has a significant effect on the structural performance of pavement (24, 25, 169). In addition, no cross anisotropy is considered in the Pavement ME Design while calculating the k-value. Recent studies have shown that granular base materials exhibit stress- dependent and cross-anisotropic behavior, which indicates the directional dependency

59 (anisotropy) of the base moduli. Newly acquired triaxial testing devices were introduced by Abu- Osei et al. and Tutumluer and Seyhan (13, 170) for determining anisotropic resilient properties of granular materials in the laboratory. Tutumluer et al. and Park and Lytton (171, 172) found that the use of nonlinear anisotropic base modulus significantly affects the stress distribution in the base layer and diminishes the horizontal tensile stresses in the bottom half of the base layer. The performance of pavement was predicted by Masad and Masad and Oh et al. (5, 39) using the nonlinear cross-anisotropic model and showed good agreement with the field measurements. Hence, there was a need to develop a modified k-value model, which involves both corrected base modulus due to cross anisotropy and the slab-base interface bond. Development of Modified k-value Model The development of modified k-value model is elaborated as below, which contains four submodels: ï· Cross-anisotropic modulus submodel for the base layer. ï· Slab-base equivalent thickness submodel. ï· Slab-base interface shear bonding submodel. ï· Modified k-value submodel. Cross-Anisotropic Modulus Submodel for Base Layer The cross-anisotropic behavior of base material is characterized by the generalized Hookeâs Law, which is shown in Eq. 4.23: 1 1 xy xx x x y x x y yxy xy x x y x E E E E E E ï® ï® ï³ ï¥ ï³ ï¥ï® ï® ï³ ïï© ï¹ï ï© ï¹ïª ïº ï© ï¹ïª ïºïª ïº ï½ ïª ïºïª ïºïª ïºï ï ï« ï»ïª ïºïª ïº ï« ï» ïª ïºï« ï» (4.23) where Ex is the horizontal modulus. Ey is the vertical modulus. Î½xy is the Poissonâs ratio to characterize the effect of vertical strain on horizontal strain. Î½xx is the Poissonâs ratio to characterize the effect of horizontal strain on the perpendicular horizontal strain. An FE program, WinJulea, is used in this submodel to determine the stress values at the mid-depth of the base layer. The estimated stress values are input in the generalized MR model, shown in Eq. 4.24, to calculate the vertical modulus in the base layer (15): 32 1 ( ) ( 1) kkv oct R a a a M k P P P ï´ï± ï½ ï« (4.24)

60 H R V R Mn M ï½ (4.25) where V RM is the resilient modulus in the vertical direction. H RM is the resilient modulus in the horizontal direction. n is the modulus ratio. The resilient moduli in the vertical and horizontal direction follows a constant ratio, n. The corrected base modulus is the golden mean of the vertical and horizontal modulus, as shown in Eq. 4.26: 1 / 2 1 / 2( ) ( )V H VR R R RM M M M nï½ ï´ ï½ (4.26) An iterative model is developed within this submodel to determine the cross-anisotropic MR of the base material. Figure 28 is a schematic flowchart for the iterative method. The input parameters include the pavement layer thicknesses and initial moduli in the FE program. Initial values of the vertical and horizontal stresses are computed from the FE model analysis. The computed initial stresses are then used to calculate the output (i.e., calculated modulus) from Eqs. 4.24, 4.25, and 4.26. Eq. 4.27 presents the convergence criterion.

61 Figure 28. Flowchart of Corrected Base Modulus due to Cross Anisotropy. ( ) ( ) ( ) 1%R Calculated R Input R Calculated M M M ï ï£ (4.27) This indicates that the base modulus is finalized whenever the difference between the input and calculated modulus is less than 1 percent of the calculated modulus. If the difference between calculated and input MR > 1 percent of ( )R C a lc u la te dM , the calculated modulus is input again in the FE program, and the iteration process is continued until the desired criterion is reached. Slab-Base Equivalent Thickness Submodel The formulation of the equivalent thickness submodel includes three steps. First, a cooperating slab and base system is introduced using the transformed-section method (Figure 29).

62 Figure 29. Illustration of Transformed-Section Method for a Cooperated Concrete Slab and Base Course System. Second, the moment of inertia of the pavement cross section is calculated based on the transformed section as follows: 2 t r s la b b a s e i iI I I A dï¤ï½ ï« ï« ï¥ (4.28) where Itr is the moment of inertia of the transformed pavement section. s la bI is the moment of inertia of the slab. ba seI is the moment of inertia of the transformed base. Ai is the area of the slab and the transformed area of the base course. Ídi is centroidal distance to each of the areas. Î´ is the interface shear bonding. Table 17 shows the calculation steps of the moment of inertia for a typical cooperated slab and base system. Table 17. Steps of Moment of Inertia Calculation for a Cooperated Slab and Base System. BodyÂ Area(A)Â y Â Ay Â z Â . .C GI Â d Â . .N AI Â (1)Â sbh Â / 2sh Â 2 / 2sb h Â i i i A y A ï¥ ï¥ Â 31 12 s bh Â 2 shz ï Â 1 2 . . 1 1C GI A dï« Â (2)Â ( )b b s Eb h E Â / 2s bh hï« Â ( / 2) ( )bs b b s Eh h b h E ï« Â 3 1 ( ) 12 b b s Eb h E Â 2 b s hh zï« ï Â 2 2 . . 2 2C GI A dï« Â Â 1 2iA A Aï½ ï«ï¥ Â Â 1 1 2 2i iA y A y A yï½ ï«ï¥ Â Â Â Â Â Third, the thickness of the equivalent section is calculated using the moment of inertia of the transformed section. The transformed section that consists of concrete and base is converted into an equivalent cross section composed of only concrete. The thickness of the equivalent section is estimated by considering the same moment of inertia for both sections:

63 3 *12tr eq Ih b ï½ (4.29) Clearly, the thickness of the equivalent cross section depends on the degree of bonding, ï¤ , in the slab-base interface. A formulation of slab-base interface shear bonding is presented in the next section. Slab-Base Interface Shear Bonding Model The interface shear bonding, ï¤ , in Eq. 4.30 is expressed by the ratio of in situ shear stress, 2( ) fz x ï±ï´ , in the base course on the slab-base interface, and the shear stress, 2v , in the base layer on the interface when full shear is transferred: 2 2 2 ( ) , ( ) f f zx zx fsv ï± ï± ï´ ï¤ ï´ï½ ï£ (4.30) where fs is the shear strength of the base course. Figure 30 shows a schematic figure of the developed shear stress in the PCC-base interface. The parameters 1v and 2v are shear stress on the interface in the slab and base layer, respectively, when full shear is transferred. The parameter 2( ) fz x ï±ï´ is the shear stress in the base course on the interface for in situ conditions, which is limited by the shear strength, fs , of the base course on the failure plane. When 2( ) fzx ï±ï´ is greater than 2v , the interface is considered as fully bonded. Depending on the ratio of the 2( ) fzx ï±ï´ and 2v , the partial bonding condition is defined in the PCC-base interface. Figure 30. Illustration of In Situ Shear Stress in the Base Course on the PCC-Base Interface Using a Mohr-Coulomb Failure Envelope. The expression of interface shear bonding, ï¤ , is given as follows; the detailed derivations are presented in Appendix E:

64 ï¨ ï© 2 2 2 2 1 3 2 2 2 2 1/22 2 2 2 2 5/2 3 ( ) ( tan ) cos sin[tan ( )] 3( )1 (1 2 )[ ] 2 ( ) ( )3 2 2 ( ) ( ) 2 s s s s s s s b b s s ss s ah h a c h h a h a h a h hh h zahP hh a h z ï³ ï¦ ï¦ ïµ ï¤ ï° ï ï© ï¹ ïª ïºï«ï« ï»ï« ï ï ï« ï« ï© ï¹ï« ï«ïª ïºï« ï»ï½ ï« ïï© ï¹ ïª ïºï«ï« ï» ï (4.31) Modified Subgrade k-value Submodel The modified k-value submodel is developed using the calculated heq and Î´ values from the previous submodels. The determination of the modified k-value is divided into the following steps: ï· The deflection patterns generated by FWD are used to determine the modified k-values. The FWD sensor deflections (0 cm, 30.48 cm [12 inch], 60.96 cm [24 inch], 91.44 cm [36 inch], and 121.92 cm [48 inch] away from the loading point) are obtained, and the basin area (BA) is calculated as: 0 1 2 1 0 [ 2( ....... ) ] 2 j j SSBA D D D D D D ï ï½ ï« ï« ï« ï« ï« ï´ (4.32) where SS is the FWD sensor spacing (30.48 cm); jD is the sensor deflection (j=0 to 4). ï· The effective relative stiffness length is calculated as follows: 2( ) ( )el a b B A c B Aï½ ï« ï´ ï« ï´ (4.33) where a, b, and c are the coefficients obtained from the field correlation, a=0.992, b=â0.2891, c=0.0284 (180). ï· Finally, the subgrade k-value is formulated in Eq. 4.34. The concept of equivalent slab thickness is applied from Eq. 4.29 to consider the effects of both the slab and base layer on the pavement subgrade k-value: 3 2 412(1 ) s eq e E h k lï® ï½ ï (4.34) where sE is the elastic modulus of the PCC; ï® is the Poissonâs ratio of PCC. Estimation of Modified k-value for LTP Pavement Sections The pavement data used in this study to estimate the modified k-value are collected from the LTPP database. A total of 505 sets of rigid pavement section data are collected, including both JPCP and CRCP. In addition to the backcalculated layer moduli, the database includes layer thicknesses and FWD deflection patterns. In the LTPP, the FWD test uses three loading

65 sequences with a target load of 40 kN (9000 lb), 53.38 kN (12000 lb), and 71.17 kN (16000 lb). In this study, only the deflection basins generated from 53.83 kN load are selected for analysis. To avoid the discrepancy in measured deflection patterns, the deflection basin tests that were performed along the mid-lane path are selected for the analysis. The mid-lane test locations are designated as J1 and C1 for JPCP and CRCP pavements, respectively. Correction of Base Modulus As shown in Eq. 4.25, the resilient moduli of base material in the horizontal and vertical directions follow a constant ratio. Based on the existing laboratory results, the n value is in the range of 0.3â0.5. In this study, the n value is assumed as 0.4 for all base materials collected from the LTPP pavement sections. Following the flow chart in Figure 28, the modulus of the base layer is corrected for all pavement sections. The procedure to determine the corrected base modulus is illustrated for an example LTPP section 27-4054 (State code-SHRP ID). The pavement section 27-4054 consists of a 0.24 m (9.4-in.) PCC layer, a 0.15 m (6-in.) unbound base layer, and a semi-infinite silty clay subgrade layer. Herein, the backcalculated modulus value is collected from the LTPP database and used as the initial input at each layer. Table 18 lists the value of collected MR coefficients k1, k2, and k3 for the base layer and the determined vertical and horizontal stress values after first iteration. Figure 31 shows the plot of mean base modulus against iteration number. As seen, final value is achieved after four iterations when it meets the convergence criteria shown in Eq. 4.27. Table 18. Collected MR Coefficients and the Simulated Stress Values at the Mid-Depth of Base Layer for LTPP Section 27-4054. k1 (kPa) k2 k3 Ï1 (kpa) Ï3 (kPa) 998.5 0.55 â0.023 16.75 9.38 Figure 31. Base Resilient Modulus Convergence with Iteration Number. Estimation of Slab-Base Interface Bonding Ratio The interface shear bonding ratio is calculated for the collected LTPP pavement sections. From the collected 505 pavement sections, 267 sections have a treated base layer beneath the PCC surface layer and 238 pavements section are constructed using unbound base course.

66 (a) (b) Figure 32. Formulation of Friction Angle from Mohr-Coulomb Failure Envelope for (a) Treated Base; (b) Unbound Base. As shown in Eq. 4.31, the expression for the degree of bonding, ï¤ , is based on the shear strength parameters c and ï¦ . In the LTPP database, the available experimental data are not sufficient to determine the shear strength parameters. Unconfined compressive strength test data are available for cement and lime treated base materials, whereas triaxial shear strength tests are conducted for unbound base materials only at 34.47 kPa (5 psi) confining pressure. Hence, empirical models are used in this study to estimate the c and ï¦ parameters from base strength and saturation properties. Stabilized Base Layer Existing studies have investigated the relationship between the cohesion parameter, c, and the physical and strength properties of stabilized base materials. Balmer and Clough et al. (176, 177) reported that cement treatment increased cohesion while internal friction angle remained constant. Thompson (178) stated that addition of lime to unbound material yields a substantial increase in cohesion and minor improvement in internal friction angle. This study uses the relationship suggested by Park and Lytton (172) to determine the cohesion of treated soils based on unconfined compressive strength: 9 .3 0 .2 9 2 cc ï³ï½ ï« (4.35) where cï³ is the unconfined compressive strength, unit: psi. Then the friction angle, ï¦ , is determined from c and the collected unconfined compressive strength test data (Figure 32a). Unconfined compressive strength test data are collected for 55 treated base materials from the LTPP database. Degree of bonding, ï¤ , on the slab-base interfaces is calculated from Eq. 4.31 using the estimated shear strength parameters.

67 Unbound Base Layer Similarly, many researchers had analyzed the shear strength test data and identified the influence of several physical properties on the shear strength parameters of unbound aggregate materials. Density and degree of saturation are the two most important parameters affecting the shear strength parameters of unbound base (179). The effect of relative density and degree of saturation on the shear strength parameters was examined by Theyse (180) and proposed the following relationship for cohesion: 12.12 2.30c 0.0107 RD Se eïï½ (4.36) where RD is the relative density, unit: %. S is the degree of saturation, unit: %. c is the cohesion, unit: kPa. Density and water content data are collected from the LTPP database and used to calculate the RD and S in Eq. 4.36. Figure 32b presents the equation to estimate the friction angle ï¦ . Shear strength and compressive strength test data are collected from 88 granular unbound base courses, which are further used to calculate ï¤ in Eq. 4.31. The black and green bars in Figure 33a denote the degree of bonding on slab-base interface using the proposed shear bonding submodel and BBF approach (181), respectively. Using the proposed bonding submodel, approximately 50 percent of the pavement sections with treated base material is estimated to be fully bonded (ï¤ =1) with the PCC slab layer and rest of the pavements have partially bonded slab-base interface. However, no partial bonding condition is adopted in the BBF approach. The bonding condition in the BBF approach is assumed based on the PCC modulus (167). The slab-base interface is considered as fully bonded if the PCC modulus is greater than 26890 MPa (3900 ksi) whereas a non-bonded interface is applied for lower PCC modulus values. Figure 33b compares the estimated degree of bonding using the proposed bonding submodel with the BBF approach for unbound base course. Most pavement sections with unbound granular base layer are estimated as partially bonded with the PCC layer whereas a full bonding condition is assigned in the BBF approach. (a)

68 (b) Figure 33. Comparison of Calculated Slab-Base Interface Degree of Bonding Ratio with the BBF Approach for (a) Treated Base; and (b) Unbound Base Layer. Comparison with Wheelpath Faulting The degree of interface bonding has a significant impact on erodibility prediction of rigid pavements (168, 169, 182, 183). In the process of erosion, the concrete slab deforms under repeated traffic loading, which causes mechanical shear stress between the slab and base. As a result, adequate interface bonding can significantly reduce erosion. Erosion is addressed through the faulting distress model (66). A better interfacial bonding model can infer better sensitivity to measured faulting. Figure 34 presents the correlation between the measured wheelpath faulting and the calculated degree of bonding. The wheelpath faulting data are collected for JPCP pavement sections from the LTPP database. Figure 34. Sensitivity of Slab-Base Degree of Bonding on Wheelpath Fault (mm).

69 The calculated degree of bonding shows better sensitivity to the observed faulting data. The faulting value decreases with a higher degree of bonding in the slab-base interface. However, the assigned interface bonding in the BBF approach shows no sensitivity to observed faulting. The developed interface bonding model can better predict the interface condition for rigid pavements. The modified k-values are estimated in the next section using the calculated interface degree of bonding. Estimation of Modified k-value The modified k-values are calculated for the collected LTPP pavement sections following the steps in modified k-value model that were outlined previously. Figure 35 plots the modified k-values versus BBF k-values for the same pavement structure collected from the LTPP. As presented, the BBF approach has higher k-values than the modified k-value model. This discrepancy has occurred due to the different interface bonding ratios. The interface between the slab and base is considered as fully bonded in the BBF approach for most of the pavement sections whereas most of them are partially bonded in the modified k-value model. Figure 35. Comparison of Modified versus LTPP k-values. An ANN model is developed in the next section to determine the modified k-values for various combinations of pavement structure and layer modulus. A wide range of pavement structural properties, which include layer thicknesses and material strength properties, resilient moduli, and slab-base interface bonding ratios, are input in the ANN model. The development of the ANN model is elaborated as follows. Development of ANN Model The ANN approach is an adaptive information processing technique, which establishes the correlations between the input variables Xi and the output variables Yj through the inter- connected neurons (i.e., weight factor, wji). The input variables Xi and the output variables Yj are usually normalized to xi and yj, respectively. The correlations developed by the ANN models between the normalized input parameters xi and the normalized output variables yj are shown in Eq. 4.37:

70 1 j n j ji i i by f w x ï½ ï¦ ï¶ ï«ï§ ï· ï¨ ï¸ ï½ ï¥ (4.37) where f is a transfer function, which normally uses a sigmoidal, Gaussian, or threshold functional form. wji and bj are the unknown weight factors and bias term, respectively. The ANN model specifically adjusts the weight factors wji and bias bj in Eq. 4.37 based on the minimum error function. In pavement engineering, the ANN approach is usually used to develop prediction models based on many data collected from experiments and numerical analysis. The ANN models have been successfully developed to predict the crack growth function (e.g., reflective cracking and top-down cracking) in asphalt concrete (144, 145). The ANN approach was also used to predict the geogrid-reinforced flexible pavement performance (118). In general, the development of ANN models includes two critical steps: 1) deflection basin estimation; and 2) construction of ANN architecture. These two steps are described below. FE Model Development The deflection basin of a specific pavement structure for FWD loading is estimated using FE program ABAQUS. Researchers developed a total of 1296 simulation cases with different combinations of pavement layer thicknesses, layer modulus, and PCC-base interface bond. Pavement structural responses are calculated under FWD loading using the FE software ABAQUS. Figure 36 shows a typical rigid pavement structure and the corresponding FE model. In this study, a large variety of pavement structures are modeled with different combinations of PCC, base thickness, base, and subgrade modulus. To properly represent the effect of the strength of the base course and its level of erosion on the interface, it is necessary to accurately characterize the slab-base interface. A Coloumb interface model is placed between the bottom of the concrete and the top of the base course. The interface model has a normal and a tangential load transfer behavior at each nodal point. By adjusting the friction coefficient and the elastic slip distance, it is possible to represent various shear strength levels of the base course on the interface, including zero strength. Varying the strength of the base course interface will generate different deflection patterns. These deflection patterns are used to generate the different levels of the subgrade k-values that are needed to represent the effect of lowered shear strength and erosion of the surface of the base course on the foundational support. The FWD sensor deflections (0, 30.48 cm, 60.96 cm, 91.44 cm, and 121.92 cm away from the loading point) are obtained from the FEM analysis, and the modified k-value is calculated using Eqs. 4.32, 4.33, and 4.34.

71 (a) (b) Figure 36. (a) Schematic Plot of a Typical Pavement Structure; (b) Axisymmetric Model of Pavement in ABAQUS. Construction of ANN Architecture Researchers constructed a three-layered ANN model in this study (Figure 37). Pavement structural properties (i.e., PCC slab and base thicknesses), strength properties (i.e., PCC, base, and subgrade moduli), and PCC-base interface bonding ratio were introduced as input parameters for the ANN model. Table 19 summarizes the range of the input parameters in the ANN model. The hidden layer is assigned 20 neurons to establish the connection between the output layer and the input layer. The number of neurons is selected based on the number of data points (146). The developed ANN model uses the sigmoidal transfer function, as shown in Eq. 4.38 (118): ï¨ ï© ï¨ ï© 1 1 expi i f I Iïª ï½ ï« ï (4.38) where Ii is the input quantity. Ï is a positive scaling constant. The parameter Ï controls the steepness between the two asymptotic values 0 and 1. The ANN model determines these weight factors wji through the two major functions: training and validating. The training data set is used to determine the trial weight factors, wji, and bias term, bj,. and the validating data set is employed to examine the accuracy of the model prediction.

72 Table 19. Selected Range of Input Parameters in ANN Training Data Set. Input parameters Level Input values PCC thickness (mm) 3 178, 254, and 348 Base thickness (mm) 3 101.6, 203.2, and 254 PCC modulus (MPa) 3 14,420, 41,400, and 82,737 Base modulus (MPa) 4 69, 690, 6894, and 25,000 Subgrade modulus (MPa) 3 34.5, 282, and 551 PCC-base interface bonding 4 0, 0.3, 0.6, and 1 Figure 37. Illustration of Three-Layered Neural Network Architecture for k-values. In this study, 80 percent of the data set is used for training and 20 percent of the data set is for validation. The training algorithm uses the Levenberg-Marquardt back propagation method to minimize the MSE. The gradient descent weight function is employed as a learning algorithm to adjust the weight factors wji. Figure 38 shows the comparison between output and target k-values from the ANN model. The output k-values shows high R2 for the training and the validation data set, which indicates a good prediction accuracy. To validate the developed ANN model, the ANN predicted k-values are compared with the values calculated with the modified k-value model. Figure 39 presents the comparison of calculated and ANN predicted k-values for LTPP pavement sections at the calculated degree of bonding. The ANN model shows high prediction accuracy with an R2 value of 0.92. This indicates that the developed ANN model is capable of accurately predicting the modified k-value at any given degree of bonding. A sensitivity analysis for the degree of bonding is performed for selected pavement sections and shown in Figure 40. Figure 40 illustrates the modified k-values increase with the increase of degree of bonding for most pavement sections. However, for pavement sections 8-0214 and 29-500, the k-values have a decreasing trend with the increase of bonding ratio. The reason behind that is the increase of slab-base bonding ratio causes a higher thickness of equivalent transformed section and basin area. The k-value increases with the increase of equivalent thickness, but it continues to decrease when the basin area increases. More sensitivity analyses are presented in Appendix F.

73 Figure 38. Target and Output k-values for Training, Validation, and Overall Data Sets for 1296 Simulation Cases. Figure 39. Comparison of Calculated versus Predicted Modified k-values.

74 Figure 40. Modified k-values at 0, 0.3, 0.6, and 1 Degree of Bonding for Selected LTPP Pavement Sections. SHEAR STRENGTH OF BASE AND SUBGRADE FOR FLEXIBLE AND RIGID PAVEMENTS The Pavement ME Design mainly considers the elastic behavior of unbound base and subgrade materials, but little attention has been paid to their shear strength. Previous studies have presented that the shear strength of underlying layer materials is related to the pavement performance through the following aspects: ï· Influence on accumulation of permanent deformation in flexible pavements. ï· Influence on extent of erosion in rigid pavements. ï· Influence on degree of bonding between the concrete slab and the base course in rigid pavements. In flexible pavements, the shear strength directly affects the amount of total rutting (179, 184, 185). For example, the dominant factor in determining permanent deformation is the relationship between the shear strength of the soil and the applied shear stress. The ratio of the shear strength to shear stress performs well in limiting the permanent deformation of granular materials against shear failure (179, 184). Shear strength affects the extent of erosion in the base course of rigid pavements, which is the key factor in the development of faulting in JPCP and punchouts in CRCP. Soil Erodibility (186) depicted that a small change in the shear strength has a considerable effect of the degree of erosion in the lower range of shear strength. Furthermore, in rigid pavements, the shear strength is directly related to the degree of bonding between the concrete slab and base course and the degree of bonding has a significant impact on the performance of rigid pavements (168, 169, 182, 183). In light of the critical role of shear strength in performance prediction, researchers incorporated in the model formulation of both flexible and rigid pavements. The general shear strength model is defined according to the Mohr-Coulomb failure envelope, which is determined from triaxial tests on laboratory molded specimens:

75 ta nncï´ ï³ ï¦ï½ ï« (4.39) where ï´ is the shear stress. c is the total cohesion. nï³ is the normal stress on the failure plane. ï¦ is the angle of internal friction. Furthermore, the impact of moisture variations on the shear strength is considered in this study. As the water content increases by a small amount, the shear strength decreases significantly depending on the magnitude of normal stress (187). Such a reduction accelerates shear failure and intensifies rutting in flexible pavements and erosion in rigid pavements. Researchers have developed a moisture-sensitive shear strength model for unbound base materials (60): ' ( ) tan 'n mc fhï´ ï³ ï± ï¦ï½ ï« ï (4.40) where 'c is the effective cohesion. 'ï¦ is the effective friction angle. ï± is the volumetric water content. f is the saturation factor. hm is the matric suction. Figure 41 illustrates Eq. 4.40, which represents the dependence of the shear strength of unbound materials on the matric suction. To make the moisture-sensitive shear strength model more applicable in the pavement design, researchers developed prediction models for the shear strength parameters c and Ïâ. In this way, the shear strength of unbound layers and subgrade can be estimated using common soil properties in the absence of triaxial test data. Figure 41. Schematic Plot of Mohrâs Circle Showing Dependence of Shear Strength on Matric Suction. Prediction Model for Unbound Base To develop the prediction models for unbound base materials, researchers collected the measured shear strength parameters of different types of base materials from Epps et al.,

76 Tutumluer et al., and Chow et al. (74, 189, 190). Additional data were collected to obtain soil physical properties such as gradation, Atterberg limit, optimum and saturated moisture content, and the SWCC of the various materials tested. Finally, regression analysis was conducted to determine the relationships between the shear strength parameters (câ and Ïâ) and the collected physical properties: ' 0.221* 4.6 * (%) 455.62 * 1262.75 ' 0.0272 * 0.638* (%) 1.487 * (%) 69.92 sat s opt sat c PI G PI MC ï± ï¦ ï± ï½ ï ï« ï« ï ï½ ï ï ï« (4.41) The goodness of fitting by the developed regression equations is shown in Figure 42 for the collected base materials. R2 value of 0.81 and 0.87 for câ and Ïâ parameters, respectively, indicate that the regression equations have a good prediction accuracy. Â (a) (b) Figure 42. Comparison of Predicted and Measured Shear Strength Model Parameters (a) câ and (b) Ïâ. Prediction Model for Subgrade Similar prediction models have been developed for the shear strength parameters of subgrade. Due to insufficient triaxial test data for subgrade, unconfined compressive strength test data were collected from LTPP database. The friction angle had been found to be solely dependent on the PI by Holtz and Kovacs (188). They developed an empirical correlation between Ïâ and PI for normally consolidated soils as shown in Eq. 4.42: 2' 0 .0014 * 0.28 * 35 .87P I P Iï¦ ï½ ï ï« (4.42) The matric suction, friction angle, porosity, and gradation of the subgrade material affect the cohesive shear strength, and this is confirmed with subsequent studies in the mechanics of unsaturated soils (e.g., 60, 74). In this study, an ANN model has been developed to predict the cohesion parameter from subgrade physical and strength properties. Figure 43 shows the architecture of the ANN model. MDD, soil porosity, effective friction angle, and suction value at OMC are selected as input parameters. The hidden layer is consisted of 10 neurons and the câ 0 50 100 150 Measured cohesion,c (kPa) 0 50 100 150 R2 = 0.81 20 25 30 35 40 45 50 55 60 65 70 Measured friction angle, (degree) 20 25 30 35 40 45 50 55 60 65 70 R2 = 0.87

77 parameter is the output of the developed model. The cohesive strength parameters were calculated from the collected 432 unconfined compressive strength test samples. Figure 43. Illustration of Three-Layer Neural Network Architecture to Predict câ Parameter. Figure 44 shows the prediction accuracy of the developed ANN model for training and validation data sets. The obtained R2 value 0.98 and 0.97 for training and validation data set demonstrate that the developed ANN model can accurately predict the câ parameter.

78 Figure 44. Target and Output câ Values for Training, Validation, and Overall Data Sets for 432 Subgrade Soils. PERMANENT DEFORMATION OF BASE AND SUBGRADE FOR FLEXIBLE AND RIGID PAVEMENTS Mechanistic-Empirical Permanent Deformation Model for Unbound Materials Permanent deformation is a major form of distress in flexible and rigid pavement. In the Pavement ME Design, the total permanent deformation is the sum of the individual layer permanent deformations (i.e., surface layer, unbound layers, and subgrade layer). The permanent deformation model for unbound base and subgrade layers used in the Pavement ME Design is shown in Eq. 4.43: ( )0( ) Ns v r e h ï¢ï²ï¥ï¤ ï¢ ï¥ ï¥ ï ï½ (4.43) where Î´ is the permanent deformation for the layer. Ïµr is the resilient strain imposed in the laboratory test. Ïµv is the average vertical resilient strain in the layer.

79 h is the thickness of the layer. N is the number of traffic repetitions. Ïµ0, Ï, and Î² are the model coefficients. Î²s is the global calibration coefficient (1.673 for granular materials and 1.35 for subgrade soils). However, the rutting model used in the Pavement ME Design is less sensitive to the modulus and the thickness of the unbound layers (4). As a result, a different permanent deformation model is proposed by Tseng and Lytton (70), shown in Eq. 4.44: ( ) 0 Ne h ï¢ï² ï¤ ï¥ ï ï½ (4.44) The Tseng-Lytton model was improved, and a new ME permanent deformation model for unbound materials was proposed, which is capable of predicting the permanent deformation behavior at different stress states using the single-stage test protocol (73). The formulation of the model is given as follows: 2( / ) ( )1 0( ) ( ) ( ) ( ) N m n m n p a a a J I KN e p p p ï¢ï² ï¡ï¥ ï¥ ï ï ï«ï«ï½ (4.45) 2 sin ' 3 (3 sin ') ï¦ï¡ ï¦ ï½ ï (4.46) ' .6 cos ' 3 (3 sin ') cK ï¦ ï¦ ï½ ï (4.47) where J2 is the second invariant of the deviatoric stress tensor. I1 is the first invariant of the stress tensor. Îµ0, Ï, Î², m, and n are model coefficients. câ and Ïâ are effective cohesion and friction angle, respectively. In this model, the two terms, J2 and (Î±I1 + K), are incorporated into the Tseng-Lytton model, which is used to reflect the influence of a stress state on the permanent deformation of unbound materials. Researchers have conducted many laboratory tests and validated the accuracy of the model proposed in Eq. 4.45. They concluded that the proposed model matches well with the measured permanent deformation curves. Gu et al. (118) studied the sensitivity of the proposed model and demonstrated that it improved the sensitivity of unbound layer significantly. Hence the model proposed by Gu and Zhang is adopted is this study to predict permanent deformation of base and subgrade soils. The model coefficients (Ïµ0, Ï, Î², m, and n) in Eq. 4.45 are determined by two approaches: ï· By fitting the model to the permanent strain measured in the repeated load permanent deformation test. ï· Using prediction models developed from a set of physical properties such as the MDD, percent fines content, aggregate gradation, etc. From these two approaches, approach 1 is time consuming, requires an extensive laboratory setup, and requires expert personnel to conduct the test. In this study a regression model is developed to predict the permanent deformation of unbound materials.

80 Regression Models for Permanent Deformation Model Coefficients for Unbound Base Layers In order to develop the prediction models of the coefficients of Eq. 4.45, researchers used the measurements from the repeated load triaxial test on 108 different types of base materials collected from various sources (74, 189â192). In addition to triaxial test data, soil physical properties such as gradation, Atterberg limit, MDD, and OMC data were also collected for the base samples. Finally, regression analysis was conducted to determine the relationships between the permanent deformation model coefficients (Ïµ0, Ï, Î², m, and n) and the collected physical properties: 0 4log 0.017* 0.00967* 0.99P MDDï¥ ï½ ï ï« (4.48) 4log( ) 0.0574* 0.0937 * 5.0756P MCï² ï½ ï ï« (4.49) 40.01483* 0.00813* 0.00136* (%)P MDD satï¢ ï½ ï ï« ï (4.50) 00.2153* 0.396m ï¥ï½ ï« (4.51) 0.22993* 0.746n mï½ ï ï (4.52) where P4 is the percent of material passing No. 4 sieve (unit: %). MDD is the maximum dry density (unit: lb/ft3). MC is the test moisture content (unit: %). Sat(%) is the degree of saturation (unit: %). Figure 45 shows the goodness of fitting by the developed regression equations for the collected base materials. The predicted coefficients showed a good match with the measured values.

81 (a) (b) (c) (d) (e) Figure 45. Comparison of Predicted and Measured Permanent Deformation Model Parameters (a) 0ï¥ , (b) ï² , (c) ï¢ , (d) m, and (e) n.

82 FAULTING OF BASE LAYER FOR RIGID PAVEMENTS Illustration of Development of Faulting An understanding of the faulting progression is beneficial to develop the faulting prediction model. Moisture infiltration and non-uniform support of base course are essential keys for the development of faulting (193). Based on the key reasons, Figure 46 schematically illustrates the development of the faulting under moving wheel loads, which can be classified into four stages to describe the process of faulting. In the first stage, the upward curling of slab edge caused by temperature and moisture variations allows water to infiltrate into the base layers through joints and fills the voids between the curled slab and the base course. In the second stage, the moving wheel forces the curling approach slab downward. At the same time, water is forced out and scours the moistened fine aggregates in the base course and ejects them from the joint. This erosion forms a void beneath the approach slab. In the third stage, when the moving wheel reaches leave slab, the approach slab rebounds upward again and leaving voids beneath it and fines are ejected from beneath the leave slab. Finally, in the fourth stage, the wheel load moves away. Both slabs rebound to the curled condition. With repeated loading, the void beneath both slabs increases in size due to both the loss of the scoured fines and plastic deformation of the moistened base course. The larger dynamic load on the leave slab creates a larger void and a lower elevation. Repeated loading causes the evolution of faulting. The difference in elevation between the upward approach slab and downward leave slab results in the faulting.

83 Figure 46. Schematically Illustration of the Development of Faulting. Development of Faulting Model in Jointed Concrete Pavement An ME model was developed to estimate the growth of faulting with time. An inflection point in the field faulting data (shown in Figure 47) is found to differentiate two stages of faulting. The inflection point, as a critical faulting depth, can be directly determined by this model and is an indicator of the beginning of erosion for the concrete pavement design. The faulting before reaching the critical depth results in the accumulated deformation of the underlying layers (Figure 47). Water infiltrates through the joint and fills the void caused by the

84 permanent deformation. Beyond this inflection point, faulting accelerates because this infiltrated water being driven by the moving traffic scours the surface of the base course and causes erosion. In summary, two faulting prediction models were developed. One is to predict the entire faulting development over time. This model of the full faulting curve shows that there is an inflection point in the faulting curve. Before reaching the inflection point, the accumulation of faulting is caused by the permanent deformation of the supporting layers. After passing the inflection point, faulting accelerates due to the action of erosion. The second model is to predict the faulting depth below the inflection point with traffic. The proposed models were proven to be considerably accurate and reliable by using LTPP data. The coefficients in the models are statistically calibrated with performance-related factors using multiple regression analysis. Figure 47. Illustration of Field Faulting Data, Including the Critical Faulting Depth. Model Development for Entire Faulting Progression The ME model was developed to characterize the faulting depth over time. The mathematical form of the model, according to the faulting depth curve as shown in Figure 47, is expressed as follows: 1/ 0 [ln( )] ee Nf N N ï¢ï² ïï¥ï½ ï (4.53) where f is the faulting depth. N is the number of days after pavement construction date. 0N is the number of days when faulting initiates. Nï¥ is the number of days to failure due to erosion. eï² and eï¢ are model coefficients. Critical Faulting Depth Model A model was developed to predict the critical faulting depth in rigid pavements using the field faulting data by collecting from the LTPP database. Figure 47 indicates this concept of critical faulting depth. The critical faulting depth is defined as the inflection point in the faulting

85 depth at which erosion begins. The point of inflection is where the curvature of the faulting depth curve changes from negative to positive, as shown in Figure 46. To determine the inflection point linking those two curved segments, take the second derivative of Eq. 4.53 and set it equal to zero to yield: 1( ) 0 e e ipN N e N ï¢ ï¢ ï« ï ï¥ï½ ï« (4.54) where ipN is the number of days at the inflection point at which the acceleration of faulting begins. The critical faulting depth, ipf , which occurs at the inflection point, is expressed as: 1/1( ) eeip e e f ï¢ï¢ï² ï¢ ï« ï½ (4.55) Use of LTPP Data for Entire Faulting Progression The data for the full faulting prediction model were collected from the field data from the LTPP database, which provides a large amount of research quality pavement performance information. The data of each jointed concrete pavement test site were analyzed and compared between the measured and the model predicted faulting data. The field data collected from the LTPP database required data cleansing since they have limitations such as small data sets, poor data quality, and bias. Firstly, the minimum number of faulting points for each LTPP section are four points that are required to calculate the model coefficients ( eï² , eï¢ , 0N , and Nï¥ ). Individual data points that were clearly outliers were eliminated. Consequently, a total of around 130 LTPP pavement sections are obtained to perform the analysis of the faulting prediction model. Based on the collected eligible LTPP faulting data, the analysis with nonlinear optimization was performed in each pavement section. The four unknown coefficients, eï² , eï¢ , 0N , and Nï¥ , are determined using the SOLVER function in ExcelÂ®. Comparison Between Measured and Predicted Faulting The comparison between the measured and predicted faulting depths is made to test how well the developed model fits the data sets. Figure 48 plots the predicted faulting to observed data, which are total analysis results from the collected data sets of around 130 LTPP sections. The majority of scatter falls on the line of equality and the value of R2 in this plot is 0.95, indicating remarkable accuracy and reliability of the developed model.

86 Figure 48. Comparison between Measured and Predicted Faulting Depth. Calibration of Model for Entire Faulting Progression To implement the proposed model, model coefficients were calibrated with performance- related factors by using multiple regression analysis. Variable Selection Both categorical and numerical variables are involved in the calibration of model coefficients. The categorical variables include climatic zones (Dry-Freeze [DF], Dry Non-Freeze [DNF], Wet-Freeze [WF], and Wet Non-Freeze [WNF]), use of load transfer device (dowel), and use of drainage features. Because the support provided by the base course is the key to the progression of faulting, the numerical variables consist of aggregate gradation in the base layer and thickness of slab and base layer. The moisture infiltration also governs the development of faulting. Thus, the moisture content and saturation are also taken into consideration. Furthermore, some environmental information is incorporated into the model calibration, including the annual average number of days with high/low temperatures and with rainfall (194). Based on the collected performance-related data, a stepwise method is adopted for variable selection. This is a method of fitting regression models by adding or subtracting predictive variables at the various steps in either direction. The popular statistical software package R is used to perform the stepwise regression to select the predominant predictive variables. The p-value is the result from the studentâs t-test to indicate the level of significance of the variables. If the p-value is less than 0.05, the corresponding variable is of significance for the response variable at a 95 percent confidence level. Thus, the select predictive variables are required to meet the criterion that the p-value is less than 0.05 to ensure a 95 percent confidence level.

87 Multiple Regression Analysis The correlation between the model coefficients and the selected performance-related data is investigated using multiple linear regression analysis. The model parameter of was found to be a constant value of 2e4. Due to the variability of performance-related LTPP data, the model coefficient, , is also incorporated into the model calibration. Table 20 presents the results of multiple linear regression analysis for the coefficients in the model. From Table 20, the p-values of the predictive variables are less than 0.05, indicating that prediction models of each coefficient are significant at a 95 percent confidence level. The coefficient models are expressed as: 0 ln( ) 0.952394 0.816222 0.164161 0.000173 0.024705 0.012736 32 0.081787 0.451114 e edowel N intensedays days C basethick drainage ï¢ ï²ï½ ï ï´ ï ï´ ï« ï´ ï« ï´ ï ï´ ï ï´ ï ï´ (4.56) 2 4N eï¥ ï½ (4.57) 0ln( ) 9.003655 0.019214 32 0.009733 0.141653 0.631479 TN days C F basethick bound ï½ ï« ï´ ï ï´ ï ï´ ï ï´ (4.58) ln( ) 0.85277 1.605784 0.994204 1.275704 -1.474294 -1.379156 0.033711 32 0.149092 0.016701 e DNF WF WNF dowel bound days C basethick wetdays ï² ï½ ï ï´ ï ï´ ï ï´ ï´ ï´ ï« ï´ ï ï´ ï« ï´ (4.59) where dowel is a dummy variable for the use of dowels (No dowel = 0, dowel = 1). intensedays is the annual average number of days with precipitation greater than 12.7 mm. basethick is the base course thickness, in. drainage is dummy variable for the use of drainage features (the section subsurface drainage information can be collected from the LTPP table of SECTOIN_DRAINGE. The code of DRAINAGE_TYPE is used to denote the dummy variable drainage. drainage = 0 if code is no subsurface drainage; drainage = 1 if code is any other drainage features such as longitudinal or transverse drains, drainage blanket, etc.). TF is the freeze-thaw cycle that is number of days in the period when the air temperature goes from less than 0Â°C to greater than 0Â°C. 32days C is the annual average number of days with the maximum temperature greater than 32Â°C. wetdays is the annual average number of days with precipitation greater than 0.25 mm. bound is a dummy variable (bound = 1 for bound base layer; bound = 0 for unbound base layer). WF is a dummy variable (wet-freeze zone = 1, else = 0). WNF is a dummy variable (wet non-freeze zone = 1, else = 0). DNF is a dummy variable (dry non-freeze zone = 1, else = 0). Note that in the DF zone the dummy variables of WF , WNF , and DNF should be equal to zero.

88 Table 20. Results of Multiple Regression Analysis for Coefficients in the First Faulting Model. Objectives Variable Parameter Estimate Standard Error t-value p-value Prediction model of ln Intercept 9.00365 0.45718 19.69383 7.11E-28 days32C 0.01921 0.00389 4.93505 6.70E-06 FT -0.00973 0.00194 -5.01109 5.08E-06 basethick -0.14165 0.05244 -2.70138 0.008964 Bound -0.63148 0.17016 -3.71107 0.000454 Objectives Variable Parameter Estimate Standard Error t-value p-value Prediction model of ln Intercept 0.85277 0.66645 1.27956 0.205976 dowel -1.47429 0.27515 -5.35824 1.64E-06 Bound -1.37916 0.28814 -4.78633 1.28E-05 days32C 0.03371 0.00708 4.76000 1.41E-05 basethick -0.14909 0.07012 -2.12621 0.037907 DNF -1.60578 0.47300 -3.39488 0.001268 WF -0.99420 0.29543 -3.36525 0.001386 WNF -1.27570 0.38121 -3.34647 0.001467 wetdays 0.01670 0.00523 3.19625 0.002289 Objectives Variable Parameter Estimate Standard Error t-value p-value Prediction model of ln Intercept 0.95239 0.31918 2.98391 0.004185 dowel -0.81622 0.21215 -3.84747 0.000304 rho -0.16416 0.03220 -5.09822 4.07E-06 No 0.00017 0.00005 3.78797 0.000368 intensedays 0.02471 0.01168 2.11502 0.038812 days32C -0.01274 0.00461 -2.76011 0.007758 basethick -0.08179 0.04858 -1.68352 0.097744 Drainage -0.45111 0.25342 -1.78008 0.080394 Figure 49 shows the results of the coefficient models by comparing the faulting parameters predicted by Eqs. 4.56, 4.58â4.59 with the raw parameters. To achieve a valid calibration model, 80 percent of selected LTPP sections were randomly taken for multiple regression and the rest were taken for validation. Figure 49 shows satisfactory and accurate correlations between faulting parameters and pavement performance-related factors.

89 (a) N0 (b) (c) Figure 49. Comparison between Measured and Predicted Coefficients in Faulting Model. Calibration Model of the Inflection Point Faulting Due to the importance of the inflection point in the proposed model, multiple regression analysis is also performed to calibrate the faulting depth. Table 21 lists the regression results for faulting at the inflection point. The faulting depth at the inflection point is expressed:

90 ln( ) 0.582841 1.515503 1.55997 1.570976 0.25955 0.013998 0.026763 32 ipf drainage dowel bound basethick wetdays days C ï½ ï ï´ ï ï´ ï ï´ ï ï´ ï« ï´ ï« ï´ (4.60) Figure 50 presents the comparison between measured and predicted faulting depth at the inflection point. The multiple regression of Eq. 4.60 is used to determine the predicted faulting depth. The regression model of Eq. 4.60 has the capability to accurately estimate the faulting depth with an R2 value of 0.89 in the validation plot of Figure 50. Figure 50. Comparison between Measured and Predicted Faulting at Inflection Point. Table 21. Results of Multiple Regression Analysis for Faulting at Inflection Point. Objectives Variable Parameter Estimate Standard Error t-value p-value Prediction model of ln Intercept 0.58284 0.77211 0.75487 0.453382 Drainage -1.51550 0.44948 -3.37166 0.001335 dowel -1.55997 0.20638 -7.55860 3.37E-10 Bound -1.57098 0.23781 -6.60602 1.34E-08 basethick -0.25955 0.07067 -3.67270 0.000525 wetdays 0.01400 0.00536 2.61277 0.011419 days32C 0.02676 0.00462 5.78941 3.02E-07 Figure 51 presents the mean critical faulting depth with or without dowels. The mean critical faulting depth is about 1.2 mm when dowels are used but roughly 2.4 mm when no dowels are used. Thus, Figure 50 shows that the use of dowels results in lower critical faulting depth.

91 Figure 51. Mean Critical Faulting Depth with or without Dowels. Faulting Modeling Based on Permanent Deformation Characterization The faulting depth accumulated before the inflection point (Figure 52) is related primarily to the permanent deformation of the underlying layers under traffic loads. To characterize the stress-dependent deformation behavior of unbound aggregates, a new ME faulting model is proposed to summarize faulting depths at different stress levels and expressed as: 0 ( ) 2 1 1 ( ) ( ) ( )i k N N i m ni i i a a J I Kf N f e P P ï¢ï² ï¡ï ï ï¥ ï½ ï« ï½ï¥ (4.61) 1 x y zI ï³ ï³ ï³ï½ ï« ï« (4.62) 2 2 2 2 2 22 1 [( ) ( ) ( ) ] 6 x y x z y z xy yz zx J ï³ ï³ ï³ ï³ ï³ ï³ ï´ ï´ ï´ï½ ï ï« ï ï« ï ï« ï« ï« (4.63) 2sin 3(3 sin ) ï¦ï¡ ï¦ ï½ ï (4.64) 6 cos 3(3 sin ) cK ï¦ ï¦ ï½ ï (4.65) where ( )if N is the total faulting depth. iN is the cumulative number of axles at axle load level i. is the number of the axle when faulting occurs. aP is atmospheric pressure, 101.305 kPa. 2iJ is the second invariant of the deviatoric stress tensor at stress level i. 1iI is the first invariant of the stress tensor at stress level i.

92 fï¥ , ï² , ï¢ , m , and n are model coefficients. xï³ , yï³ , and zï³ are normal stresses. xyï´ , yzï´ , and zxï´ are shear stresses. In this model, the two terms of 2J and 1I Kï¡ ï« are used to reflect the influence of stress state. Figure 52. Plastic Deformation Curve before the Inflection Point. The variables, as known in the model including predictive traffic variables ( iN ), power coefficients ( m and n ), stress state terms ( 2J and 1I ), cohesion c , and friction angle ï¦ , are required to perform nonlinear optimization analysis to determine the unknown coefficients (ï², ï¢ , 0N , and fï¥ ). Therefore, the determination of the variables is required by taking following steps: 1. The power coefficients of m and n are determined by Eqs 4.66 and 4.67 from the developed permanent deformation model. 2. The cohesion c and friction angle ï¦ are obtained by Eqs. 4.68 and 4.69 from the developed shear strength model. 3. The stress state terms of 2J and 1I are determined from the analytical elastic method that is elaborated in Appendix G. 4. The predictive traffic variable of iN represents the cumulative number of axles at the axle load level i. To determine the traffic variable, the load distributions are required to be rearranged into a tabulation of axle load levels by number of axles. Based on this established tabulation, the calculated faulting at various axle load levels were added up to the total faulting depth:

93 4 0.017 0.00967 0.990.2153 10 0.396P MDDm ï´ ï ï´ ï«ï½ ï´ ï« (4.66) 0.22993 0.746n mï½ ï ï´ ï (4.67) 20.22134 4.6 455.62 1262.75sat sc PI w Gï½ ï ï« ï« ï (4.68) 20.0014 0.2796 35.87PI PIï¦ ï½ ï ï« (4.69) where 4P is the percent of the base course passing the #4 sieve. MDD is the max dry density. w is the water content. satw is the saturated water content ( wsat d w n ï§ ï§ ï½ ï´ ). n is porosity. wï§ is the unit weight of water. dï§ is the dry unit weight. sG is the specific gravity of base course. PI is the plasticity index. Determination of Stress State in the Faulting Model The faulting model that incorporates the effects of a load spectrum on the permanent deformation of the base course requires a simple, consistent, and reasonably accurate method of estimating the state of stress in the middle of the base course. The method adopted for this faulting model is from Timoshenko (195) for estimating the stress state in a half space in which the only material property that is needed is the Poissonâs ratio. The tire load for both single and dual tires is represented as a point load on the concrete pavement surface. The use of these equations permits the superposition of the stresses caused by both dual tires. Appendix G introduces the calculations of stress in the base course for single and dual tires. Categorization of Traffic Loads The traffic input for the faulting model is the number of axles corresponding to different load levels. The weigh-in-motion (WIM) data can be used to determine the traffic input. The WIM tables present the number of axles by types of axles within vehicle classes at various load levels. The WIM data can be directly collected from the LTPP database. If the data are unavailable in the database, it can be estimated with the annual average daily truck traffic (AADTT) and axle load distributions. The axle distributions shown in WIM tables are required to be categorized into the groups (single and dual tire) to permit the determination of the stress state. Appendix H elaborates on the categories of axle distribution and the approach to estimate the number of axles with AADTT.

94 Comparison between Measured and Predicted Faulting Based on Permanent Deformation Characterization The validation of the model for the load-related faulting is based on the LTPP field faulting data. First, the field data were collected for validation of the first faulting model to characterize the entire life of faulting and determine the critical faulting depth at the inflection point. Then the field data before reaching the determined critical depth were retained for the validation of the second faulting model. The coefficients ( ï² , ï¢ , 0N , and fï¥) in the second load-related faulting model were determined by nonlinear optimization analysis. Figure 52 presents the comparison between the measured and predicted faulting before reaching the critical faulting depth at the inflection point. Figure 53 shows that the overwhelming majority of the scatter falls close to the identity line. Moreover, the value of R2 is 0.97 indicating a high degree of accuracy of the proposed second faulting model. Figure 53. Comparison between Measured versus Predicted Faulting before Critical Depth. Calibration of Model for Faulting on Deformation Characterization The calibration coefficients for the second load-related faulting model are similar to the first faulting model. Similarly, both the categorical and numerical variables are included in the calibration coefficients for the second faulting model and were collected from the LTPP database. The stepwise multiple regression analysis was performed to investigate the significant predictive variables with a p-value below 0.05 at a 95 percent confidence level. Table 22 presents the results of the multiple regression analysis for the coefficients in the second load- related faulting model. The model parameter of was solved to a constant value of 2e7. The regression equations for the coefficients are as follows:

95 0 429.470561 15.72206 180.18123 196.160669 0 831.66223 23.586441 32 N FI FT days C dowel days C ï½ ï ï ï´ ï ï´ ï« ï´ ï« ï´ ï« ï´ (4.70) 0.028917 0.000628 0.006929 0.016353 0.031231 0.02095 0.007663 0 FI FT DNF WF WNF days C ï¢ ï½ ï ï´ ï ï´ ï« ï´ ï ï´ ï« ï´ ï« ï´ (4.71) 0.184024 0.009054 0.060236 0.021429 f basethick bound dowel ï¥ ï½ ï ï´ ï ï´ ï ï´ (4.72) 2 7eï² ï½ (4.73) where 0days C is the annual average number of days with the temperature lower than 0Â°C. FI is the calculated freezing index for year (it can be collected from LTPP table of CLM_VWS_TEMP_ANNUAL). Table 22. Results of Multiple Regression Analysis for Coefficients in the Load-Related Faulting Model. Objectives Variable Parameter Estimate Standard Error t-value p-value Prediction model of N Intercept -429.471 592.975 -0.72426 0.473858 FI -15.7221 2.424347 -6.48507 2.03E-07 FT -180.181 26.09081 -6.90593 5.88E-08 days0C 196.1607 29.10503 6.739752 9.59E-08 dowel 831.6622 235.89 3.525637 0.001231 days32C 23.58644 4.154377 5.677491 2.26E-06 Objectives Variable Parameter Estimate Standard Error t-value p-value Prediction model of Î² Intercept 0.028917 0.028854 1.002181 0.322766 FI -0.00063 8.91E-05 -7.04431 2.47E-08 FT -0.00693 0.001018 -6.80987 5.09E-08 DNF 0.016353 0.031002 0.527485 0.601006 WF -0.03123 0.012386 -2.52151 0.016122 WNF 0.02095 0.019958 1.049707 0.300661 days0C 0.007663 0.001129 6.787173 5.46E-08 Objectives Variable Parameter Estimate Standard Error t-value p-value Prediction model of f Intercept 0.184024 0.019227 9.571042 2.51E-10 basethick -0.00905 0.003044 -2.97453 0.005982 Bound -0.06024 0.010982 -5.48506 7.38E-06 dowel -0.02143 0.011367 -1.88516 0.069825

96 Figure 54 presents the results of the comparison between the coefficients predicted by Eqs. 4.70â4.73. The correlations between those coefficients and pavement performance-related factors are shown in Figure 54 together with satisfactory R2 values. Figure 54. Comparison between Measured and Predicted Coefficients in the Permanent Deformation Faulting Model.