Research on Multi-Parameter Error Model of Backcalculated Modulus Using Abaqus Finite Element Batch Modeling Based on Python Language

Abstract

The error in modulus backcalculation is a crucial component in validating the rationality and reliability of results for engineering applications. The objective of this study is to identify the theoretical limitations associated with backcalculated modulus errors under typical parameter uncertainties and to determine the primary factors contributing to these errors. Firstly, using the actual measurements or data from the Long-Term Pavement Performance (LTPP) project, the statistical distributions of errors for typical parameters in the modulus backcalculation model were determined. Subsequently, a factor level table for orthogonal experimental design was developed, leading to the construction of 81 orthogonal design experimental schemes and their corresponding theoretical pavement structure models based on the actual error distributions. The deflection responses of 81 theoretical pavement structure models were then computed using an ABAQUS finite element batch analysis method devised in Python. Furthermore, a multi-parameter error model for modulus was established using multiple linear regression and variance analysis. Finally, the theoretical limitations of modulus errors under actual errors were analyzed. The results show that the errors of thickness, load amplitude and load frequency follow a normal distribution, while the distribution of backcalculated modulus errors follows an approximate mixed Gaussian distribution. When the errors of multiple parameters are combined randomly, the modulus errors range from −100% to 595%, and the probability of the modulus errors being less than 15% is highest in the asphalt surface layer, followed by the subgrade, and then the base and subbase layers. Within the same error range, the modulus error is random. However, with different error ranges, the overall level of modulus error increases in proportion to the size of those ranges. Compared to factors such as thickness, load amplitude, and load frequency, the errors in deflections have a highly contribution rate on the modulus errors exceeding 99%.

1. Introduction

The modulus back-calculation technique is an important method for rapidly and non-destructively obtaining the mechanical parameters of pavement structures using road surface displacement response data. Since the 1970s, modulus backcalculation has achieved fruitful development, with over ten backcalculation programs and even more than fifty subcategories having been publicly reported [1,2,3,4]. However, the controversy surrounding this technology continues to this day. In fact, there is no consensus regarding the best or most reliable and accurate backcalculation program [5]. One of the primary controversies surrounding the use of backcalculated modulus methods lies in the significant errors associated with these techniques. Several studies have highlighted these issues, with some researchers going as far as to consider the magnitude of these errors to be unacceptable for practical applications [6,7,8,9]. These concerns underscore the critical need for more accurate and reliable approaches to modulus determination.

Studies have shown that the deviation between simplified theoretical models and actual models is considered to be the fundamental cause of errors in backcalculated modulus [5,10,11]. The simplification of models introduces multiple sources of error, specifically including errors in the mechanical calculation model of the pavement structure, thickness of the pavement layers, load amplitude, load frequency, temperature, surface displacement response, as well as systematic and random errors in testing [12,13]. Irwin et al. found that mere random error is sufficient to cause more than double the error in backcalculated modulus [14]. Jooste, F. et al. concluded that the backcalculation error due to small variations in layer thicknesses was at least 5% [15].

To address the issue of error in backcalculated moduli, a large number of researchers have carried out correction studies on model deviations [16,17], including aspects such as pavement structure elasticity, viscoelasticity, linear and nonlinear mechanical models, temperature, load amplitude, load frequency, and deflection error, yielding some meaningful results [18,19]. However, existing studies primarily utilize theoretical calculations and statistical methods to analyze the impact of individual parameters on the backcalculated modulus error. The theoretical parameters are based on experience, and statistical conclusions have limitations in their scope of application. Some studies establish empirical relationships using backcalculated moduli [20]; the reliability of the backcalculation programs and the corresponding results is questionable, casting doubt on the research conclusions. Moreover, existing studies have not clearly defined the coupled influence of multiple sources of error in backcalculated moduli. They are unable to systematically address the limitations of the error in the backcalculated modulus caused by various parameters’ errors under actual measurement levels.

This study aims to construct a model for multi-parameter errors in pavement structural modulus and, in conjunction with the actual levels of error in typical parameters, to determine the theoretical limitations and key influencing factors of backcalculated modulus errors.

2. Typical Parameter Composition of Modulus Backcalculation Model

According to the general principles of modulus backcalculation [5], the model consists of three parts: the measured model, the theoretical calculation model, and the mathematical model, as shown in Figure 1.

For the measured model, the actual load applied to the actual pavement structure generates stress waves, which propagate outward from the load center, causing actual vertical displacement response at the pavement surface. The theoretical calculation model exists in relation to the measured model, with parameters such as the combination of structural layers, thickness, and modulus assumed to be known. Assumed load can be derived from actual load collected by sensors in the measured model. In the theoretical calculation model, the displacement response monitored according to the actual sensor placement. The mathematical model serves as the bridge between the measured and theoretical calculation models. The main parameters involved in the modulus backcalculation model including thickness and modulus of each structural layers, load, and deflection response.

3. Statistical Distribution of Actual Measurement Errors in Model Parameters

Studying the statistical distribution of the actual measurement errors of typical parameters in modulus backcalculation models using field measurement data or LTPP data.

3.1. Statistical Distribution of Errors in Pavement Structural Layer Thickness

Three-dimensional ground penetrating radar was employed to inspect the pavement structural layer thicknesses of six typical roads in Guangdong Province in China and tallies the discrepancies between the measured thicknesses of each structural layer per road and their corresponding designed thicknesses. The error in pavement structural layer thickness is defined as Equation (1).

