Development and validation of a pavement temperature profile prediction model in a mechanistic-empirical design framework

An accurate temperature prediction tool is an important part of any mechanistic-empirical (M-E) pavement design and performance prediction method. In this paper, a one-dimensional finite control volume (FCV) model is introduced that predicts the temperature within a pavement structure as a function of time and depth. The main input data required for the model are continuous time series of air temperature for conductive heat transfer, solar radiation for radiative heat transfer, and wind speed for convective heat transfer. The heat balance equation for each control volume of the FCV model is solved using an implicit scheme. To validate the numerical model, comparisons were made to measured temperature data from four test sections in Sweden located in regions with different climatic conditions. A good agreement was obtained between the calculated and measured temperature values within the asphalt layer, and temperature in the granular layers with the values of the coefficient of determination R 2 ranging from 0.866 to 0.979. The model is therefore suitable to be implemented as a pavement temperature prediction tool in M-E design.


Introduction
The infrastructure system is one of the main investments every modern society must make.It is therefore of great importance to base decisions on a well-founded design method and to have a good overview of the maintenance needed during a system's service life in order to minimize both construction and maintenance costs.To date, the main methods used worldwide to analyse and design pavement structures have relied on empirical correlations with past performance.However, these methods have severe limitations regarding performance prediction as background factors are often not very well understood as they are not based on the principles of engineering mechanics [1].New Mechanistic-Empirical (M-E) pavement design methods are therefore being developed in different countries in the world with the main purpose of adequately predicting pavement response and performance.The new M-E pavement design methods are capable of adapting to changing conditions in loading and environmental variables, capable in considering novel materials in the design process, while offering improved performance prediction [1,2].To further develop M-E methods, validation and calibration is required for the design inputs, distress analyses, and performance predictions.The validation process requires data that may be collected in the laboratory, by accelerated pavement testing, or by performing full-scale testing on roads [3].
Environmental factors such as temperature, moisture, and freezingthawing cycles are known to heavily influence the behavior of flexible pavements [4].The bearing capacity and performance of the asphalt layer are highly temperature-dependent while the performance of the granular materials in the base course and subbase are generally reduced at increasing moisture content [4][5][6][7].The influence of temperature on the performance of asphalt is generally attributed to the viscoelastic nature of the material [8].In cold regions, surface temperatures below 0 • C in combination with moisture from precipitation and a low amount of sunlight cause freezing of the surface.After the freezing period is over, the thawing of the granular layers causes an increase in moisture in the granular layers and lowers the bearing capacity of flexible pavements [9].The availability of data related to pavement temperature facilitates decision making processes related to winter maintenance [10], investigating the variations in year-to-year performance [11], and investigating the influence of changing climatic factors on pavement performance [12].Temperature-related information is also useful for binder selection in pavement design [13], as well as for evaluating the structural capacity results obtained by in-situ nondestructive tests such as falling weight deflectometer (FWD) [14].Furthermore, pavement temperature data can be used for the development and validation of new models and design tools.
Several models that aim to predict the temperature distribution in the cross-section of pavements are available in literature.Based on their approach the models can be classified as empirical, analytical, and numerical.Empirical pavement temperature prediction models consist of simple linear or nonlinear equations typically developed to perform statistical analyses of existing datasets such as to predict the minimum and maximum temperature.Examples of empirical prediction models found in literature include linear equations as described in the works of Bosscher et al. [15], Raad et al. [16], and Diefenderfer et al. [17], nonlinear equations as the ones described in the works of Huber et al. [18], Ovik et al. [19], and Park et al [20], and artificial neural network based approaches as described by Shao [21].A notable example of the implementation of empirical models in M-E pavement design are Long-Term Pavement Performance (LTPP) program's equations for maximum and minimum pavement temperatures [22].Their main drawback is that they do not consider physical processes directly.Thus their accuracy is established only for the original dataset by which they are developed [23].
Analytical methods aim to address the shortcomings of empirical methods by offering a derived closed form analytical solution to the heat transfer partial differential equation (PDE).A notable example of analytical methods used in pavement temperature prediction can be found in the work of Wang [24] where the analytical solution is derived based on Laplace transform in which several steady state solutions were superimposed.Analytical methods are capable of accurately predicting temperature changes in pavements; however, it is difficult to implement them into M-E pavement design due to the complexity of the boundary conditions [23].
Numerical models for pavement temperature prediction address the shortcomings of empirical and analytical methods by offering reasonable accuracy, possibility to consider complex boundary condition and flexibility in case of changes.When using numerical methods, the heat transfer partial differential equation (PDE) within a discretized domain of computation.Based on the discretization type of the computation domain, numerical methods can be categorized into finite difference (FD), finite element (FE), and finite control volume (FCV) methods.
FD methods are commonly used for pavement temperature prediction.They use differences based on Taylor series expansions to approximate the derivatives of the heat transfer PDE.The first example of FD methods being implemented to predict pavement temperature is the work of Dempsey and Thomas [25].The model was later included in the Enhanced Integrated Climatic Model (EICM) [26].EICM considers the radiation component of heat transfer using a regression equation and has some restrictions in considering discontinuities in thermal properties as a function of depth which may lead to temperature prediction inaccuracies [23].Other implementations of FD models in pavement temperature predictions have been made to extend the application to 2D problems [27], use different time discretization schemes such as explicit, implicit, or Crank-Nicolson [28], and consider a variety of boundary conditions [29].FD models are computationally efficient [23] but encounter difficulties when dealing with complex geometries, and 2D and 3D problems [30].
FE methods discretize the computational domain into finite elements to solve the heat transfer PDE.Examples of FE based pavement temperature prediction models can be found in the works of Ho and Romero [31], and Minhoto et al. [30].FE methods are capable of considering complex geometries, and can be easily coupled with mechanical response analysis of pavements which makes it easier to investigate the effect of temperature on mechanical properties [30,32].However, the technique is not numerically efficient which makes them not attractive when large computations are needed such as in the case of real-time surface temperature predictions.
An alternative approach to FD and FE methods is the FCV method.
When using this method, the media is first discretized into FCV ś, the energy balance is made for each FCV, and then the heat transfer PDE is solved for each control volume.This technique has been extensively used in building physics simulations and hygrothermal modeling of roofs [33] as well as to a limited extent in pavements [23].The method is computationally efficient, can consider complex geometry and boundary conditions, and can consider the variation of the thermal properties in the pavement layers since it is an integration based method [23].A comprehensive review of all the three methods in pavement temperature predictions problems is given in Chen et al (2019) [34].Currently, a new M-E pavement design and performance prediction tool is being developed in Sweden.Temperature prediction within the pavement cross-section is an essential module of the Swedish M-E pavement design tool under development.It is important for the temperature prediction model to be able to provide an accurate and fast prediction of the temperature mainly in the asphalt layer but also at larger depths to aid in estimating the frost penetration depth.
This paper presents a 1-D FCV model that predicts pavement temperature based on simple and easily available historical data of meteorological inputs.The base inputs for the model are air temperature, global radiation, and wind speed.The model is aimed to be used in temperature-related aspects of pavement design processes.It uses an implicit scheme to provide predicted temperature values at any given time interval.This was done with the aim of making it possible to correlate the pavement temperature to mechanical response and performance predictions.Temperature data obtained from in-situ field measurements over multiple years has been used for the validation of the model.

