Research Article  Open Access
Yuzhi Zhang, Jianghui Bei, Pei Li, Xiaojie Liang, "Numerical Simulation of the ThermalHydroMechanical Characteristics of HighSpeed Railway Roadbeds in Seasonally Frozen Regions", Advances in Civil Engineering, vol. 2020, Article ID 8849754, 14 pages, 2020. https://doi.org/10.1155/2020/8849754
Numerical Simulation of the ThermalHydroMechanical Characteristics of HighSpeed Railway Roadbeds in Seasonally Frozen Regions
Abstract
A multiphysics mathematical model of highspeed railway (HSR) roadbeds is necessary to facilitate a good level of understanding of the frost heaving mechanism. Based on the classical hydrodynamic model and fundamental thermoelasticity theories, we propose a thermohydro coupled model, based on the soilwater characteristic curve and solidliquid ratio as the relation equations, with the effects of the icewater phase change and water migration due to temperature change considered. With the linear expansion coefficient related to the temperature and the mass of ice content in roadbeds as the relation equation, we establish a macroscopic thermalhydromechanical model for unsaturated soil to calculate the roadbed deformations. Based upon the field data of a typical cross section of the HarbinDalian HSR roadbed, the variation of the thermalhydromechanical characteristics is simulated and studied. The results demonstrate that the increase of water content in the roadbed’s central line mainly appears in soil layers at depths less than 1.2 m and most icecontaining soil layers are at depths less than 0.6 m. Under the driving force of thermal and hydraulic migration, the vertical displacement of the west shoulder is increased to 18 mm. Then the settled maximum surface unevenness reaches 16 mm between the shoulder and centre line.
1. Introduction
The coarsegrained fillings and wellgraded broken stone widely used in HSR roadbeds are generally regarded as materials with little or no frost heaving. However, field monitoring of the HarbinDalian HSR roadbed in seasonally regions of China reveals that frost heaving is commonly present [1, 2]. Hence, a complete understanding of the mechanism of freezethaw deformations of HSR roadbeds in cold regions is essential. Currently, the mechanism of freezethaw deformations in roadbeds is most often studied inside a laboratory by considering the frost heaving behaviour of unsaturated coarsegrained soil with respect to macroscopic factors such as particle size and grade [3, 4], but the heat and water migration in HSR roadbeds and the variation of their thermalhydromechanical characteristics have seldom been studied.
Due to limitations of the testing conditions and long duration, neither laboratory tests nor field monitoring can completely clarify the freezethaw deformations of coarsegrained fillings. Therefore, a thermalhydromechanical numerical model is useful and facilitates an understanding of the mechanism of freezethaw deformations of coarsegrained fillings. Typical mechanical models of frozen soil already presented in the literature include a hydrodynamic model [5, 6], a segregation potential model [7, 8], a rigid ice model [9–12], and a thermodynamic model [13]. Based on the classical theory of hydrodynamics, various thermalhydromechanical coupled models that consider the effects of volumetric strain caused by phase transformations between ice and water that occur in situ and the water migration in soil deformation have been proposed [14–16]. For example, Lai et al. [17] proposed a thermalhydromechanical coupled model of saturated soil in which the ice content function is introduced and the Clapeyron equation is used as the equilibrium condition of the ice/water phase transformation. Similarly, thermalhydromechanical coupled models have been established based on the partial expression of frozen soil components [18, 19]. Also, based on the similarity to frost heaving due to the ice/water phase transformation and the thermal linear expansion of materials, unsaturated thermalhydromechanical frost heaving models have been proposed [20, 21].
In previous studies, thermalhydromechanical models were constructed by considering the Clapeyron equation generally used to describe phase equilibrium conditions. Also, the applicability of the Clapeyron equation under unsteadystate conditions is controversial, and moreover, it cannot distinguish between different soil characteristics such as particle composition, water content, and saturation degree, which are very important for practical engineering. The segregation criterion in the thermalhydromechanical models is often empirical and cannot reach a unified judgment. A practical and realistic thermalhydromechanical coupled model has not been put forward yet, since the existing thermalhydromechanical models are often very complicated and some physical parameters are difficult to be measured. At the same time, previous studies explored moisture migration, frost heaving mechanism, and frost heaving variation of frozen soils in different ways, but most of the existing calculation models have only been verified through laboratory tests. Moreover, few of these consider the variation of the thermalhydromechanical characteristics of HSR roadbeds in seasonally frozen regions, especially those with sunnyshady slopes.
Therefore, based on the hydrodynamic model and the basic theory of elasticity, a unified soilwater characteristic curve (SWCC) is adopted to represent the soil water potential in unsaturated soil. The permeability coefficient related to the volumetric water content is introduced, and the parameter changes due to the water content in unsaturated soil are considered. The specific solidliquid ratio is used as the coupling equation, and thus, the soil freezing law is looked at from the macroscopic perspective, instead of the microscopic study on the segregation criterion and the ice lens initiation mechanism. Meanwhile, the influence of the transformation of ice into water and water migration on the temperature field is considered, so as to realize the full hydrothermal coupling. The linear expansion coefficient, which is related to the soil temperature and mass of water content, is introduced into the mechanics of porous media equilibrium equation as the coupling equation between the hydrothermal state and the stressstrain relationship. The step function of hydrothermal characteristic parameters for frozen and unfrozen soil is used to accelerate the convergence. Subsequently, the thermalhydromechanical mathematical model for unsaturated soil is constructed. Also, based on the field data of a typical section of the HarbinDalian HSR roadbed, a model considering the sunnyshady slope effect is developed and verified. Additionally, the variation of thermalhydromechanical characteristics of HSR roadbeds in seasonally frozen regions is investigated, as well as the variations of the temperature field, moisture field, and deformation field.
The organization of the remaining sections of this article is as follows. In Section 2, the development of the thermalhydromechanical model is presented, followed by an application of the model to the HSR roadbed in Section 3. In this latter section, the achieved results are presented and discussed. In Section 4, a conclusion is made.
2. Development of the ThermalHydroMechanical Model
In this section, the development of the thermalhydromechanical model is presented.
2.1. Basic Assumptions
In deriving the model, the following basic assumptions were considered:(1)The moisture migration mainly occurs at the frost front and unfrozen area, while that in the frozen area is negligible(2)Moisture migration is only in the form of liquid water. In other words, migration of vapors (i.e., sublimation and condensation of moisture) and moisture evaporation at roadbed boundaries are excluded(3)Moisture migration can be described by Darcy’s law(4)Soil particles, ice, and water are incompressible substances with constant densities(5)The change of solute is neglected, namely, the osmotic pressure caused by concentration difference is neglected(6)Soil is an isotropic linear elastic continuum material that is completely consolidated and no account of creep is taken into consideration(7)Once the mass of ice content exceeds a critical level, the linear expansion coefficient of soil follows a linear relationship with respect to the mass content of ice
2.2. Governing Equations
2.2.1. Moisture Migration
According to the principle of mass conservation and Darcy’s law, the governing equation of moisture for unsaturated soil is given as follows [5]:where is the volumetric content of unfrozen water; t is the time (s); is the ice density (kg/m^{3}); is the water density (kg/m^{3}); is the volumetric content of ice; is the permeability coefficient (m/s); and is the soil water potential.
The soil water potential can be calculated by considering the generalized Clapeyron equation and the hydrodynamic model. Since the Clapeyron equation does not include any soil properties, it is only applicable to a state of phase equilibrium and can hardly measure the icewater phase pressure. Hence, the soil water potential is calculated by the hydrodynamic model in this study. The soil water potential of unsaturated soil mainly consists of the matrix potential and the gravity potential. By introducing the diffusion coefficient , equation (1) can be rewritten as follows:where is the water diffusion coefficient in soil (m^{2}/s); is the soil permeability coefficient related with ; C_{m} is the specific water capacity (m^{−1}); and is the matric potential (m). In most cases, the relation of matrix potential and is described by the unsaturated soilwater characteristic curve (SWCC).
2.2.2. Temperature Field
Since the temperature field distribution is affected by both the moisture field and ice/water phase transformation, the heat transfer induced by migrating water in the hydrodynamic model cannot be ignored [6]. Considering the effects of migrating water and the ice/water phase transformation on heat transfer in soil, the control equation for frozen soil temperature is given bywhere C is the soil volumetric heat capacity (kJ·m^{−3}·°C ^{−1}); T is the temperature (°C); is the volumetric heat capacity of water (kJ·m^{−3}·°C^{−1}); λ is the thermal conductivity (W·m^{−1}·°C^{−1}); and L is the latent heat of the ice (kJ/kg).
2.2.3. StressStrain Relationship
In order to simplify the model and according to the basic assumptions, the effects of soil gravity and expansion induced by ice/water phase transformation on soil deformation are considered, but the effects of train loads on soil deformation are neglected. The stressstrain control equation of frozen soils is given bywhere is the elasticity matrix, is the strain matrix, is the thermal strain matrix, and is the volumetric force matrix:where E is the elastic modulus (Pa) and is Poisson’s ratio.
As for the thermal expansion of common materials, linear expansion coefficients are introduced to describe thermal strains in frozen soils induced by the ice/water phase transformation frost heaving:where is the total strain matrix, , , and are components of the normal strain of thermal expansion in x, y, and z directions, respectively. For isotropic frozen soils, , where is the linear thermal expansion coefficient (1/°C) and it is related to the mass content of ice and temperature; is the reference temperature (°C). are components of the shear strain of thermal expansion in x, y, and z directions, respectively. For isotropic frozen soils, .
2.2.4. Coupling Equations
From equations (2) and (4), two differential equations and three unknown basic variables (, , and T) are identified for coupling the moisture and temperature fields, and so one more coupling equation is required to close and solve the system of equations. Indeed, liquid water cannot completely transform into ice during freezing, and consequently, unfrozen water is present around the soil particles. The dynamic equilibrium of the volumetric content of unfrozen water and temperature can be characterized by the soil freezing characteristic curve (SFC) or other empirical equations. Based on the model for the prediction of unfrozen water [22], the solidliquid ratio is introduced to describe the relation of the volumetric content of ice and volumetric content of unfrozen water and temperature, thus establishing a coupling equation. The icecontaining state of soil can be described by considering the temperature and the content of unfrozen water as shown in equations (8) and (9) to effectively improve the efficiency of the calculation:where is the initial gravimetric moisture content (); is the mass content of ice (); is the mass content of unfrozen water (); is the soil dry density; is the soil temperature (°C); is an empirical constant related to soil texture [22]; and is the freezing temperature of water (°C). Also,where is the solidliquid ratio and can be described by and T. A complete hydrothermal coupling analysis can be achieved by considering equations (2), (4), and (9). Based on the solution of the moisture and temperature fields, the thermalhydromechanical coupling analysis can be realized by combining equation (5) with the linear thermal expansion coefficient () as the coupling equation.
2.2.5. Model Parameters
The ice content calculated by the hydrothermal model is combined with the hydrodynamic frost heaving model. The determined frost heaving coefficient based on the ice content was substituted by the thermal expansion coefficient in COMSOL (a Multiphysics software). The linear thermal expansion coefficient of clay loess in the foundation was set to be −0.0015, while the linear thermal expansion coefficients of wellgraded broken stone and coarsegrained fillings in roadbeds [20] were determined by the following equation:where is the linear thermal expansion coefficient of the wellgraded broken stone layer (1/K); is the linear expansion coefficient (1/K) of coarsegrained fillings; and is the mass content of ice.
To facilitate the model’s numerical convergence, the hydrothermal parameters of soil [22] were selected by considering step functions so that these are continuous and differentiable:where is the volumetric thermal capacity of the soil (kJ·kg^{−1}·°C^{−1}); is the specific heat capacity of soil skeleton (kJ·kg^{−1}·°C^{−1}); is the specific heat capacity of water (kJ·kg^{−1}·°C^{−1}); is the specific heat capacity of ice (kJ·kg^{−1}·°C^{−1}); is the specific heat capacity of frozen soil skeleton (kJ·kg^{−1}·°C^{−1}); is the specific heat capacity of thaw soil skeleton (kJ·kg^{−1}·°C^{−1}); H (, ) is the Heaviside step function, as shown in Figure 1; and is the 1/2 phase transition interval. Thus, the soil thermal conductivity λ (W·m^{−1}·°C^{−1}) is considered to vary according to the following equation:where is the frozen soil thermal conductivity (W·m^{−1}·°C^{−1}) and is the thaw soil thermal conductivity (W·m^{−1}·°C^{−1}). Also, the soil permeability K (m/s) is considered to vary according to the following equation:where is the frozen soil permeability (m/s) and is the thaw soil permeability (m/s).
Based on the model proposed by Van Genuchten [23], which is suitable for the calculation of the permeability coefficient of unsaturated soil, and the specific water capacity C_{m} are calculated by the following equation:where is the saturated permeability coefficient (m/s); is the relative permeability; o, m, and k are empirical constants related to soil properties [23]; and h is the pressure head (m), namely, (m). Andwhere is the saturated volumetric content of water and is the volumetric content of residual water. Also,where I is the impedance coefficient. According to Taylor’s model [6, 14], .
In summary, the complete hydrothermal coupling analysis can be achieved by considering equations (2), (4), (9), and (11)–(16), and then, the thermalhydromechanical analysis can be obtained by solving simultaneously equations (5) and (10).
3. Application of the Developed Model
3.1. Application to the Roadbed
Currently, in frost heaving models, the boundary temperatures are usually simplified in practical applications by a symmetric boundary condition. Considering a section of the HarbinDalian HSR as an object, the thermalhydromechanical relation of the HSR roadbeds in seasonally frozen regions is investigated using field data and with the sunnyshady slope effect taken into consideration, the effect of which was previously reported by Zhang et al. [24] in the northsouth oriented part of the HarbinDalian HSR.
Steadystate monitored field data were obtained, and the results were used for a transient analysis of the initial temperature field. Once the initial temperature field is established, the thermalhydromechanical coupling is divided into two steps. First, the hydrothermal coupling analysis of the roadbeds is performed by secondary development of the mathematical module in COMSOL, to obtain the distributions and variations of moisture and temperature fields. Second, stress and deformation fields are calculated by the geomechanics module in COMSOL, which is based on the proposed model and the results of the hydrothermal coupling analysis.
As shown in Figure 2(a), the target section is located near the Shuangcheng Station in Harbin and its altitude is 166.4 m above mean sea level [25]. The left part of Figure 2(b) shows the composition of each roadbed component [25]. The width and the height of the roadbed are 13.6 m and 5.433 m, respectively, and the roadbed upper layer is a 0.4 m layer of wellgraded broken stone. A composite geomembrane insulating layer is laid just below the upper layer with a 0.05 m medium and coarse sand cushion below and above, respectively. The part under the insulating layer consists of a 1.0 m antifreezing layer with coarsegrained fillings (High quality and good filling materials in China railway roadbed design code [26]), while the lowest layer is composed of normal coarsegrained fillings. An antifrost heaving berm is set at the foot of the embankment slope. The foundation soil is clayey loess, silty clay, and sand, reinforced by CFG (cement flyash gravel) piles arranged in equilateral triangles. No surface water spillover was observed at the test site. The groundwater type is Quaternary porous phreatic water, which is mainly recharged by atmospheric precipitation and exposed to severe seasonal variation. The water depth is relatively large (between 13 m and 18.5 m), and the groundwater level has an amplitude approximately between 1 m and 2 m.
(a)
(b)
The right part of Figure 2(b) shows the arrangement of the roadbed measuring points considered [25]. Slope foots and temperature measuring holes are only arranged on the west part of the roadbed and other measuring points are arranged symmetrically with respect to the central line of the roadbed. Temperature measuring points are arranged at intervals of 0.5 m between 0.8 m and 3.8 m at intervals of 1 m between 3.8 m and 10.8 m below the roadbed surface. Five measuring lines are arranged at the east shoulder, west shoulder, central line, west slope foot, and 20 m away from the west slope foot, respectively, and each measuring line consists of 14 measuring points.
With the thermistor as the main component, the temperature sensor has a measuring range of −40°C to +60°C. Its accuracy is ±0.03°C in the range of −20°C to +20°C.
At the target section, wellgraded broken stone and tracks were laid in August 2010 and SeptemberOctober 2010, respectively. The temperature was measured every five days from August 01, 2010, to August 01, 2014. The ice density and water density were 917 and 1000 kg m^{−3}. The ice thermal conductivity and the water thermal conductivity were 2.22 and 0.58 W m^{−1}·°C^{−1}, respectively. The latent heat of the ice L was 334.56 kJ kg^{−3}. And the specific heat capacity of ice and the specific heat capacity of water were 2.09 and 4.18 kJ kg^{−1}·°C^{−1}. Tables 1 and 2 summarize the basic characteristic parameters of the roadbed and foundation layers [22, 27].