$${\epsilon }_{h_{ri}}=h_{ri}-h_{di}$$

(1)

where $h_{ri}$ is the measured thickness of layer $i$ in the pavement structure, and $h_{di}$ is the designed thicknesses of layer $i$ in the pavement structure.

Taking the error in asphalt layer thickness as an example, the distribution of measured errors in asphalt layer thickness for each project is shown in Figure 2.

From Figure 2, it can be observed that the histograms of asphalt layer thickness errors for various projects are approximately “bell-shaped”, and the majority of the project’s thickness error distribution histograms have a high degree of overlap with their optimally fitted normal distribution curves. The distribution of asphalt layer thickness errors can be considered to be approximately normal. Similarly, the thickness error distributions of the cement stabilized base layer and crushed stone layer also follow an approximate normal distribution. The confidence intervals for the parameters of the approximate normal distribution of layer thickness errors are presented in

Table 1. Confidence intervals for the parameters of the approximate normal distribution of thickness errors in various structural layers.

Layer	$95\% \mathbf{Confidence} \mathbf{Interval} \mathbf{for} \mathbf{Parameter} \mathit{\mu }$	$95\% \mathbf{Confidence} \mathbf{Interval} \mathbf{for} \mathbf{Parameter} \mathit{\sigma }$
Asphalt layer	[−0.50 cm, 1.03 cm]	[0.35 cm, 3.44 cm]
Cement stabilized base layer	[−1.48 cm, 1.11 cm]	[0.58 cm, 1.11 cm]
Crushed stone layer	[−1.45 cm, 2.95 cm]	[1.07 cm, 1.98 cm]

.

3.2. Statistical Distribution of Errors in Pavement Structural Layer Modulus

Given that there is currently no publicly available and fully accepted modulus backcalculation programs, this paper investigates the distribution of in situ modulus errors in pavement structural layers based on data from LTPP project made publicly available by the Federal Highway Administration. The modulus error of the pavement structural layer is defined as Equation (2).

$${\epsilon }_{E_i}=E_{backcalculatedi}-E_{di}$$

(2)

where $E_{backcalculatedi}$ is the backcalculated modulus of layer $i$ in the pavement structure. $E_{di}$ is the designed modulus of layer $i$ in the pavement structure.

Using the data provided by the LTPP database, such as asphalt viscosity, gradation of asphalt mixtures, void content, and effective asphalt content, the asphalt layer modulus can be calculated according to the Witczak model [21,22]. Predicted asphalt moduli from Witczak model are shown in

Table 2. Parameters of Witczak model and the predicted asphalt layer modulus.

SHRP ID	Layer	Pass Rate/%	Cumulative Sieve Residue/%			${\mathit{V}}_{\mathit{a}}$/%	${\mathit{V}}_{\mathit{b}\mathit{e}\mathit{f}\mathit{f}}$/%	$\mathit{\eta }$/106 Poise		${\mathit{T}}_{\mathit{r}}$/°C	${\mathit{E}}_{\mathit{d}}$/ MPa
		${\mathbf{p}}_{0.075}$	${\mathbf{p}}_{4.75}$	${\mathbf{p}}_{9.5}$	${\mathbf{p}}_{19}$			60 °C	135 °C
31-0116	HMAC in base layer	2.2	50.0	32.0	5.0	8.35	4.80	3376.0	4.9	34.2	3663
39-0105	AC in surface layer	5.9	48.0	12.0	0.0	4.67	4.27	4949.0	5.1	28.5	15,663
48-0117	AC in surface layer	6.7	56.0	29.0	2.0	6.03	3.70	2286.0	3.7	34.0	6275

. The laboratory-measured elastic moduli of typical pavement structures in the LTPP database, such as PCC slab, semi-rigid base layers, granular base layers, and subgrade, were taken as the design modulus for the corresponding structural layers. Laboratory-measured elastic moduli of PCC slab, semi-rigid base layers, granular base layers, and subgrade are shown in

Table 3. Laboratory-measured elastic modulus of PCC slab, semi-rigid base layers, granular base layers, and subgrade, MPa.

SHRP ID	31-0116	48-0114	48-0121	48-7165	87-2811
Type of pavement	Flexible	Semi-rigid	Semi-rigid	Composite	Composite
PCC slab	/	/	/	44,818	40,508
Semi-rigid base layers	/	159.1	148.4	/	/
Granular base layers	178.6	175.7	174.1	/	/
Subgrade	97.1	/	54.5	95.3	62.0

.

Based on the comparison between the backcalculated modulus and the corresponding design modulus for specific pavement structures in the LTPP database, the error distributions of backcalculated moduli for different pavement structures and different structural layers were obtained; see

Table 4. Distribution parameters of backcalculated modulus errors.

Layer	Number of Valid Data	SHRP ID	Distribution		Distribution Parameters
			Shape Test	Chi-Square Test	${\mathit{\mu }}_1$	${\mathit{\sigma }}_1$	${\mathit{w}}_1$	${\mathit{\mu }}_2$	${\mathit{\sigma }}_2$	${\mathit{w}}_2$	${\mathit{\mu }}_3$	${\mathit{\sigma }}_3$	${\mathit{w}}_3$
Asphalt layer	56	31-0116	Bimodal	passed	−2371.5	455.8	0.47	1193.5	1583.7	0.53	/	/	/
	41	39-0105	Trimodal	passed	6661.0	2053.8	0.42	−8296.9	1869.2	0.39	−2372.6	2236.9	0.18
	32	48-0117	Unimodal	failed	−508.0	2525.9	1.00	/	/	/	/	/	/