Description of the model
Based on the amount of heat transfer at the surface due to the interaction of the pavement with the surrounding environment, the energy balance is written for the upper boundary condition.After the upper boundary condition is set based on the environmental conditions, the structure is discretized, and each pavement layer is further divided into smaller thin finite control volumes.The energy balance equation is written for each finite control volume with energy being transferred between neighboring control volumes.A constant time step Δt of one hour is used for the time domain discretization.Fig. 1 illustrates the environmental variables considered at the upper boundary condition and the discretization of the structure into FCV ś.Since the model defined here is one-dimensional, it is not possible to assume heat transfer at the sides of the pavement.Nevertheless, the temperature changes with time in the horizontal direction are negligible to the changes in the vertical direction.
The energy balance for the upper boundary conditions was set based on i) conduction due to the temperature difference between the pavement surface and the surrounding air, ii) convection due to the wind movement at the pavement surface, iii) solar shortwave radiation, and iv) incoming and outgoing longwave radiation.The consideration of each component of the heat flow at the surface in the modelling process is described below:

Solar shortwave radiation
The sun emits high-frequency shortwave radiation due to the high temperature at its surface.Out of the total radiation transmitted from the sun to the earth, only a part of it called the diffuse incident radiation manages to reach the earth.Direct shortwave radiation is defined as the radiation coming from the sun that manages to reach the surface of the earth.The balance between direct and diffuse shortwave radiation depends on the cloudiness, meaning that in the case of lower cloudiness, the direct shortwave radiation is higher.The heat flux at the surface coming from shortwave solar radiation denoted as q sol,sw is the net solar radiation absorbed by the pavement surface, and it is given by: D. Saliko et al. where: α, is the albedo, a unitless reflection coefficient.
Q solar is the global radiation incident obtained from meteorological data or generated using empirical models.
Albedo is defined as the ratio of the reflected radiation on the surface and the incident radiation.The parameter depends on the condition of the pavement surface.The typical range of values for the albedo is 0.04-0.18for asphalt [35,36] depending on the level of wear, approximately 0.8 for freshly fallen snow, around 0.4-0.6 for old snow, and in the range of 0.3-0.4 for ice [35,37].A constant albedo value of 0.07 was used in this model.The effect of snow cover was neglected, assuming that snow cover was removed by the winter maintenance operations.