Considering the effects of temperature on the elastic modulus and Poisson’s ratio, the mechanical parameters of roadbeds are determined by the method reported elsewhere (see Table 3) [24].

By fitting a function of the ground temperature obtained from field monitoring data, as shown in Figure 3, the Dirichlet temperature boundary conditions considering the sunnyshady slopes effect for different locations are designed [27]:
The initial phase angle ω is determined according to the initial data obtained from the monitored area. Based on the boundary layer principle (boundary layer thickness was considered as 0.5 m) [24], the average ground temperatures at the east shoulder, central line, and west shoulder of the roadbed were determined to be 9.55°C, 7.71°C, and 7.80°C, respectively. is the temperature amplitude of the boundary layer bottom, is the time from the initial moment, measured in hours, and is the vibration period, adopted as 8760 hours. is the thermal gradient of the boundary layer bottom temperature increased by global warming, adopted as 0.04°C/8760 h [28], n is the outer normal direction of the boundary, and is the temperature at the lower boundary of the calculated region.
The dimensions of the model are shown in Figure 3 [27]. The height of clayey loess, silty clay, and sand was determined to be 11 m, 11 m, and 3 m, respectively, by drilling.
East shoulder to central line (the depth from the roadbed surface z range = 0–5.433 m):
Central line to west shoulder (the depth from the roadbed surface z range = 0–5.433 m):where x is the distance from the east shoulder (0–13.8 m).
To discretize the variations of the moisture and temperature fields of roadbeds in different seasons, the total time step and substep of the hydrothermal coupling transient analysis considered are 365 days and 0.5 days, respectively. The starting date of the transient analysis is August 1 (temperature peak moment).
3.2. Results and Discussion
3.2.1. Model Verification
Figure 4 shows a comparison between the simulated frost depth (SFR) and the monitored frost depth (MFR) of the west shoulder and a comparison between the simulated vertical displacement (SVD) of the roadbed and the monitored vertical displacement (MVD) within 2.7 m from the adjacent roadbed section [29]. As it is evident, simulated and monitored frost depths follow consistent trends. Simulation of the deformation demonstrates a reasonable match to the monitored deformation basically. The region with frost depth less than 1.2 m is exposed to severe soil frost heaving. When the frost depth exceeded 1.2 m, frost heaving developed slowly even though the frost depth increased further. The duration of thawing subsidence was short. As shown in Figure 4, field monitoring of the adjacent roadbed section indicated 2 mm residual deformation after thawing subsidence of the roadbed [29], while previous studies claimed that roadbed deformations almost returned to the initial condition after thawing subsidence. In some cases, a small amount of settlements were observed after thawing subsidence [27].
Therefore, from the results, it can be concluded that the proposed numerical model is effective and applicable for a comprehensive analysis of the temperature field, moisture field, and deformation field of the roadbed.
3.2.2. Analysis of the Simulation Results
Figure 5 illustrates the temperature field distribution of the section at certain moments in time. As observed, the freezing process is unidirectional and lasted for a long period, while the thawing process is bidirectional and has a short period. Due to the sunnyshady slope effect, the temperature fields of the roadbed are asymmetrically distributed in the transverse direction during the freezingthawing process. The maximum frost depth occurred at the west shoulder (right shoulder, sunny shoulder), followed by the east shoulder (left shoulder, shady shoulder) and that at the central line. The freezing sequence is west side slope, east side slope, and top surface. During thawing, frozen soil nuclei are observed in the roadbed and the thawing sequence is west side slope, east side slope, and roadbed top surface. Additionally, thawing of shoulders occurred earlier than that of the slope foot.
(a)
(b)
(c)
(d)
(e)
Figure 6 illustrates the distribution of the mass content of ice at the roadbed section. Due to the presence of the temperature field, icecontaining soil layers are observed at depths less than 1.2 m. At the same moment, the distribution of ice content of the roadbed section is uneven; the ice content of the west part is relatively large during freezing and melted rapidly during thawing.
Figure 7 shows time history curves of the temperature field and moisture fields of the central line of the roadbed. As observed, the icecontaining soil layer at the central line has a thickness of about 1.2 m, mainly within the range 0 to 0.6 m. Moisture migration is observed in both positive and negative temperature periods, but most especially in the negative temperature period. The moisture migration is concentrated in the region at depths less than 2.5 m. Here, the moisture distribution in positive and negative temperature periods is mainly affected by selfequilibrium of the moisture and temperature fields, respectively. During freezing, the mass content of moisture in the roadbed increased with depths less than 1.2 m decreased with depths between 1.2 m and 2.5 m and remained almost constant with depths greater than 2.5 m. Here, the mass content of moisture increased near the top surface of the roadbed and within the layer at a depth of 1 m. The maximum variation of the mass content of moisture is 2%, which was about 20% higher than that before freezing. In other words, moisture migration has a significant effect. The critical depth of moisture variation is 1.2 m. At this depth, the variation range of temperature is −0.7∼0°C (with the phase transition temperature zone) during freezing. Also, a frost front is presented, resulting in the migration of moisture to the lower part of the roadbed from an unfrozen region to a frozen one.
(a)
(b)
(c)
(d)
The distribution of the deformation of roadbeds is investigated based on the above analysis. In this work, the vertical displacement of the roadbed surface was the major concern. Figure 8 shows the time history curve, temperature curve, frost depth curve of the roadbed deformations. Figure 9 shows vertical displacements of roadbed section at different moments. Figure 10 shows the horizontal and vertical displacement fields of the roadbed section on March 10, 2014. As it is evident, vertical displacement of the roadbed surface followed the sequence of west shoulder, east shoulder, and central line, while the vertical displacement followed the sequence of west shoulder which is higher that of the east shoulder which is higher than that at the central line. The maximum vertical displacement of the roadbed surface is approximately 11∼18 mm (15 mm is allowed according to Chinese railway roadbed design code). At the early stage of freezing, the vertical displacement increased with the frost depth. Once the frost depth exceeded a critical level, the vertical displacement did not increase further with the frost depth. The roadbed vertical displacements are asymmetrically distributed in the transverse direction and also an uneven deformation of the roadbed surface appeared. In this study, the left and right shoulders are exposed to ambient conditions, while the central line is covered by the track slab. Hence, roadbed deformations at the shoulders and the central line are uneven in the early negative temperature period (up to 14 mm) and the difference subsequently decreases.
(a)
(b)
3.2.3. Discussion
For frost heaving of roadbeds in cold regions, changes of the temperature, moisture, and deformation are the cause, the medium, and the result, respectively. Roadbed deformations are macroscopic effects of hydrothermal migrations in roadbed soils.
Field monitoring revealed that frost heaving of the HarbinDalian HSR roadbed is present from the surface to a frost depth of 1.2 m. Here, frost heaving of the wellgraded broken stone dominated the frost heaving of the roadbed to a frost depth of 0.5 m [27]. Additionally, the moisture content varied drastically within the depth of 0.6 m.
The results of this study demonstrated that moisture migration at the freezing front is observed in the roadbed of the HarbinDalian HSR under the driving force of the temperature field. Here, the moisture content increased significantly in the range between 0 m and 1.2 m, and the icecontaining soil layers are concentrated in this region especially in the range between 0 and 0.6 m. Therefore, phase transition and moisture migration in HarbinDalian HSR roadbed driven by an uneven distribution of the temperature field are the main causes of roadbed deformations.
The transverse variation of the ground temperature of roadbeds can be attributed to the differences in the solar radiation intensity and also the surface covering layer [24, 30]. These variations of the ground temperature lead to an unbalanced development of the frost depth time and value, thus resulting in an unbalanced moisture migration. In this case, uneven transverse deformations are inevitable [30, 31]. Field monitoring data revealed differences in deformations of the east and west shoulders of the HarbinDalian HSR roadbed [32]. Simulation results demonstrated that the distribution of deformation field of roadbed in the transverse direction is asymmetric and uneven due to the effects of temperature and mass content of ice. Indeed, the deformation is maximized at the shoulders. The deformation of a specific part of roadbed is affected by a variation of the ground temperature and moisture at this part and of the surrounding soil. Hence, roadbed deformations are highly complicated. Frost heaving was first observed at the west shoulder, since freezing was first observed over there. However, deformation was initiated before freezing, since the freezing of the side slope occurred earlier than roadbed freezing. The uneven transverse deformation of the roadbed results in severe frost heaving at the shoulders and side slope. In fact, slope foots should be protected as weak spots due to severe ice accumulation results in horizontal displacement.
The results demonstrated that a reduction of the roadbed freezethaw deformations should be achieved based on the analysis of thermohydro migration characteristics of roadbeds, instead of a thermal insulation layer on the roadbed surface. Elimination of the uneven transverse deformations should be obtained by minimizing the difference of the ground temperature in the transverse direction, which can be achieved by asymmetric transitional layout of roadbed soils [30] or by an asymmetric antifrost heaving berm [33]. Based on a theoretical study, this study provides references to the design and maintenance of HSR roadbeds in seasonally frozen regions.
4. Conclusions
With SWCC, solidliquid ratio, and the linear thermal expansion coefficient as coupling equations, considering the influence of the transformation of ice into water and water migration on the temperature field, a thermalhydromechanical model applicable for unsaturated soils was proposed. The hydrothermal characteristic parameters of frozen and unfrozen soils were expressed by step functions. Then the thermalhydromechanical model is introduced in this study to explore the thermal field, water migration, and deformation variations of the HSR roadbed.
From the work described in this article, it is concluded that(1)The physical property fields of the roadbed varied significantly in soil layers with depths less than 2.5 m. Redistribution of the moisture field was triggered by a variation of the temperature field. Due to the significant effect of moisture migration, the mass content of moisture increased within the roadbed at depths greater than 1.2 m and the maximum variation of the mass content of moisture was 2 %. Most icecontaining soil layers were at depths less than 0.6 m.(2)Affected by the boundary of the sunnyshady slopes, distributions of the temperature field, moisture field, and deformation of the roadbed in the transverse direction were asymmetric and uneven. The maximum frost depth was greatest at the west shoulder (right shoulder, shady shoulder) followed by the east shoulder (left shoulder, sunny shoulder) and least at the central line. The maximum uneven settlement between the shoulder and the central line was 16 mm, while that between the centre of the two tracks was 4 mm. Here, the maximum vertical displacements of the west shoulder, east shoulder, and central line were 18 mm, 16 mm, and 12 mm, respectively.(3)The roadbed deformations are affected by the distribution and variation of the temperature and moisture fields. Therefore, reduction of freezethaw deformations and uneven deformations of roadbeds should be achieved by controlling the roadbed thermohydro migration characteristics.
Abbreviations
:  The volumetric content of unfrozen water 
t:  Time (s) 
:  The ice density (kg/m^{3}) 
:  The water density (kg/m^{3}) 
:  The volumetric content of ice 
K ():  The soil permeability coefficient related with (m/s) 
:  The soil water potential 
D ():  The water diffusion coefficient in soil (m^{2}/s) 
C_{m}:  The specific water capacity (m^{−1}) 
:  The matric potential (m) 
C:  The soil volumetric heat capacity (kJ·m^{−3}·°C^{−1}) 
T:  Temperature (°C) 
:  The volumetric heat capacity of water (kJ·m^{−3}·°C^{−1}) 
:  The thermal conductivity (W·m^{−1}·°C^{−1}) 
L:  The latent heat of the ice (kJ/kg) 
[D]:  The elasticity matrix 
[]:  The strain matrix 
:  The total strain matrix 
[]:  The thermal strain matrix 
[]:  The volumetric force 
E:  The elastic modulus (Pa) 
:  Poisson’s ratio 
:  Components of the normal strain of thermal expansion in x, y, and z directions 
α:  Linear thermal expansion coefficient (1/°C) 
:  Reference temperature (°C) 
:  Components of the shear strain of thermal expansion in x, y, and z directions 
:  The initial gravimetric moisture content 
:  The mass content of ice 
:  The mass content of unfrozen water 
:  The soil dry density 
b:  An empirical constant related to soil texture 
:  The freezing temperature of water (°C) 
:  The solidliquid ratio 
:  The linear thermal expansion coefficient of the wellgraded broken stone layer (1/K) 
:  The linear expansion coefficient (1/K) of coarsegrained fillings 
:  The specific heat capacity of soil skeleton (kJ·kg^{−1}·°C^{−1}) 
:  The specific heat capacity of water (kJ·kg^{−1}·°C^{−1}) 
:  The specific heat capacity of ice (kJ·kg^{−1}·°C^{−1}) 
:  The specific heat capacity of frozen soil skeleton (kJ·kg^{−1}·°C^{−1}) 
:  The specific heat capacity of thaw soil skeleton (kJ·kg^{−1}·°C^{−1}) 
H (, ):  The Heaviside step function 
:  The frozen soil thermal conductivity (W·m^{−1}·°C^{−1}) 
:  The thaw soil thermal conductivity (W·m^{−1}·°C^{−1}) 
:  The ice thermal conductivity (W·m^{−1}·°C^{−1}) 
:  The water thermal conductivity (W·m^{−1}·°C^{−1}) 
:  The soil permeability (m/s) 
:  The frozen soil permeability (m/s) 
:  The thaw soil permeability (m/s) 
:  The saturated permeability coefficient (m/s) 
:  The relative permeability 
o, m, k:  Empirical constants related to soil properties 
h:  The pressure head (m) 
:  The saturated volumetric content of water 
:  The volumetric content of residual water 
I:  The impedance coefficient 
:  The average ground temperatures 
:  The temperature amplitude of the boundary layer bottom 
P:  The vibration period, adopted as 8760 hours 
:  The thermal gradient of the boundary layer bottom temperature increased by global warming, adopted as 0.04°C /8760 h 
n:  Outer normal direction of the boundary 
:  Temperature at the lower boundary of the calculated region 
ω:  The initial phase angle. 
Data Availability
The data used to support the findings of this study are available from the corresponding author upon request.
Conflicts of Interest
The authors declare that they have no conflicts of interest in the submission of this manuscript. They do not have any commercial or associative interest that represents conflicts of interest in connection with the work submitted.
Acknowledgments
This research was funded by the National Natural Science Foundation of China (NSFC) under grant no. 51708369, National Key Research and Development Program of Ministry of Science and Technology under grant no. 2019YFC150960402, Natural Science Foundation of Hebei Province under grant no. E2017210110, Science Foundation of Chinese Postdoctoral under grant no. 2018M633607, Science and Technology Innovation Project of National Energy Investment Group under grant no. SHGF1541, and the Special Support Project of Hebei Province Basic Research Team under grant no. 31100803. The authors appreciate the help from Niu Fujun for supplying the deformation data monitored in the field.
References
 D. Cai, “Test on frost heaving spatialtemporal distribution of highspeed railway subgrade in seasonallyly soil region,” China Railway Science, vol. 37, pp. 16–21, 2016, in Chinese. View at: Google Scholar
 G. Shi, S. Zhao, and X. M. Li, “The frost heaving deformation of highspeed railway subgrades in cold regions: monitoring and analyzing,” Journal of Glaciology and Geocryology, vol. 36, pp. 360–368, 2014. View at: Google Scholar
 T. Wang and Z. Yue, “Influence of fines content on frost heaving properties of coarsegrained soil,” Rock and Soil Mechanics, vol. 34, pp. 359–364, 2013, in Chinese. View at: Google Scholar
 Q. Wang, J. Liu, Y. Tian, and X. Zhu, “A study of orthogonal design tests on frostheaving characteristics of graded crushed rock,” Rock and Soil Mechanics, vol. 36, pp. 2825–2830, 2015, in Chinese. View at: Google Scholar
 R. L. Harlan, “Analysis of coupled heatfluid transport in partially frozen soil,” Water Resources Research, vol. 9, no. 5, pp. 1314–1323, 1973. View at: Publisher Site  Google Scholar
 G. S. Taylor and J. N. Luthin, “A model for coupled heat and moisture transfer during soil freezing,” Canadian Geotechnical Journal, vol. 15, no. 4, pp. 548–555, 1978. View at: Publisher Site  Google Scholar
 J.M. Konrad and N. R. Morgenstern, “A mechanistic theory of ice lens formation in finegrained soils,” Canadian Geotechnical Journal, vol. 17, no. 4, pp. 473–486, 1980. View at: Publisher Site  Google Scholar
 J. M. Konrad and N. Morgenstern, “The segregation potential of a freezing soil,” Canadian Geotechnical Journal, vol. 18, pp. 482–491, 2011. View at: Google Scholar
 K. Hu, “Development of separated ice model coupled heat and moisture transfer in freezing soils,” 2011, Ph.D. thesis, China University of Mining and Technology, Xuzhou, China, p. 130, in Chinese. View at: Google Scholar
 R. D. Miller, “Freezing and heaving of saturated and unsaturated soils,” Highway Research Record, vol. 393, pp. 1–11, 1972. View at: Google Scholar
 K. O’Neill and R. D. Miller, “Exploration of a rigid ice model of frost heave,” Water Resources Research, vol. 21, pp. 281–296, 1985. View at: Google Scholar
 Y. Zhou, “Study on frost heave model and frost heave control of frozen soils,” Journal of China University of Mining and Technology, p. 113, 2009, in Chinese. View at: Google Scholar
 C. Duquennoi, M. Fremond, and M. Levy, “Modelling of thermal soil behavior,” VTT Symposium, vol. 95, pp. 895–915, 1989. View at: Google Scholar
 X. Mao and B. Ma, Studies on the Stability of Permafrost Subgrade Based on Coupled Water and Heat Transfer, China Communications Press, Beijing, China, 2011, in Chinese.
 J. Wu and T. Han, “Numerical research on the coupled process of the moistureheatstress fields in saturated soil during freezing,” Engineering Mechanics, vol. 26, pp. 246–251, 2009, in Chinese. View at: Google Scholar
 Q. Xu, G. Peng, N. Li, and Z. Liu, “Numerical method of phasechange field of temperature coupled with moisture, stress in frozen soil,” Journal of Tongji University, vol. 33, pp. 1281–1285, 2005, in Chinese. View at: Google Scholar
 Y. Lai, W. Pei, M. Zhang, and J. Zhou, “Study on theory model of hydrothermalmechanical interaction process in saturated freezing silty soil,” International Journal of Heat and Mass Transfer, vol. 78, pp. 805–819, 2014. View at: Publisher Site  Google Scholar
 F. Ming and D. Li, “Modeling and experimental investigation of one dimension coupled moisture and heat in unsaturated freezing soil,” Journal of Central South University (Science and Technology Edition), vol. 45, pp. 889–894, 2014, in Chinese. View at: Google Scholar
 J. Zhou and D. Li, “Numerical analysis of coupled water, heat and stress in saturated freezing soil,” Cold Regions Science and Technology, vol. 72, pp. 43–49, 2012. View at: Publisher Site  Google Scholar
 Q. Bai, “Determination of boundary layer parameters and a preliminary research on hydrothermal stability of subgrade in cold region,” 2016, M.S. thesis, Beijing Jiaotong University, Beijing, China, p. 73, in Chinese. View at: Google Scholar
 B. Tai, J. Liu, X. Li, Z. Yue, and Y. Shen, “Numerical model of frost heaving and anti frost heaving measures of high speed railway subgrade in cold region,” China Railway Science, vol. 38, pp. 1–9, 2017, in Chinese. View at: Google Scholar
 X. Xu, J. Wang, and L. Zhang, Physics of Frozen Soils, Science Press, Beijing, China, 2010, in Chinese.
 M. T. Van Genuchten, “A closedform equation for predicting the hydraulic conductivity of unsaturated soils,” Soil Science Society of America Journal, vol. 44, no. 5, pp. 892–898, 1980. View at: Publisher Site  Google Scholar
 Y. Zhang, Y. Du, and B. Sun, “Temperature distribution in roadbed of highspeed railway in seasonallyly frozen regions,” Chinese Journal of Rock Mechanics and Engineering, vol. 33, pp. 1286–1192, 2014, in Chinese. View at: Google Scholar
 Y. Zhang, Y. Du, and B. Sun, “Temperature distribution analysis of highspeed railway roadbed in seasonally frozen regions based on empirical model,” Cold Regions Science and Technology, vol. 114, 2015. View at: Publisher Site  Google Scholar
 China State Railway Administration (CSRA), Code for Design of Railway Earth Structure (TB 10001–2016), China Railway Publishing House, Beijing, China, 2014.
 Y. Zhang, Y. Du, B. Sun, S. Zhang, and J. Han, “Roadbed deformation of highspeed railway due to freezingthawing process in seasonallyly frozen regions,” Chinese Journal of Rock Mechanics and Engineering, vol. 33, pp. 175–182, 2014, in Chinese. View at: Google Scholar
 H. Jin, Q. Yu, L. Lü et al., “Degradation of permafrost in the Xing’anling mountains, northeastern China,” Permafrost and Periglacial Processes, vol. 18, no. 3, pp. 245–258, 2007. View at: Publisher Site  Google Scholar
 F. Niu, A. Li, J. Luo et al., “Soil moisture, ground temperatures, and deformation of a highspeed railway embankment in northeast China,” Cold Regions Science and Technology, vol. 133, pp. 7–14, 2016. View at: Google Scholar
 H. Liu, F. Niu, and H. Guan, “An engineering evaluation index of thermal asymmetry in subgrade and its optimal design in cold regions,” Cold Regions Science and Technology, vol. 137, pp. 1–6, 2017. View at: Publisher Site  Google Scholar
 Y. Zhang, Y. Du, B. Sun, S. Zhang, and J. Han, “Analysis on deformation characteristics of roadbed transition section of highspeed railway in seasonally frozen region,” China Railway Science, vol. 37, pp. 39–45, 2016, in Chinese. View at: Google Scholar
 H. Liu, F. Niu, Y. Niu, J. Xu, and T. Wang, “Effect of structures and sunnyshady slopes on thermal characteristics of subgrade along the HarbinDalian Passenger Dedicated Line in Northeast China,” Cold Regions Science and Technology, vol. 123, pp. 14–21, 2016. View at: Publisher Site  Google Scholar
 A. Wen, Y. Zhang, and W. Zhao, “The frost characteristics influence of highspeed railway roadbed with different antifrost berms,” Journal of Shijiazhuang Tiedao University, vol. 32, pp. 27–33, 2019, in Chinese. View at: Google Scholar
Copyright
Copyright © 2020 Yuzhi Zhang et al. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.