PCC slab	24	48-7165	Bimodal	passed	1500.9	8100.8	0.77	−8164.9	702.9	0.23	/	/	/
	48	87-2811	Unimodal	passed	−17,930.8	6726.1	1.00	/	/	/	/	/	/
Granular base	32	31-0116	Unimodal	passed	−89.8	26.5	1.00	/	/	/	/	/	/
	56	48-0114	Unimodal	failed	301.5	84.3	1.00	/	/	/	/	/	/
	56	48-0121	Bimodal	failed	326.6	63.6	0.87	651.6	91.2	0.13	/	/	/
Semi-rigid base	56	48-0114	Bimodal	failed	338.2	60.1	0.53	526.6	116.8	0.47	/	/	/
	56	48-0121	Bimodal	failed	313.3	28.0	0.69	538.5	71.6	0.31	/	/	/
Subgrade	32	31-0116	Unimodal	passed	71.7	20.3	1.00	/	/	/	/	/	/
	56	48-0121	Unimodal	failed	266.9	18.3	1.00	/	/	/	/	/	/
	24	48-7165	Unimodal	passed	54.4	10.7	1.00	/	/	/	/	/	/
	36	87-2811	Bimodal	failed	60.3	73.4	1.00	/	/	/	/	/	/

.

From Table 4, it can be seen that the error distribution of backcalculated moduli in pavement structural layers follows an approximate mixed Gaussian distribution [23,24], and its mathematical model is shown in Equation (3).

$$f({\epsilon }_{E_i})=\sum _{j=1}^k{w}_j{\displaystyle \frac{1}{\sqrt{2\pi {{\sigma }_j}^2}}}{e}^{-{\displaystyle \frac{{\left(x-{\mu }_j\right)}^2}{2{{\sigma }_j}^2}}}$$

(3)

where $f\left({\epsilon }_{E_i}\right)$ is the probability density function of the error distribution for the backcalculated modulus of pavement structural layer $i$. $k$ is the number of Gaussian distributions that make up the mixture. $w_j$, ${\mu }_j$, and ${\sigma }_j$, respectively, represent the weight, mean, and standard deviation of each Gaussian distribution.

From Table 4, it can be observed that the error distributions of backcalculated modulus for asphalt layers and semi-rigid bases tend to favor bimodal or trimodal distributions, whereas the error distributions of backcalculated modulus for granular base layers and subgrade lean more towards a unimodal distribution. The more complex the structural layer itself and its environment, the more pronounced the multi-modality of the modulus error distribution for that structural layer.

3.3. Load Characteristics and the Statistical Distribution

Combining the characteristics of the measured FWD (CFWD-10T produced in China) load time-history curve and previous research experience [25,26], this study focuses on two parameters: load amplitude and load frequency. Taking a typical semi-rigid asphalt pavement as the research object, a 50 kN load is used to testing at a fixed point repeatedly 90–100 times. The load amplitude for each measurement is recorded, and the load frequency calculation model proposed by Leiva-Villacorta et al. is used to calculate the load frequency for each measurement [27]. The approximate normal distribution of load amplitude has a mean of 50.2 kN and a standard deviation of 1.15 kN. The approximate normal distribution of load frequency has a mean of 34.5 Hz and a standard deviation of 0.87 Hz.

4. Construction of Finite Element Batch Processing Model Based on Python Language

4.1. Standard Model of Pavement Structure

A standard four-layer pavement structure model was constructed, whose geometric parameters and mechanical parameters are shown in

Table 5. Parameters of the standard pavement structure.

Parameters	Asphalt Layer			Cement Stabilized Base Layer			Cement Stabilized Subbase Layer		Subgrade
Thickness/cm	18			40			20		-
Modulus/MPa	12,500			8500			5000		70
Poisson’s ratio	0.3			0.25			0.25		0.40
Damping coefficient	0.9			0.04			0.04		0.06
Density/kg/m3	2400			2300			2200		1800
Viscoelastic parameters (Prony series)
tau_i_Prony	10−5	10−4	10−3	10−2	10−1	1	10	100	1000
g_i_Prony	0.0689	0.0954	0.1917	0.2808	0.2177	0.0877	0.028	0.0089	0.0042

.

The length, width, and depth of the pavement model are taken as 20 m, 20 m, and 6 m, respectively. The model boundaries are set as fully constrained, with the layers set to be completely continuous. The load amplitude and load frequency are represented by a semi-sine function, see Equation (4).

$$F\left(t\right)=F_{max}\mathit{sin}w\left(t+\phi \right)$$

(4)

where $F\left(t\right)$ is the load-time history data. The load application area is taken as 0.2659 m $\times$ 0.2659 m, the standard load amplitude is set at 50 kN, and the standard load frequency is 34.5 Hz. $w=2\pi F_{freq}$. $\phi$ is the phase of the semi-sinusoidal load–time history curve, which is determined according to each actual measurement; see Figure 3.

The positions of the surface displacement sensors are arranged according to

Table 6. The layout scheme for the actual sensors or the displacement value extraction points during finite element analysis.

Sensors	D1	D2	D3	D4	D5	D6	D7	D8	D9
Distance from load center	0	20	30	60	90	120	150	180	210

.

4.2. Batch Parameterized Modeling Method