Incoming longwave radiation
A part of the solar radiation is absorbed by the atmosphere and emitted to the earth as longwave radiation.The incoming longwave radiation, q sol,in absorbed by the surface of the pavement is calculated in a similar manner to previous studies [36,38] as shown below: where: ε a is the absorption coefficient for the incoming longwave radiation.
T air is the air temperature in degrees Celcius, • C. Cloud coverage affects the incoming long wave radiation and therefore the energy balance at the pavement surface.Information about local cloud cover is generally difficult to obtain.Existing radiation models are usually empirical and do not consider cloud cover properly [34].In case metrological data regarding radiation is not available empirical model might be used however adjust for cloud cover is of extra importance in such cases [39].

Outgoing longwave radiation
The outgoing longwave radiation emitted from the surface of the earth to the atmosphere is calculated by using the Stefan-Boltzmann law, given by: where: ε, is the emissivity coefficient of the surface for the outgoing longwave radiation.
T s is the pavement surface temperature in degrees Celcius, • C. The surface emissivity coefficient depends on the type of surface.The emissivity coefficient of dry asphalt typically ranges between 0.85 and 0.90 [36,38].In the case of snow cover, the emissivity coefficient range is 0.82-0.99 and may reach up to 1.0 in the case of fresh snow [37].The surface temperature is required for the estimation of the outgoing longwave radiation, therefore, the surface temperature calculated at the previous time step was used.Since the time steps are small, the differences in calculated heat flux due to outgoing longwave radiation are small.
A previous study by Mirzanamadi et al. [40] was used as a reference for the values of emissivity and absorption coefficients of the surface.The emissivity and absorption coefficient were selected as 0.89 and 0.78 respectively.

Convection
The convective heat flux q conv between air and pavement surface that is caused by the flow or movement of air has been calculated as: where: h c is the convection coefficient, W/(m 2 K).Multiple empirical equations exist for the calculation of the convection coefficient h c .A comparison of convection coefficient equations can be found in the work of Dehdezi [41].Depending on the purpose of the application, these differences may be considered as negligible.The Vehrencamp model is used for the calculation of the convection coefficient h c which is expressed as a function of air temperature, surface temperature, and wind speed [25,42]: where a = 1.4 and d = 0.5 are empirical model parameters [23,42].U is the wind speed in m/s.

Conduction/Diffusion
Assuming constant thermal conductivity and isotropic layers, the heat conduction into each finite control volume can be expressed by the one-dimensional heat conduction equation as follows including temperature T = T(z,t): where: ρ is the density, kg/m 3 .
c p is the heat capacity, J/(kg K). k is the conductivity, W/(mK).z is the depth, m.An implicit time integration method was employed for the temperature calculations throughout the pavement cross-section depth.The implicit method was selected with the intention to have stable solutions independent of the selected time step.The selected time step for the calculation was 1 h.To initialize the simulation, the mean annual air temperature was used as an initial temperature value for all the control volumes.A smaller time step was used at the start of the computations and later increased to 1 h.The temperature is calculated at the mid-point of each FCV except for the first layer, see Fig. 1.The discretization of Equation ( 6) using an implicit scheme for a FCV leads to: where the superscript i denotes the discretization in time, and the subscript n the discretization in space.
The boundary between two different layers was considered through the differences in thermal conductivity, specific heat capacity and density between two neighbouring control volumes.
Modifications were made for the control volume at the surface boundary, as the temperature is calculated at the outer border of the control volume in order to calculate the surface temperature.Thus the temperature of the first control volume T 1 is equal to the surface temperature T s .The right side of the Equation (7) has been modified to account for conduction, convection modelled as a resistance at the surface, and radiation modelled as heat fluxes at the surface.The heat balance for the first control volume located at the surface is written as: For the bottom boundary condition, a constant temperature boundary condition with zero heat flow has been prescribed at a sufficiently large depth as recommended by Doré and Zubeck [4].The depth of 4 m has been used here Hermansson [38].To account for the constant temperature boundary condition, Equation (7) has been modified for the last finite control volume as expressed in Equation ( 9) below: where T G is the pre-set constant ground temperature equal to the mean annual air temperature.The index m indicates the last finite control volume.
FCV with a smaller thickness were used for the part of the structure closer to the surface.This was done to capture the rapidly changing temperature due to the upper boundary condition.For the deeper parts of the structure, a coarser discretization with thicker finite control volumes was used.This was done since fewer temperature changes are expected at larger depths.
The terms r n , h n , and ho n have been used in the computations to facilitate the notation of repeating terms: The terms r n , h n , and ho n were introduced with the aim of simplifying the subsequent Equations (13)(14)(15)(16)(17)(18)(19)(20).Using the abovementioned notations and simplifications, Equations ( 7), (8), and (9) can be written as: where a n , b n , and c n are given by: For the first control volume that is located at the surface, d n is given by: + q sol,sw + q sol,out + q sol,in + q conv (17) For the internal control volumes, d n is given by: For the control volume located at the bottom d m the boundary condition is prescribed by: The use of the implicit method leads to a system of linear equations in the form of: The system of linear equations is solved to obtain the calculated temperatures for each finite control volume.Since matrix A is tridiagonal, the Thomas algorithm which is a simplified form of Gaussian elimination has been used to solve Equation (20), without inversing the square matrix A.

Description of the test sites
For the validation of the FCV approach, temperature data from 4 different test sites in Sweden were used.The test sites were selected based on the cross-section of the structures and on their geographical location with the aim to validate the test model for a variety of pavement cross-sections throughout the entire country.
The first selected structure was E18, located in the E18 motorway between the cities of Västerås and Enköping at the coordinate 59 • 38 ′ 0.2 ′′ N, 16 • 51 ′ 14.0 ′′ E. The road on which the structure is located is a four-lane motorway, with two lanes in each direction and a median strip between the lanes.The test site is located on level terrain, in an open unshaded area.The pavement structure consists of a 20 cm asphalt layer, placed on top of an 8 cm granular base course, and a 100 cm subbase layer.
The test site Ullevileden was selected as it was the southernmost available test section, located at the coordinate 58 The geographical locations of the test sites and the corresponding pavement cross-sections are illustrated in Fig. 2.
The thermal properties for the materials in each layer of the pavement structures were obtained partially from previously performed tests and partially from literature.For the pavement located at test site E18, the values for the thermal conductivity of the asphalt layer (AC surface-, binder-and base course) were obtained through laboratory testing by Mirzanamadi et al. (2018) [43].For the unbound granular layers and subgrade soils, the values for thermal conductivity k, specific heat capacity c, and density ρ were obtained from the literature [35,44,45].For  the pavement structures located at the other test sites, a similar literature-based approach was used in selecting the thermal properties of the layers.The range of values of thermal properties for the unbound layers are taken from Andersland and Ladanyi (2013) [35], that is specific heat capacity around 800 -900 J/(kg⋅K) and a thermal conductivity in the range 1300 -1700 W/(m⋅K).The layer thicknesses and the selected thermal properties for each layer for all the test sites are given in Table 1.The thermal properties were selected without any optimization process in order to improve the fit.Due to limitations in the scope of the paper, simplifications had to be made in the selection of the thermal properties for each layer.The thermal conductivity of the asphalt layer was assumed as constant.The composition of the asphalt mixture consisting of the aggregates and filler used and the type of bitumen is known to affect the conductivity of the mixture [46,47].Furthermore, the thermal conductivity of both the asphalt mixture and the granular layers are variable due to the changes in external factors such as temperature, moisture, and freeze-thaw cycles [35,47].Obtaining suitable values for the thermal properties of each layer's material was one of the main challenges of the modelling process.
Air temperature were available on-site at all test sites.Wind speed data was available on-site at three stations (Ullevileden, Svappavaara and E-18) and at Långträsk the wind speed data was available at a meteorological stations operated by SMHI (The Swedish Meteorological and Hydrological Institute) located at 20 km distance from the site.Regarding radiation data a mesoscale radiation model called STRÅNG developed by SMHI for the Nordic countries was used.The model uses real-time radiation data measured at multiple locations by pyrheliometers combined with real-time meteorological data to generate values of global, direct CIE-weighted, and photosynthetically active radiation [48,49].An example of recorded hourly values of the air temperature, global incident solar radiation, and wind speed used as input for the modelling process for the test site E18 for a period of 3 years from 1 January 2019 until 31 December 2021 is shown in Fig. 3.
For all test sites, built-in temperature recording instrumentation was available on-site.Depending on the instrumentation setup available on each test site, either the surface temperature, the temperature at a certain depth inside the asphalt layer, or the temperature at a certain depth within the granular layers or a combination of them was available.Infrared radiometers were used to record the surface temperature.The temperature inside the asphalt layer was recorded using thermocouples.For the granular layers, the temperature and the frost penetration depth were recorded using a frost rod called Tjälstav 2004 [50].