4.2.1. Construction of Error Parameter Sequence

Based on the distribution of errors in different structural layer thicknesses, different structural layer moduli, load amplitude, and load frequencies, a data sequence with errors for each specific parameter was constructed according to Equation (5).

$$C_{{\epsilon }_r}=C+{\epsilon }_r\left({\epsilon }_r~f\left(·\right)\right)$$

(5)

where $C_{{\epsilon }_r}$ is the data sequence that contained error factors for the specific parameters. $C$ is the parameter for the standard model. To simplify the study, the number of pavement structural layers was fixed at four layers, and the influence of parameter errors such as Poisson’s ratio, damping, density, and viscoelastic parameters of the pavement structural layers was ignored. Thus, there were, in total, nine parameters including thickness of each layer ($h_1$, $h_2$, $h_3$), modulus of each layer ($E_1$, $E_2$, $E_3$, $E_4$), load amplitude ($F_{max}$), and load frequency ($F_{freq}$). The values in ${\epsilon }_r$ were independently sampled errors following the error distribution function $f\left(·\right)$.

4.2.2. Batch Parameterization Modeling Process

Step 1: Creation of parameterized scripts.

The action commands of ABAQUS are stored in the abaqus.rpy file, which is located in the working directory and can be updated in real time, recording all action commands including view operations. Additionally, when saving the CAE model of ABAQUS, the .jnl file can also record the operational commands during the modeling process, but it does not include view operations and the code is more concise. The abaqus.rpy file and the .jnl file are the basic scripts for batch parameterized modeling in ABAQUS.

Step 2: Establishment of a parameterized model.

Based on the parameterized script files abaqus.rpy and .jnl, one can directly utilize Python commands to establish models containing design parameters within the ABAQUS finite element analysis software and write them into an .inp file.

Step 3: Model parameter modification.

Batch modification of model parameters in the .inp file through Python, resulting in $N$ number of .inp files.

Step 4: Batch computation and result extraction.

The $N$ number of batched .inp files are sequentially imported into the ABAQUS finite element software for computational analysis, resulting in $N$ number of .odb files. Typically, a .bat file is written to include jobs for the $N$ number of .inp files, and the .bat file is placed in the Abaqus working directory for execution. Abaqus can then sequentially complete the computations for the $N$ number of jobs.

4.2.3. Effects of Batch Parameterization Modeling Using Python Language

(1)

Accuracy. Taking the standard model in Table 5 as the research object, the differences between the deflections obtained by the batch analysis method and those obtained by manual interaction modeling and manual extraction are compared. The results obtained by the Python method are completely consistent with those obtained by the manual interaction method.

(2)

Application efficiency. Taking the standard model in Table 5 as the research object, directly modeling and analyzing in ABAQUS software (version 2021) using manual interaction takes about 20 min per model. If considering errors during operation, rework, and rest time, with an effective daily work time of 12 h, completing 5000 models would take over 140 days. In contrast, using the Python batch parameterization modeling method, running on a computer equipped with an AMD Ryzen 9 5900X 12-Core Processor at 3.70 GHz, each model takes only about 1 min, and it can achieve 24 h of uninterrupted automatic calculation with zero rework. The actual total time to complete the modeling and analysis of 5000 models is only about 3.5 days, which can improve the time efficiency by more than 40 times.

5. Multi-Parameter Error Model

5.1. Orthogonal Experimental Design

According to above researches, there are nine typical parameters closely related to the modulus. In a full factorial experimental design, even if each factor is set at only 9 levels, the total number of experiments exceeds 9⁹. Using a batch parametric modeling approach based on the Python, under the current computational resource configuration, the time required to complete all simulation experiments exceeds 747 years. To analyze the impact of multiple parameters’ errors, an orthogonal experimental design using the standard orthogonal design table L₈₁(9⁹) was adopted due to its efficiency in minimizing the number of experiments required while maximizing the exploration of potential interactions among variables [28]. Based on the distribution of errors for different parameters mentioned in Section 3, nine parameter sequences considering the corresponding errors were randomly constructed respectively. The orthogonal design factor level table is shown in

Table 7. Orthogonal design factor level table.

Level	${\mathit{E}}_1$/MPa	${\mathit{h}}_1$/m	${\mathit{E}}_2$/MPa	${\mathit{h}}_2$/m	${\mathit{E}}_3$/MPa	${\mathit{h}}_3$/m	${\mathit{E}}_4$/MPa	${\mathit{F}}_{\mathit{m}\mathit{a}\mathit{x}}$/kN	${\mathit{F}}_{\mathit{f}\mathit{r}\mathit{e}\mathit{q}}$/Hz
1	11,278	0.182	6454	0.385	4444	0.184	44	47.8	33.6
2	11,753	0.192	7499	0.392	4744	0.189	56	48.0	34.1
3	11,987	0.199	7594	0.399	4812	0.199	59	49.0	34.2
4	12,033	0.207	7979	0.401	4916	0.205	68	49.3	34.3
5	12,338	0.208	8195	0.407	4956	0.207	76	50.0	35.1
6	12,812	0.209	8353	0.408	5028	0.209	80	50.4	35.2
7	14,035	0.211	8449	0.409	5074	0.217	87	51.3	35.5
8	14,081	0.212	9102	0.413	5242	0.223	87	51.8	35.7
9	14,458	0.215	10314	0.417	5353	0.231	107	51.9	37.0

.