E18 Test site
For the initial validation of the FCV model, test site E18 was selected.The test site has acted as a permanent road research station, collecting multiple types of data continuously.To evaluate the quality of fit, a comparison was made between the surface temperature value measured using an infrared radiometer and the corresponding values calculated by the FCV model for the years 2019-2021 as shown in Fig. 4.
Further comparisons were made to validate the FCV model's capability to predict the temperature in the granular layers.The calculated temperature values by the FCV model were compared against measured values recorded using a frost rod consisting of 13 sensors placed every 25 cm on a 3 m long rod starting from a depth of 20 cm down to a depth of 320 cm.The recorded temperature values in the pavement cross section from 1 January 2019 until 15 June 2020 are shown in Fig. 5 a).Some minor loss of data occurred due to malfunctions of the instrumentation.This is visible in Fig. 5 a) as gaps.The calculated temperature values for the same period are shown in Fig. 5 b).
From Fig. 4 and Fig. 5, it can be observed that the comparisons between the measured and calculated values exhibit a good fit.Thus, it is shown that the FCV model is able to predict the temperature throughout the whole year.For the surface temperature, minor deviations were observed mainly for days with lower temperatures during the months of January 2020, and January 2021.Looking closer at the meteorological data, it was noticed that the main deviations occurred during the dates when snowfall occurred.In the case of fresh snowfall, another layer with different thermal properties is added on top of the pavement, instead of direct exposure to air temperature.Furthermore, snowfall affects the emissivity of the surface which leads to changes in albedo and surface emissivity, thus affecting the radiative heat transfer.The prediction accuracy typically improved a few hours after snowfall which corresponds to the occurrence of winter maintenance procedures.Thus, the snowfall may be attributed to the minor accuracies in predicting the surface temperature during the winter months.The effect of snowfall was not considered in the model for simplicity reasons and because it is difficult to predict when the winter maintenance procedures will occur at a specific location.The temperature prediction in the granular layers was observed to be more accurate as shown in Fig. 5.
To further evaluate the prediction accuracy of the model, the hourly  recorded temperature values were plotted against the predicted temperature values.The coefficient of determination R 2 was used as a statistical measure of the prediction accuracy.For both the surface temperature and the temperature in the granular layer, a high R 2 value was calculated, 0.958 for n = 26305 data points for the surface temperature, and 0.959 for n = 341951 data points for the temperature in the granular layers.The plots and the R 2 coefficient are shown in Fig. 6.
To further evaluate the capability of the model to predict pavement behavior, the same comparison between the measured and calculated surface temperature was made for the average daily temperature, maximum daily temperature, and minimum daily temperature.This information is useful in pavement design for binder selection, decisions

Table 2
The average difference and standard deviation between the predicted and measured surface temperatures for each month at test site E18.  on maintenance procedures, prediction of distresses, etc.The surface temperature was selected as the variable for comparison since the largest temperature variations occur at that specific location.The comparisons are given in Fig. 7.The R 2 values were calculated as 0.991 for the average daily surface temperature, 0.982 for the maximum daily surface temperature, and 0.989 for the minimum daily surface temperature indicating that the model is capable of capturing both the daily surface temperature extrema and average values.For further comparison, the average difference between the measured and calculated temperatures, and the standard deviation values were calculated for each month of each year.The results are shown in Table 2.
From Table 2, the largest temperature deviations between the measured and calculated temperature values were noticed at the months of May, June and July.The reason for this is the difficulties in considering the radiative heat flux at the surface.The radiative heat flux at the surface is affected by multiple site-specific factors such as shading, and cloudiness thus making it difficult to have a precise estimate of it.

Ullevileden test site
At the Ullevileden test site, hourly temperature measurement data for the asphalt temperature recorded by thermocouples at 2 cm, 7.5 cm, and 13 cm depth was available.The recorded temperature values were compared against the calculated temperature values by the FCV model at each corresponding depth.The comparisons were made for two full years from 1 January 2019 to 31 December 2020 as illustrated in Fig. 8. Furthermore, the prediction accuracy was evaluated by plotting the calculated temperature values against the corresponding measured temperature values and calculating the R 2 coefficient for each depth as shown in Fig. 9.
The prediction accuracy was the lowest for the thermocouple located at 2 cm depth and highest for the thermocouple at 13 cm depth.The reason for this is the proximity of the thermocouple at 2 cm depth to the pavement surface where the variations in heat flux are the largest.This can be observed in Fig. 8 a) where larger differences are observed for the thermocouple at 2 cm depth compared to the other two thermocouples in Fig. 8 b) and c).Fig. 9 illustrates this better with the data points in Fig. 9 c) for the thermocouple at 13 cm depth being closer to the quality line with an R 2 = 0.979compared to the data points for the thermocouple at 2 cm depth with an R 2 = 0.940 as shown in Fig. 9 a).The overall quality of fit is nevertheless high for all the three depths which means that the model provides an accurate temperature prediction for the asphalt layer.

Långträsk test site
Långträsk test site was selected in order to investigate the suitability of the model in predicting the temperature distribution in low-volume roads typically consisting of thin surfaced asphalt pavements.At the Långträsk test site, the data available consisted of asphalt temperature data recorded by a single thermocouple placed at 2.5 cm depth within the asphalt layer and temperature data in the granular layers recorded by a 200 cm long frost rod with 41 built-in temperature sensors with a spacing of 5 cm.The comparisons between the measured and calculated values for the temperature in the asphalt layer recorded by the thermocouples, and in the granular layers by the frost rod are given in Figs.10-12.
The model provides an adequate overall prediction of the asphalt temperature that captures the temperature variation with time.The coefficient of determination R 2 is calculated as 0.866 which is relatively high.In some specific cases, the model encountered difficulties in predicting the maximum surface temperature values by overestimating the extrema in the temperature prediction during the summer months as shown in Fig. 10 a).The reason for this is attributed to the lack of information related to site-specific conditions.In this case, the test site is surrounded by tall trees that provide shading during certain times of the day depending on the position of the sun.This leads to an overestimation of the heat gain due to solar radiation which resulted in higher predicted temperatures compared to the measured values.Minor difficulties were encountered also in predicting the temperature in the granular layers during the freezing and thawing periods.The model does not consider the contribution of latent heat of fusion due to the lack of information regarding the moisture content in the granular layers.
By making a visual comparison between the measured and calculated temperatures, minor deviations were noticed at larger depths specifically in temperatures close to 0 • C as shown in Figs.11 and 12.The maximum frost depth in this location exceeded the 200 cm length of the frost rod, therefore it was not possible to measure the maximum frost penetration depth.It was however possible to check and compare the starting and ending date of the frost period and the frost penetration rate.Due to the latent heat component not being considered, the model predicted the frost to occur a few days before the actual frost was measured.A steeper, more quick frost penetration rate was calculated compared to the measured rate.Furthermore, minor inaccuracies were observed during the thawing period.The reason for this might be the change in thermal properties of the granular layers due to the increase in moisture content within the structure.Implementing a coupled temperature-moisture content model that includes the effect of latent heat of fusion would increase the prediction accuracy for these types of structures.