According to the standard orthogonal design table L₈₁(9⁹), an experimental plan was designed. A batch parameterized modeling analysis was carried out, conducting a total of 81 simulation experiments.

5.2. Multivariate Regression and Analysis of Variance

A multivariate linear regression method is used to construct a relationship model between the errors in the modulus of each structural layer and the errors in the thickness of each layer, the errors of load amplitude and load frequency, and the errors of deflections, as shown in Equation (6).

$$E_{\xi err}=a_{\xi }+\sum _{\xi =1}^3{x}_{\xi \xi }{h}_{\xi err}+\sum _{j=1}^2{y}_{\xi j}{F}_{jerr}+\sum _{i=1}^9{z}_{\xi i}{D}_{ierr},\xi =1,...,4$$

(6)

where $E_{\xi err}$ is the error of modulus of layer $\xi$ in the pavement structure. $h_{\xi err}$ is the error of thickness of layer $\xi$ in the pavement structure. $F_{jerr}$ is the error in the $j$-th load parameter. $D_{ierr}$ is the error in the $i$-th deflection sensor. $a_{\xi }$, $x_{\xi \xi }$, $y_{\xi j}$, and $z_{\xi i}$ are the fitting parameters of Equation (6).

Table 8. Regression model parameters and ANOVA significance test results of layer modulus.

Parameters	${\mathit{E}}_{\mathit{e}\mathit{r}\mathit{r}1}$	p Value	${\mathit{E}}_{\mathit{e}\mathit{r}\mathit{r}2}$	p Value	${\mathit{E}}_{\mathit{e}\mathit{r}\mathit{r}3}$	p Value	${\mathit{E}}_{\mathit{e}\mathit{r}\mathit{r}4}$	p Value
$a_{\xi }$	0.007	0.032	0.008	0.189 1	−0.004	0.560 1	0.027	0.013
$x_{\xi 1}$	−1.381	0.000	−0.414	0.012	−0.549	0.010	0.077	0.788 1
$x_{\xi 2}$	0.251	0.028	−0.825	0.000	−2.308	0.000	−1.185	0.003
$x_{\xi 3}$	0.110	0.023	−0.276	0.003	−1.010	0.000	−0.459	0.005
$y_{\xi 1}$	1.160	0.000	1.339	0.000	0.857	0.001	2.639	0.000
$y_{\xi 2}$	−0.447	0.001	−0.824	0.001	−0.019	0.949 1	−3.222	0.000
$z_{\xi 1}$	−41.650	0.000	58.500	0.000	−63.300	0.001	27.600	0.275 1
$z_{\xi 2}$	−0.100	0.996 1	−49.700	0.072 1	24.500	0.490 1	9.300	0.849 1
$z_{\xi 3}$	17.900	0.249 1	−23.900	0.411 1	101.300	0.010	−90.700	0.088 1
$z_{\xi 4}$	28.600	0.039	−25.500	0.316 1	−28.100	0.396 1	86.000	0.065 1
$z_{\xi 5}$	27.600	0.126 1	−13.300	0.689 1	−17.100	0.694 1	−68.600	0.257 1
$z_{\xi 6}$	−61.600	0.013	133.900	0.004	−92.100	0.121 1	122.600	0.136 1
$z_{\xi 7}$	43.200	0.048	−94.700	0.022	76.200	0.149 1	−52.700	0.467 1
$z_{\xi 8}$	−19.500	0.116 1	15.600	0.499 1	−17.400	0.561 1	−112.500	0.009
$z_{\xi 9}$	4.380	0.536 1	−2.400	0.855 1	15.300	0.378 1	76.100	0.003
$R^2$	0.979		0.967		0.656		0.968

¹ p value is bigger than 0.5 which indicates the variance test is not significant.

presents the regression fitting parameters and their ANOVA significance test results.

From Table 8, it can be seen that compared to thickness, load amplitude, and load frequency, the error in deflections has a more significant impact on modulus errors. Meanwhile, the multivariate linear regression models of asphalt surface layer modulus error, base layer modulus error, and subgrade modulus error with $h_{\xi err}$, $F_{jerr}$, and $D_{ierr}$ all have high significance, exceeding 0.96, while the regression model of subbase layer modulus error with $h_{\xi err}$, $F_{jerr}$, and $D_{ierr}$ has relatively lower significance.

To address the multi-parameter model issue of subbase layer modulus error, a relationship model is first constructed between the error in deflections and the errors in modulus, thickness, load parameters, etc., as shown in Equation (7).

$$D_{ierr}=a_i+\sum _{\xi =1}^3{x}_{\xi i}{h}_{\xi err}+\sum _{j=1}^2{y}_{ji}{F}_{jerr}+\sum _{\xi =1}^4{z}_{\xi i}{E}_{\xi err}, i=1,...,9$$

(7)

where the meanings of $D_{ierr}$, $h_{\xi err}$, $F_{jerr}$, and $E_{\xi err}$ are the same as in the previous Equation (6). $a_i$, $x_{\xi i}$, $y_{ji}$, and $z_{\xi i}$ are the fitting parameters of Equation (7), see

Table 9. Regression model parameters and ANOVA significance test results of deflection value.