Svappavaara test site
To evaluate the suitability of the FCV model in predicting temperature changes for pavements located in cold regions, the Svappavaara test site located approximately 100 km north of the Arctic Circle was selected.For the Svappavaara test site, comparisons were made only between the measured asphalt temperature recorded by a thermocouple at 10 cm depth, and the corresponding measured temperature.The comparison was made for a 3-year period from 1 January 2015 to 31 December 2017.Temperature data for the granular layers was not available for the test site.The comparisons are illustrated in Fig. 13.
The comparisons in Fig. 13 illustrated that the FCV model is capable of predicting the temperature variations in a relatively adequate manner also in pavements located in cold regions with a calculated R 2 value of 0.947.Since the structure was located in an open area, the model was able to accurately predict the temperature during the summer months unlike in the case of the Långträsk test site.The model however calculated lower temperatures compared to the measured temperatures during the winter months, in particular for temperatures lower than − 20 • C due to the missing latent heat of fusion component and difficulties in considering the surface conditions during snowfall.Considering the latent heat of fusion in the calculations would offer significant improvements in calculating temperatures during the winter months for pavements located in cold regions.

Discussion
The 1-D FCV method is used here to predict temperature in pavement structures.The method requires an energy equilibrium at the surface based on air temperature, wind speed, and incoming and outgoing solar radiation.When solar radiation data was not available the mesoscale radiation model STRÅNG was used to generate the missing data.Further, a constant temperature at the bottom boundary at 5 m depth was exploited equal to the site annual average temperature.The used layers thermal properties were based on laboratory testing when available but otherwise literature data were employed.
To validate the model, a test section in which meteorological data and temperature measurement data were available was used.The validation was performed for a period of 3 years for hourly time steps.A generally good agreement was obtained for the surface temperature with a calculated coefficient of determination R 2 of 0.958 and for the granular layers with an R 2 value of 0.959.The average temperature difference per month between the measured and calculated temperature values was between 0.94 • C and 3.1 • C. Furthermore, the suitability of the model to be used as a temperature prediction tool in an M-E pavement design framework was investigated.This was done by comparing the hourly temperature predicted by the model to measured hourly temperature values at three other locations.Overall, a good agreement was obtained for all the test structures between the calculated and measured values of temperature within the asphalt layer, and temperature in the granular layers with R 2 values ranging from 0.866 to 0.979.Some seasonal variation in accuracy was observed and was the highest inaccuracy observed during spring and early summer months.The reason is probably related to the radiative heat flux at the surface such as shading and cloudiness that are difficult to estimate.
Some limitations of the model were identified while making the comparisons.The model does not account for site-specific variations such as shading and snow cover which in some cases affected the predicted temperature close to the surface.The model does further not include the contribution of the latent heat of fusion in the heat balance which leads to some inaccuracies in predicting the temperature in the granular layers, particularly for pavements located in frost prone regions.Since the additional heat contribution by the latent heat of fusion is not considered, lower temperatures were calculated in comparison with the measurements at the depths located in the vicinity of the freezing front.Furthermore, the changes in thermal conductivity, specific heat capacity, and density of the layers with moisture changes and freeze-thaw cycles were not considered.Considering the latent heat of fusion requires additional information related to the moisture content in the granular layers which was not available for this study.Since the bottom boundary condition is modelled using a constant temperature the application is limited to pavements located on top of subgrade soils and is not suitable for predicting temperature variations of pavements located f. ex. on bridge decks subjected to cooling from below.
Incorporating and coupling a moisture prediction model to the FCV temperature prediction model could lead to improvements in terms of accuracy.It is possible to incorporate moisture prediction into the existing FCV framework and couple them using a predictor-corrector type of algorithm.However, a larger amount of inputs would be required for the simulation such as the soil-water retention curves for each layer, precipitation data, and measured moisture data for the validation which exceed the scope of this study.Another possible improvement would be to extend the model for 2D applications to consider the variations at the shoulders of the road and 3D applications to predict temperature variations at slopes and embankments.