${\mathit{D}}_{\mathit{i}\mathit{e}\mathit{r}\mathit{r}}$	${\mathit{D}}_{1\mathit{e}\mathit{r}\mathit{r}}$	${\mathit{D}}_{2\mathit{e}\mathit{r}\mathit{r}}$	${\mathit{D}}_{3\mathit{e}\mathit{r}\mathit{r}}$	${\mathit{D}}_{4\mathit{e}\mathit{r}\mathit{r}}$	${\mathit{D}}_{5\mathit{e}\mathit{r}\mathit{r}}$	${\mathit{D}}_{6\mathit{e}\mathit{r}\mathit{r}}$	${\mathit{D}}_{7\mathit{e}\mathit{r}\mathit{r}}$	${\mathit{D}}_{8\mathit{e}\mathit{r}\mathit{r}}$	${\mathit{D}}_{9\mathit{e}\mathit{r}\mathit{r}}$
$a_i$	0.007	0.007	0.007	0.007	0.007	0.007	0.008	0.008	0.008
$x_{1i}$	−0.115	−0.104	−0.098	−0.085	−0.075	−0.071	−0.068	−0.066	−0.064
$x_{2i}$	−0.354	−0.334	−0.328	−0.309	−0.285	−0.268	−0.253	−0.238	−0.227
$x_{3i}$	−0.217	−0.218	−0.217	−0.207	−0.190	−0.170	−0.153	−0.138	−0.126
$y_{1i}$	−0.634	−0.643	−0.648	−0.653	−0.641	−0.616	−0.581	−0.542	−0.499
$y_{2i}$	−0.171	−0.173	−0.174	−0.176	−0.173	−0.167	−0.156	−0.145	−0.133
$z_{1i}$	−0.244	−0.247	−0.249	−0.251	−0.247	−0.237	−0.223	−0.206	−0.188
$z_{2i}$	−0.169	−0.174	−0.177	−0.189	−0.205	−0.222	−0.237	−0.253	−0.267
$z_{3i}$	0.967	0.967	0.967	0.968	0.970	0.972	0.973	0.970	0.971
$z_{4i}$	−0.623	−0.628	−0.633	−0.648	−0.671	−0.696	−0.722	−0.742	−0.758
$R^2$	0.977	0.977	0.977	0.975	0.974	0.973	0.972	0.971	0.969

.

As can be seen in Table 9, there is a significant correlation between the error in deflections and multiple parameters such as modulus, thickness, and load parameters, with the regression model all exceeding 0.969.

According to Equation (7) and the related regression parameters in Table 9, combined with Equation (6), the relationship model for the error in subbase layer modulus is obtained, as shown in Equation (8).