Conclusions
In this study, a one-dimensional FCV model capable of predicting the temperature changes within pavement structures has been developed.The model uses available historical meteorological input data; that is air temperature, solar radiation, and wind speed to predict the temperature variation on an hourly basis.An implicit time integration scheme was used to provide a stable solution regardless of the selected time step.Other required inputs are surface related properties such as surface emissivity, and absorption coefficient.Finally the layer thicknesses and the thermal properties (thermal conductivity, specific heat capacity, and density) for each pavement layer are a necessary input parameters.
The model was validated against test sections were meteorological data and temperature measurement data were available.Generally good agreement was observed between measured and calculated pavement temperatures.The highest inaccuracy was observed at the surface and decreased with depth.The coefficient of determination was R 2 was 0.87 or higher in all observed cases.
Including information related to site specific parameters such as  shading due to adjacent obstacles (i.e.trees or buildings) or cloud coverage could improve the temperature predictions especially during spring and early summer months.This is however usually information that is not available.Incorporating and coupling a moisture prediction model to the FCV temperature prediction model could lead to improvements in terms of accuracy.It is possible to incorporate moisture prediction into the existing FCV framework and couple them using a predictor-corrector type of algorithm.However, a larger amount of inputs would be required for the simulation such as the soil-water retention curves for each layer, precipitation data, and measured moisture data for the validation which exceed the scope of this study.
Another possible improvement would be to extend the model for 2D applications to consider the variations at the shoulders of the road and 3D applications to predict temperature variations at slopes and embankments.
The results of the study indicate that despite the abovementioned simplifications, the model is capable of predicting the temperature variations in the cross-section of pavements located in different climatic areas with reasonable accuracy.Furthermore, the model is suitable for incorporation as a part of a temperature calculations tool in an M-E pavement design method.

Declaration of Competing Interest
The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

Fig. 1 .
Fig. 1.Illustration of the heat transfer between the environment and a multilayer system and the modelling of the heat flux at the upper boundary.
• 26 ′ 05.7 ′′ N 15 • 36 ′ 07.1 ′′ E on a suburban road in an industrial area in Linköping.The cross-section of the pavement structure consists of a 15 cm asphalt layer, an 8 cm granular base course, a 47 cm subbase layer, and a 170 cm light fill.The structure is located in a non-shaded area, with steep shoulders making it a well-drained structure.The third test site Långträsk was selected since it is a thin pavement structure located on a low-volume road.The structure consists of only 4 cm of asphalt placed on a 15 cm base course.It is located at the coordinate 65 • 19 ′ 05.3 ′′ N 20 • 17 ′ 48.7 ′′ E in the northern part of Sweden along road 515 close to the village Långträsk.The structure is located in a shaded area with tall trees on both sides of the road.The side shoulders are not steep and snow cover occurs during the winter.The last section Svappavaara located close to the village Svappavaara along the road E45 at the coordinate 67 • 38 ′ 15.1 ′′ N 21 • 04 ′ 44.3 ′′ E was selected as it was the northernmost available test site.It consists of a 20 cm asphalt layer, a 10 cm base course, a 30 cm subbase layer, and a 70 cm granular material fill.

Fig. 3 .
Fig. 3.The recorded hourly values of meteorological variables a) air temperature.b) solar radiation and c) wind speed at the test site E18.

Fig. 5 .
Fig. 5. Comparison between a) the measured temperature in the granular layers by the frost rod and b) the calculated temperature using the FCV model at test site E18.Gaps in measurements are due to malfunction in the equipment.

Fig. 6 .
Fig. 6.Evaluation of the prediction accuracy at test site E18 for a 3-year period for a) the surface temperature recorded by infrared radiometer and b) temperature in the granular layers recorded by the frost rod.

Fig. 7 .
Fig. 7. Evaluation of the prediction accuracy at test site E18 for a 3-year period for a) the daily average surface temperature, b) the maximum daily surface temperature, and c) the minimum daily surface temperature.

Fig. 8 .
Fig. 8. Thermocouple measured vs calculated temperature in the asphalt layer at Ullevileden test site at a) 2 cm, b) 7.5 cm, and c) 13 cm depth.

Fig. 9 .
Fig. 9. Comparison between the temperature measured in the AC layer using thermocouples and the FCV model calculated temperature for Ullevileden test site at a) 2 cm depth, b) 7.5 cm depth, and c) 13 cm depth.

Fig. 10 .
Fig. 10.Thermocouple measured vs calculated temperature in the asphalt layer at Långträsk test site plotted a) on a time axis, and b) in comparison to each other.

Fig. 11 .
Fig. 11.Measured vs calculated temperature in the granular layers at Långträsk test site at a) 17.5 cm, b) 52.5 cm, and c) 162.5 cm depth.

Fig. 12 .
Fig. 12. a) The measured temperature in the granular layers by the frost rod and b) the calculated temperature using the FCV model at test site Långträsk.Data from 01.09.2018 to 31.06.2019 is shown.

Fig. 13 .
Fig. 13.Thermocouple measured vs calculated temperature in the asphalt layer at Svappavaara test site plotted a) on a time axis, and b) in comparison to each other.

Table 1
Thermal properties of the pavement cross-section layers at test sites E18, Ullevileden, Långträsk, and Svappavaara.