$$\left\{\begin{array}{c}{E}_{3err}={\displaystyle \frac{1}{z_{3i}}}\left(D_{ierr}-a_i-\sum _{\xi =1}^3{x}_{\xi i}{h}_{\xi err}-\sum _{j=1}^2{y}_{ji}{F}_{jerr}-\sum _{\xi }^{\mathrm{1,2},4}{z}_{\xi i}{E}_{\xi err}\right)\\ E_{\xi err}\ge -100\%\end{array}\right., i=1,...,9$$

(8)

where the error of modulus, $E_{\xi err}$, must not be less than −100%, as this would result in a modulus value of less than 0 MPa, which is impractical. $E_{\xi err}\ge -100\%$ is applicable to all mathematical models related to modulus errors.

In order to determine which parameters most significantly affect the modulus errors, a sensitivity analysis was conducted using Equation (9) and (10):

$${Contri}_i={\displaystyle \frac{{{{\beta }_i}^{\prime }}^2}{{\sum }_{j=1}^n{{{\beta }_j}^{\prime }}^2}}\times 100\%$$

(9)

$${{\beta }_i}^{\prime }={\beta }_i\times {\displaystyle \frac{SD\left(X_i\right)}{SD\left(E_{err\xi }\right)}}, \xi =\mathrm{1,2},\mathrm{3,4}$$

(10)

where ${Contri}_i$ is the relative contribution rate of the error in the $i$-th parameter of the multiple linear regression model to the modulus error of the $\xi$-th structural layer. ${{\beta }_i}^{\prime }$ is the standardized regression coefficient of the $i$-th parameter error. ${\beta }_i$ is the regression coefficient of the $i$-th parameter error. $SD\left(X_i\right)$ is the standard deviation of the $i$-th parameter error. $SD\left(E_{err\xi }\right)$ is the standard deviation of the modulus error of the $\xi$-th structural layer. The contribution rates of errors in thickness, load peak and deflection to the modulus errors are shown in Figure 4.

It can be found that the errors in thickness, load amplitude, and load frequency have almost no influence on the modulus error, while the contribution rate of deflection errors dominates overwhelmingly exceeding 99%. The deflection error has the greatest impact on the modulus error of the asphalt layer, followed by the subbase layer, base layer, and subgrade. And the modulus error of asphalt surface layer is significantly related to the deflection errors at D1, D4, D6, and D7, while the modulus error of base layer is more significantly related to the deflection errors at D1, D6, and D7. The modulus error of the subbase layer has a significant relationship with the deflection errors at D1 and D3 and the error in the subgrade modulus is significantly related to the deflection errors at D8 and D9. This can help engineers recognize the importance of correcting deflection errors and also assist them in quickly adjusting the modulus errors of various layers based on different deflection errors.

6. Theoretical Limitations of Modulus Error

6.1. Theoretical and Actual Error Levels of Typical Parameters

At the theoretical level, the accuracy of the load sensor is 1% ± 0.1 kN [29,30], the accuracy of the displacement sensor is ±2% [31], and the theoretical accuracy of ground penetrating radar thickness detection is about ±3% [32].

When repeated measurements are taken, the coefficient of variation range for different displacement sensor’s data is between 1.2% and 3.2%, with an average value of approximately 2.0%. However, when testing is not conducted at fixed points, the variability level of each displacement sensor’s data can be as high as 30%. The coefficient of variation range for load amplitude in multiple repeated measurements is between 1.1% and 2.3%, with an average value of about 1.56%. The coefficient of variation range for load frequency is between 2.51% and 4.13%, with an average value of approximately 3.06%. After eliminating the influence of moisture content, the dielectric constant of the asphalt pavement is relatively stable, and the deviation of thickness measurement results can be controlled within 5.02%.

6.2. Modulus Errors Under the Level of Multiple Parameter Errors

6.2.1. Theoretical Modulus Errors Under Specific Combinations of Multiple Parameter Errors

When the errors in thickness, load parameters, deflections all zero, according to Equations (6) or (8), the theoretical modulus error of the asphalt surface layer is approximately 0.70%, the theoretical modulus error of the base layer is about 0%, the theoretical modulus error of the subbase layer is about −0.75%, and the theoretical modulus error of the subgrade is about 2.7%.

When the errors in thickness, load parameters, deflections all take the upper limit of the measured errors, according to Equations (6) or (8), the theoretical modulus error of the asphalt surface layer approximately ranges from −41.6% to 42.9%, the theoretical modulus error of the base layer ranges from −74.9% to 76.5%, the theoretical modulus error of the subbase layer ranges from −95.3% to 94.4%, and the theoretical modulus error of the subgrade ranges from −15.7% to 21.2%.

6.2.2. Theoretical Modulus Error Under General Combinations of Multiple Parameter Errors

(1)

Monte Carlo reliability method

In practice, the different parameters of the modulus backcalculation model are usually collected separately by multiple sensors. Therefore, the error of the actual measured data from different sensors are not uniform and fixed, and the combination of multiple parameter errors is random. It is difficult to exhaust all random combinations using general methods. The Monte Carlo method is a method that uses random numbers or pseudo-random numbers to guide the solution of reliability problems through probability statistical theory [33,34]. The Monte Carlo method relies on the law of large numbers in probability theory and the central limit theorem, obtaining an approximate solution to the problem through a large number of random experiments. The basic steps are as follows:

Step 1: Define the limit state equation.

According to the key content of modulus error research, the limit state equation can be defined as Equation (11).

$$P=P\left(\left|E_{\xi err}\right|\le {err}_{cri}\right)$$

(11)

where $E_{\xi err}$ is the error of modulus of the $\xi$-th layer. $\left|E_{\xi err}\right|$ is the limit state function obtained by Equation (6) or Equation (8). When $\left|E_{\xi err}\right|\le {err}_{cri}$, it is considered that the limit state is in a safe condition, while when $\left|E_{\xi err}\right|>{err}_{cri}$, it is considered that the limit state is in a failure condition.

Based on the Monte Carlo statistical simulation experiment method, under each random multi-parameter error combination, an $\left|E_{\xi err}\right|$ can be calculated according to Equation (6) or Equation (8). Each $\left|E_{\xi err}\right|$ is substituted into Equation (11) to test whether it satisfies $\left|E_{\xi err}\right|\le {err}_{cri}$. The number of simulation experiments that satisfy $\left|E_{\xi err}\right|\le {err}_{cri}$ are recorded. If, in a sufficient number of $N$ simulation experiments, there are $m$ simulation experiments with results in a safe state, then the safety probability, $P_s$, of the limit state equation can be calculated according to Equation (12).

$$P_s=P\left(\left|E_{\xi err}\right|\le {err}_{cri}\right)={\displaystyle \frac{m}{N}}$$

(12)

(2)

Level and probability of modulus error under random combination of multiple parameter errors

To determine the impact of different parameter error combinations on the modulus error of each structural layer, the typical error ranges of different parameters were identified based on the research in Section 6.1, as shown in

Table 10. Typical error range of different parameters.

Parameters	Thickness	Load Amplitude	Load Frequency	Deflections
Error range	−5.02~5.02%	−2.30~2.30%	−4.13~4.13%	−3.2~3.2%

. For each parameter involved in Table 10, 10,000 error data points were randomly generated within its error range, as shown in Figure 5.

Allowing the error combinations of the nine parameters in Figure 5 to be randomly combined and substituted into the multiple linear regression Equation (6) or Equation (8), the theoretical errors of the modulus for each structural layer were calculated. The results are shown in

Table 11. Analysis of theoretical error in modulus estimates for varying structural layers.

Modulus Error	${\mathit{E}}_{1\mathit{e}\mathit{r}\mathit{r}}$	${\mathit{E}}_{2\mathit{e}\mathit{r}\mathit{r}}$	${\mathit{E}}_{3\mathit{e}\mathit{r}\mathit{r}}$	${\mathit{E}}_{4\mathit{e}\mathit{r}\mathit{r}}$
$\mathrm{Valid} \mathrm{samples} (E_{\xi err}>-100\%$)	7292	6093	10,000	6556
Error range	−100~297%	−100~595%	−40.6~402%	−100~459%
Average of errors	69.4%	191.9%	180.7%	136.2%
$\mathrm{Number} \mathrm{of} \mathrm{samples} \mathrm{meeting} \left|E_{\xi err}\right|<15\%$	692	333	421	438
Safety probability of limit state	9.5%	5.5%	4.2%	6.7%

and Figure 6.

From Table 11 and Figure 6, it can be seen that the theoretical modulus error caused by random combinations of multiple parameter errors ranges from −100% to 595%. The combination of multiple parameter errors has the potential to reduce the level of modulus error through coupled influence. Under the coupled influence of multiple parameter errors, the theoretical modulus error of the asphalt surface layer is the smallest, followed by the subgrade modulus and the subbase layer modulus, while the theoretical modulus error of the base layer is the highest. According to the statistical results of modulus error for different structural layers determined by the Monte Carlo method, the probability of modulus errors being less than ±15% is the highest for the asphalt surface layer, followed by the subgrade, base layer, and subbase layer.

It is noted that when the level of asphalt layer modulus error varies, there is no significant difference in the distribution of random combinations of different parameter errors. This implies that at the same level of modulus error, the errors of different parameters can always be fully distributed within the entire range of the corresponding parameter errors. Under the same parameter error range, the relationship between the modulus error of each pavement structure layer and the multi-parameter error combination is relatively random, with no clear corresponding relationship.

(3)

Impact of multi-parameter error range on modulus error

To study the impact of different parameter error range on modulus error, the error range of different parameters in Table 10 were used as standards. A set of new multi-parameter error ranges were obtained by setting the original range to 70%, 50%, and 30%, respectively. These new ranges were then substituted into Equations (6) and (8) to calculate the mean modulus errors for the pavement structural layers, as shown in

Table 12. Mean modulus error of pavement structural layers under different error range.

Layer	${\mathit{E}}_{1\mathit{e}\mathit{r}\mathit{r}}$	${\mathit{E}}_{2\mathit{e}\mathit{r}\mathit{r}}$	${\mathit{E}}_{3\mathit{e}\mathit{r}\mathit{r}}$	${\mathit{E}}_{4\mathit{e}\mathit{r}\mathit{r}}$
Range in Table 10 × 100%	69.4%	191.9%	180.7%	136.2%
Range in Table 10 × 70%	49.9%	138.8%	139.4%	91.6%
Range in Table 10 × 50%	24.1%	84.4%	99.4%	52.1%
Range in Table 10 × 30%	4.5%	33.0%	59.4%	16.6%

. Additionally, the safety probabilities for the limit state where the modulus error less than ±15% for each layer were determined and are presented in

Table 13. Safety probabilities for the limit state where the modulus error less than 15%.

Layer	${\mathit{E}}_{1\mathit{e}\mathit{r}\mathit{r}}$	${\mathit{E}}_{2\mathit{e}\mathit{r}\mathit{r}}$	${\mathit{E}}_{3\mathit{e}\mathit{r}\mathit{r}}$	${\mathit{E}}_{4\mathit{e}\mathit{r}\mathit{r}}$
Range in Table 10 × 100%	9.5%	5.5%	4.2%	6.7%
Range in Table 10 × 70%	11.6%	7.2%	4.6%	8.5%
Range in Table 10 × 50%	14.6%	9.2%	6.6%	10.7%
Range in Table 10 × 30%	21.5%	13.1%	11.3%	15.4%

.

From Table 12 and Table 13, it can be seen that the smaller the multi-parameter error range, the smaller the mean theoretical modulus error of different structural layers, and the higher the safety probability for the limit state where $\left|E_{\xi err}\right|<15\%$ is satisfied. Combined with Table 11 and Figure 6, it is not difficult to infer that under the same error range (measurement accuracy control standards), the modulus error exhibits randomness. However, under different error ranges, the overall level of modulus error is directly proportional to the size of the error ranges. This indicates for the application of modulus backcalculation techniques that controlling the measurement errors of multiple parameters or improving the accuracy of sensors has a significant effect on controlling the error in backcalculated modulus.

7. Conclusions

The main conclusions are as follows:

(1)

The errors in layer thickness, load amplitude, and load frequency follow a normal distribution, while the distribution of errors in backcalculated moduli follows approximate mixed Gaussian distribution exhibit characteristics such as unimodal, bell shape bimodal and multimodal. This provides the real parameters needed for large-scale, high-volume pavement structure simulation calculations, which will effectively improve the accuracy of pavement structure simulation calculations.

(2)

The computational results of the ABAQUS finite element batch modeling and processing method based on Python are consistent with those of human–computer interaction, and the efficiency is increased by more than 40 times. This will provide a new computational tool for the large amount of pavement structure simulation analysis.

(3)

The theoretical modulus error caused by random combinations of multiple parameter errors ranges from −100% to 595% and exhibits randomness. The probability of modulus errors being less than ±15% is the highest for the asphalt surface layer at 9.5%, followed by the subgrade at 6.7%, base layer at 5.5%, and subbase layer at 4.2%. Under the same error range (measurement accuracy control standards), the modulus error exhibits randomness. However, under different error ranges, the overall level of modulus error is directly proportional to the size of the error ranges.

(4)

Compared to thickness, load amplitude, and load frequency, the deflection error has a highly contribution rate on the modulus errors exceeding 99%, which indicates that when applying modulus backcalculation techniques, special attention should be paid to controlling the errors in deflection measurement.

Future Works

The research in this paper has clearly demonstrated the significant impact of surface displacement response on modulus error. In the next phase, the characteristics of surface displacements and their influence including cracks, external loading, and seasonal effects will be investigated. And the different types of FWD will be verified.