- Email: [email protected]
Roadbed temperature study based on earth-atmosphere coupled system in permafrost regions of the Qinghai-Tibet Plateau
Cold Regions Science and Technology 86 (2013) 167–176
Contents lists available at SciVerse ScienceDirect
Cold Regions Science and Technology journal homepage: www.elsevier.com/locate/coldregions
Roadbed temperature study based on earth-atmosphere coupled system in permafrost regions of the Qinghai-Tibet Plateau Zhi-Yun Liu a,⁎, Jian-Bing Chen b, Long Jin b, Yu-Jie Zhang a, Chao Lei a a b
College of Geology Engineering & Geomatics, Chang'an University, Xi'an, Shaanxi 710054, China CCCC First Highway Consultants Co. Ltd, Xi'an, Shaanxi 710054, China
a r t i c l e
i n f o
Article history: Received 5 July 2012 Accepted 10 October 2012 Keywords: Permafrost roadbed Temperature field Earth-atmosphere coupled system Coupled heat transfer Parameter analysis
a b s t r a c t In this work, a two-dimensional unsteady model for analyzing the roadbed temperature in permafrost regions of the Qinghai-Tibet Plateau (QTP) is developed and the heat transfer processes of comprehensive earth-atmosphere coupled system are simulated by the numerical method. In the computation, the surrounding air environment of ground surface is brought into the research model and the corresponding unsteady and non-uniform influence factors of earth-atmosphere coupled system are also considered to obtain reasonable simulation results. Meanwhile, the model is validated with monitored data at a depth of 0.5 m in different roadbed positions. Based on the above model, the spatial and temporal temperature distributions of the roadbed surface in Wudaoliang region are determined. Then the influences of the system parameters, such as the ambient temperature, wind speed, roadbed height and route trends, on the local temperature distribution are investigated. The results show that, temperature of surface regions exhibits a periodic sinusoid trend and always higher than the ambient air temperature. For different surface regions, the temperature of asphalt pavement keeps the maximum value and the natural ground surface has an abrupt step change in temperature at beginning of April and end of October. The variation of ambient air temperature has little influence on annual average temperature of permafrost roadbed surface. The deviations between the surfaces' annual average temperature increase and their average values are less than 1%. Annual average temperature of various surfaces decreases as the wind speed and the roadbed height increase. Furthermore, the annual average temperature of leeward slope (also sunny side) along EW direction is the highest one among all three route directions and such temperature decreased systematically from SW45° direction to SN direction. However, the annual temperature of upwind slope (also shady side) is just contrary to the leeward slope. © 2012 Elsevier B.V. All rights reserved.
Contents 1. 2.
3. 4.
Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Modeling of the permafrost roadbed . . . . . . . . . . . . . . . . . . . . . . 2.1. The physical and mathematical model of earth-atmosphere coupled system 2.2. Boundary condition and other computation settings . . . . . . . . . . . 2.2.1. Velocity and direction of the wind . . . . . . . . . . . . . . . . 2.2.2. Environmental temperature . . . . . . . . . . . . . . . . . . . 2.2.3. Source term of the couple wall . . . . . . . . . . . . . . . . . 2.2.4. Special source term settings of the slope couple wall . . . . . . . 2.3. Numerical procedures . . . . . . . . . . . . . . . . . . . . . . . . . . Validation of the numerical procedure . . . . . . . . . . . . . . . . . . . . . Results and discussions . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4.1. Numerical results of the permafrost roadbed . . . . . . . . . . . . . . . 4.2. Effects of parameter on the coupled heat transfer . . . . . . . . . . . . . 4.2.1. Annual average air temperature . . . . . . . . . . . . . . . . . 4.2.2. Annual average wind speed . . . . . . . . . . . . . . . . . . . 4.2.3. Roadbed height . . . . . . . . . . . . . . . . . . . . . . . . 4.2.4. Route direction . . . . . . . . . . . . . . . . . . . . . . . .
⁎ Corresponding author. Tel.: +86 29 82339356. E-mail address: [email protected] (Z.-Y. Liu). 0165-232X/$ – see front matter © 2012 Elsevier B.V. All rights reserved. http://dx.doi.org/10.1016/j.coldregions.2012.10.005
. . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . .
168 168 168 169 170 170 170 171 171 171 173 173 174 174 174 174 175
168
Z.-Y. Liu et al. / Cold Regions Science and Technology 86 (2013) 167–176
5. Conclusions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Acknowledgments . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . References . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
1. Introduction Temperature distributions of the Qinghai-Tibet Highway (QTH) roadbed in permafrost regions are influenced by many factors, such as climatic conditions, physical parameters and construction scheme. So the temperature field of a highway roadbed would be different from that of the natural ground and such temperature difference shall further influence the stability of underlying permafrost (Niu et al., 2002; Wu et al., 2005). Meanwhile, the temperature field of the highway roadbed is also the prior condition to study the interaction process between the QTH and the permafrost. It not only directly influences the heat and mass transfer process within the roadbed, but also influences morphology and the freeze–thaw process of the underlying permafrost. Therefore, it is necessary to carry out numerical studies about the temperature field of permafrost roadbed (Sheng et al., 2003; Wang et al., 2005; Wu et al., 2001). The heat transfer process of permafrost roadbed includes two aspects: i) unsteady heat transfer process between roadbed and external environment; and ii) unsteady heat transfer process between roadbed and the underlying permafrost. Unlike the main conductive heat transfer process between roadbed and the underlying permafrost, the heat transfer process between roadbed and ambient environment consists of different and coupled heat transfer types, which includes: the solar radiation absorption of highway surface, the long-wave radiation emitting from the highway surface to surrounding environment and the forced-convection between the air and the highway surface. Furthermore, the different pavement structures make the change of absorption coefficient and evaporation capacity, and the different route trends influence the absorption coefficient of the slope surface. And also, some influencing factors (such as air temperature, wind velocity and direction, solar radiation density, absorption coefficient, evaporation capacity, etc.) change with time, which all aggravate the complexity of heat transfer process numerical computation work on the solving of temperature field of permafrost roadbed (Liu et al., 2002; Wang et al., 2008; Zhang et al., 2004). Due to the above reasons, it is very hard to obtain the temperature field of permafrost roadbed directly by the numerical method and this problem has became a challenging research subject in permafrost engineering research (Lai et al., 2002; Zhang et al., 2004). In recent years, many significant numerical and experimental research works had been conducted on the temperature field of permafrost roadbed. Niu et al. (2008) investigated the surface temperature of pipe-ventilated roadbed by placing temperature measuring points and summarized its annual variation rules, which provided the boundary conditions for solving temperature field of permafrost roadbed based on boundary layer theory. Mi et al. (2002) and Su and Ning (2004) simulated the temperature field of pipe-ventilated roadbed using the constant temperature boundary condition and the results showed that the duct-ventilated structure can effectively cool the roadbed. Zhang et al. (2005) and Lai et al. (2004) compared the temperature distributions of different roadbed structures and the results showed that the ripped-rock embankment would be the better choice for the high temperature permafrost region. However, the above research works on the temperature field of permafrost roadbed are on the premise of knowing the detailed information of the surface temperature distribution, which will need a long time field-observation or a suitable empirical formula. In this way, the influences of a large number of coefficients such as air temperature, wind direction and air speed, etc., are neglected. Consequently, in the numerical results there exists a large deviation from the actual situation and with rigid regional limitations. So in this work, to reveal the spatial and temporal temperature distribution characteristics of
175 176 176
permafrost roadbed surface and to optimize the design of highway in permafrost region, the complex heat transfer processes of the comprehensive earth-atmosphere coupled system have been simulated. In the simulation, the surrounding air environment of the ground surface is brought into the research model and the environmental influencing factors, such as the solar radiation, air temperature, wind speed and direction, surface vegetation and evaporation, are also considered. Based on the above model, the spatial and temporal temperature distributions of the roadbed surface in Wudaoliang region are determined. Then the influences of the system parameters, such as the ambient temperature, wind speed, roadbed height and route trends, on the local temperature distribution are investigated. 2. Modeling of the permafrost roadbed 2.1. The physical and mathematical model of earth-atmosphere coupled system Owing to direct heat exchanging processes between roadbed and underlying natural permafrost, the temperature distribution conditions of underlying natural permafrost have great influence on the temperature field of roadbed since the engineering construction is finished. To obtain a more reasonable calculation result, the temperature distribution of underlying natural permafrost that is on the time node of construction work finished should be imported into the computation module as the initial field. So in this work, there are two corresponding physical models for numerical calculation: i) the initial field model with natural permafrost and air environment; and ii) the permafrost roadbed model with comprehensive earth-atmosphere coupled system. The schematic of the above two models is shown in Fig. 1. As shown in the figure, the underlying natural permafrost is composed of silty clay, gravelly clay and fully weathered mudstone and the thickness of each layer is 2 m, 5 m and 23 m from the top down. The highway engineering is built by subgrade soil and asphalt pavement. The width and thickness of pavement structure is designed as 7 m and 0.3 m respectively. In order to decrease the influence of top environment, the height of the air layer is chosen as 20 m. The longitudinal length of the initial field mode is defined as 54 m. To avoid the influence of entrance effect, the length of bilateral natural permafrost of the roadbed is defined as 60 m. The width of bilateral road shoulder is defined as 1.5 m. Meanwhile, the height and slope ratio of the roadbed is defined as 3 m and 1:1.5 respectively, and the route direction highway is chosen as NS direction. The Wudaoliang region is located at the 4000 m height of the Qinghai-Tibet Plateau, so the reference pressure is set as 60,000 Pa in the numerical simulation. Furthermore, considering that the variation range of the ambient air temperature is within 20 centigrade, the air properties are defined as a fixed value in the numerical simulation. The detailed physical properties of air used in this work are listed as following: ρ=0.7441 kg/m3, Cp =1005 J/(kg⋅K), μ=1.751×10−5kg/(m⋅s). And because the rough ground surface can only be defined as the smooth surface in the calculation software, so the conductivity of air (λ) is amplified as 0.0586 W/(m·K) to obtain a reasonable results. The density of pavement structure, subgrade soil and different natural permafrost layers is also defined as a fixed value. But their specific heat capacity and conductivity are set as the temperature correlation function, except for the fixed definition of pavement structure. The specific value of the above mentioned physical properties are listed in Table 1 (Xu et al., 2001).
Z.-Y. Liu et al. / Cold Regions Science and Technology 86 (2013) 167–176
169
(b) Permafrost roadbed model
(a) Initial field model
Fig. 1. Schematics of the two computation model.
A two-dimensional unsteady model has been developed in the present work to predict the spatial and temporal temperature distribution characteristics of Qinghai-Tibet plateau permafrost roadbed surface. The heat transfer processes between the ground surface and ambient environmental are coupled in the boundary condition settings and the source terms of energy equations. The vector form of the continuity, momentum and energy equations are expressed as: Continuity equation: → div U ¼ 0
ð1Þ
Momentum equation: → ∂ðρui Þ ∂p þ divðη gradui Þ þ div ρui U ¼ − ∂t ∂xi
ð2Þ
2.2. Boundary condition and other computation settings In this paper, the surrounding air environment is brought into the research model and computation model is modified from the closed system (based on the boundary layer theory) to the opening system with air inlet and outlet. Therefore, more influencing factors are considered in the numerical simulation and the corresponding boundary conditions become more complex. The detail boundary condition settings are listed as the following: i) both the left and right sides of the air layer are defined as the velocity inlet boundary conditions. The wind enters the computation area vertically to the boundary and its speed, direction and temperature change annually during a year. Except for changing with time, the wind speed is also set as the function of the height from the ground. All the boundary condition parameter settings of the left and right sides are same, but only the wind direction is set as the contrast value to keep the quantity conservation of system. And the turbulence density and turbulence length scale can be calculated by: 1 I ¼ 0:16 ReDH 8
Energy equation: ! → ∂ðρT Þ k gradT þ ST þ div ρ U T ¼ div cp ∂t
ð3Þ
ð4Þ
l ¼ 0:07DH
ii) top of the air layer is defined as the wall boundary condition and its temperature varies annually during a year; iii) the bilateral side of the
Table 1 Physical property of the soil layers. Physical properties ρ (kg ⋅ m−) Cp (J ⋅ kg−1 ⋅ K−1)
λ (W ⋅ m−1 ⋅ K1)
T (K)
Pavement
Roadbed
Silty clay
Gravelly clay
Mudstone
233.15 253.15 263.15 268.15 271.15 272.15 272.65 272.95 273.05 273.15 313.15 233.15 272.75 272.91 273.07 273.23 313.15
2177.78 815.78 815.78 815.78 815.78 815.78 815.78 815.78 815.78 815.78 815.78 815.78 1.55 1.55 1.55 1.55 1.55 1.55
2183.6 8.38E + 02 8.43E + 02 8.52E + 02 8.70E + 02 1.14E + 03 2.06E + 03 4.12E + 03 8.24E + 03 1.65E + 04 1.02E + 03 1.02E + 03 1.5 1.5 1.49 1.48 1.4 1.4
2080 1.05E + 03 1.06E + 03 1.20E + 03 1.44E + 03 2.88E + 03 6.01E + 03 1.20E + 04 2.40E + 04 4.94E + 04 1.17E + 03 1.17E + 03 1.8 1.8 1.78 1.73 1.5 1.5
2340 8.78E + 02 8.97E + 02 9.83E + 02 1.07E + 03 1.92E + 03 3.85E + 03 7.69E + 03 1.50E + 04 2.95E + 04 1.22E + 03 1.22E + 03 2.2 2.2 2.13 1.94 1.0 1.0
2310 2.64E + 06 2.70E + 06 2.80E + 06 3.00E + 06 4.00E + 06 8.00E + 06 1.70E + 07 3.30E + 07 6.56E + 07 2.97E + 06 2.97E + 06 2.50 2.50 2.47 2.39 2.0 2.0
170
Z.-Y. Liu et al. / Cold Regions Science and Technology 86 (2013) 167–176
Table 2 Average wind velocity of different months in the Wudaoliang region. Month
Jan.
Feb.
Mar.
Apr.
May
June
July
Aug.
Sep.
Oct.
Nov.
Dec.
v (m ⋅ s−1)
5.69
5.83
5.92
4.79
4.37
4.11
3.72
3.37
3.38
3.58
4.29
5.23
permafrost layers are all defined as the adiabatic boundary condition; iv) the bottom fully weathered mudstone layer is defined as the heat flux boundary condition with the specific value of 0.06 W/m2 (Zhang et al., 2005); v) the natural permafrost and road engineering surface are defined as the 0.1 m thick couple heat transfer wall. The above mentioned complex heat transfer processes between the ground surface and ambient environment are coupled in the source term computation of couple wall. In the following sections, the specific value of air temperature, wind speed and direction of different times as well as the computation method of couple wall's source terms will be given. 2.2.1. Velocity and direction of the wind According to the research, the wind velocity of atmospheric environment is the function of the height from the ground (Liu et al., 2009): v ¼ vref ðh=href Þ
ð5Þ
where, vref is the wind speed of reference height; href is the height reference spot and chosen as 3 m in this work. Meanwhile, the wind speed and direction also change periodically for years. The average wind speed of different months in the Wudaoliang region is listed in Table 2. To simplify the numerical computation, the wind direction is selected as the maximum wind frequency during a certain time. And based on the meteorological data, the wind direction is defined as from left to right in the time range of Feb. 20 to Aug. 20 and from right to left in the other time of year. 2.2.2. Environmental temperature The ambient air temperature changes periodically during a year and its calculation function is set as the sinusoid in this paper: 2πt 2π þ 3600 24 30 12 3
6
4
3
2
þ 0:4266d −0:1875d−1:455d þ 135:8
ð8Þ
where d is the period order of ten days during a year (since Aug. 20) and its value range from 1 to 36 (assuming 30 days a month). Because the local ground surface temperature is different for the surface position and changes with times, so the value of Qrad also changes with time and local surface temperature: 4 4 Q rad ¼ εσ T local −T sky
ð9Þ
where ε is the emissivity of ground surface and defined as 0.9; Tlocal is the temperature of local ground surface region and can be derived through the numerical calculation for every time step; Tsky is the sky temperature and can be calculated as (Duffie and Beckman, 1996): 1:5
T sky ¼ 0:0552T air
ð10Þ
ð6Þ Qeva can be calculated as:
In this equation, t is the influence time since the construction work of highway is finished (Aug. 20); T0 is the mean annual temperature of Wudaoliang region and its specific value is defined as −5.3 °C based on the local meteorological data; g(t) is the annual temperature ascending speed and defined as 0.022 °C/year (Wu et al., 2005). Because it is assumed that all the parameters keep a constant value or vary periodically in the initial field model, so the term of g(t) is not included in the air temperature computation equation in such model. 2.2.3. Source term of the couple wall The source term of couple wall is the net energy absorbed by ground surface and it can be obtained by analyzing the heat transfer process between roadbed and ambient environment: Q sou ¼ αQ sol −Q rad −Q eva
7
Q sol ¼ 3:214e−7d −4:111e−5d þ 0:002025d −0:04629d
0:14
T air ðt Þ ¼ T 0 þ g ðt Þ þ 11:3 sin
road shoulder are both 0.72 and pavement is 0.9. However, the absorption coefficient of the natural ground surface varies with time for the ground vegetations. Considering that Wudaoliang belongs to the typical desert steppe region, the absorption coefficient of the natural ground surface is simplified as changing with seasons and its specific value is listed in Table 3. The solar radiation intensity casting on the ground surface includes the direct normal solar radiation and diffuse solar radiation and its value changes periodically during a year in generally. For lack of comprehensive meteorological data, the solar radiation intensity data used in the numerical computation are based on the whole Qinghai-Tibet Plateau solar radiation data (Zhong, 1989) and obtained by the following time correlation function:
ð7Þ
where α is the absorption coefficient of the ground surface; Qsol is the total solar radiation projected on the ground surface; Qrad is the long-wave radiation heat loss from ground surface to environment; Qeva is the heat loss of surface evaporation. Because the specific values of the above terms in Eq. (7) fluctuate for different surface types and some of them change with time, so the source terms of the couple wall should be determined according to the surface type. The detailed computation methods of every term in Eq. (7) are accounted below. The absorption coefficients of the slope, road shoulder and pavement surface are defined as the fixed values, where the slope and
Q eva ¼ UG
ð11Þ
where G is the latent heat of vaporization of water and defined a fixed value of 2500 kJ/kg; U is the month surface evaporation quantity and can be calculated as (Wang, 2000): U ¼ 0:3082⋅uw −0:849
ð12Þ
and uw is the month water surface evaporation of the Wudaoliang region, which can be obtained by inquiring the hydrology data (see Table 4). Meanwhile, considering the climatic condition of the Wudaoliang region, Qeva has a certain value only from April to November in this paper for the ground surface is exposed to the ambient environment and reduces to Table 3 Absorbing coefficient of the natural ground surface under different seasons. Time range
Spring
Summer
Autumn
Winter
Alpine steppe
0.644
0.719
0.691
0.668
Table 4 Water surface evaporation of different months. Month
Apr.
May
June
July
Aug.
Sep.
Oct.
Nov.
uw (mm)
114
156
158
167
161
112
79
68
Z.-Y. Liu et al. / Cold Regions Science and Technology 86 (2013) 167–176 Table 5 Slope coefficient of different months and different route directions. Month
Jan. Feb. Mar. Apr. May June July Aug. Sep. Oct. Nov. Dec.
NS
SW45°
EW
Upwind
Leeward
Upwind
Leeward
Upwind
Leeward
0.958 0.946 0.934 0.931 0.928 0.922 0.928 0.931 0.934 0.946 0.958 0.97
0.958 0.946 0.934 0.931 0.928 0.922 0.928 0.931 0.934 0.946 0.958 0.97
0.58 0.658 0.76 0.82 0.886 0.904 0.886 0.82 0.76 0.658 0.58 0.554
1.3 1.18 1.066 1.0 0.94 0.916 0.94 1.0 1.066 1.18 1.3 1.366
0.418 0.52 0.67 0.79 0.88 0.898 0.88 0.79 0.67 0.52 0.418 0.4
1.414 1.24 1.126 1.03 0.992 0.898 0.992 1.03 1.126 1.24 1.414 1.498
2.2.4. Special source term settings of the slope couple wall Compared to the surfaces of the natural ground, road shoulder and pavement, the slope surface has different solar radiation absorption quantities for the different route directions and its angle factor does not keep to 1.0 for the slope inclination. Therefore, to reflect the effects of the route direction and angle factor, the calculation equation of the slope couple wall's source term has been modified as: ð13Þ
where k and X is the slope coefficient and angle factor respectively. Because the direct solar radiation accounts for 60% and diffuse solar radiation is 40% in the total solar radiation of Qinghai-Tibet Plateau (Qiao et al., 2004), so the integral slope coefficient used in present work is defined as: κ ¼ 0:6κ dir þ 0:4κ dif
ð14Þ
where κdif is the slope coefficient of diffuse solar radiation and defined as 1.0 for its non-directional characteristic; κdir is the slope coefficient of direct solar radiation and also different from the sunny slope and shady slope as well (Wang, 2000):
κ dir−shady ¼ ð cot α s − cot βj cosðδ−θÞjÞ sin α s κ dir−sun ¼ ð cot α s − cot βj cosðδ−θÞjÞ sin α s
above equations, it can be inferred that the slope coefficient not only changes with time but also varies with the route direction. To simplify the numerical computation, the slope coefficient is assumed as the month of piecewise function and the detailed slope coefficients of different months (from June to Dec.) and different route directions (NS, SW45° and EW) used in this work are listed in the Table 5. The angle factor, X, can be obtained through the radiation heat transfer related theory formulas (Tao, 2001). Under the physical model used in this work, the specific value of the angle coefficient is calculated as 0.852. 2.3. Numerical procedures
0 in other times for the snow covering. And the detailed time of snow melting and freezing is judged by whether the ground surface temperature is larger than 0 °C or not. Furthermore, owing to the structure of the asphalt pavement, the heat loss carried away by surface vaporization is not included in its source term computation.
Q sou ¼ αkQ sol −XQ rad −Q eva
ð15Þ
where αs is subgrade slope angle; β is sun elevation angle; θ is sun azimuth angle; and δ is the angle between east and north. Under the
The geometry and quadrilateral mesh of the permafrost roadbed are created by GAMBIT 2.3.16 software, as shown in Fig. 2. To get better numerical results, local grid refinement is applied on the roadbed slope, road shoulder and pavement and the final ratio between minimal grid scale and maximal grid scale is 1/15. The governing equations for the forced convection and heat conduction are solved in the FLUENT 6.3 software package. The 2-D, unsteady state and turbulence model with implicit solver is adopted in the simulation. In the boundary condition settings, except for the insulation boundary condition of the bilateral side of permafrost layers, other complex boundary conditions and source term of the couple wall are imported into the FLUENT by UDF program. The properties of working fluid are chosen as the air at 273.15 K and the physical properties of different soils are taken as the step function of temperature. Meanwhile, in order to embody the influence of periodically changing environment on the coupled heat exchange of whole computation model, ambient temperature, wind speed, solar radiation and evaporation capacity are set as the annual periodical change function. The SIMPLE algorithm has been used for pressure velocity coupling. Standard scheme is adopted for the discretization of the pressure term and the second order upwind scheme is applied to discretize other terms. A convergence criterion of 10 −6 is imposed on the residuals of the momentum, energy and 10−5 on continuity the equations. The time step is defined as 12 h and the inner iterations of each time step are chosen as 200. Besides, the calculation time of initial field model is defined as 100 years and then its calculation results are imported into the permafrost roadbed model as the initial field of non-steady state computation. 3. Validation of the numerical procedure In order to validate the present numerical procedure, the comparative studies have been carried out for the temperature monitoring data and numerical calculation result of No.7 Tangbei permafrost roadbed in Wudaoliang. In the numerical computation, the prototype of the computation model is based on the geological section map of the Wudaoliang region in “Frozen Soil Engineering Geology Section for
20
Y /m
10 0 -10 -20 -30 -30
-20
-10
0
10
20
30
40
171
50
60
70
80
90
X /m Fig. 2. Computation grid of the permafrost roadbed.
100
172
Z.-Y. Liu et al. / Cold Regions Science and Technology 86 (2013) 167–176
4.0 0.0
Y /m
-4.0 -8.0
-12.0 -16.0
30
20
10
0.0
10
20
30
X /m Fig. 3. Geological section map of highway in Wudaoliang region.
(a) Natural ground surface
(b) Pavement
6
15 Numerical results
Experimental results
2 0 -2 -4 -6 -8
-10 -12
9 6 3 0 -3 -6 -9
-12 2008/2/20
2008/8/21
2009/2/20
2009/8/20
2010/3/18
-15
2010/9/21
Date / yy-mm-dd
2008/8/21
2009/2/20
2009/8/20
2010/3/5
2010/9/5
--
(d) Leeward road shoulder Temperature of upwind road shoulder / ºC
12 Experimental results
2008/2/18
Date / yy-mm-dd
(c) Upwind road shoulder Temperature of leeward road shoulder / ºC
Numerical results
12
4
Temperature of pavement / ºC
Temperature of natural surface / ºC
Experimental results
Numerical results
9 6 3 0 -3 -6 -9
15 12
Experimental results
Numerical results
9 6 3 0 -3 -6 -9
-12
-12
-15
2008/2/18 2008/8/21 2009/2/20 2009/8/20 2010/3/5 2010/9/5
Date / yy-mm-dd
--
2008/2/18 2008/8/21 2009/2/20 2009/8/20 2010/3/5 2010/9/5
Date / yy-mm-dd
Fig. 4. The temperature distributions of 0.5 m depth soil under different positions and times.
--
Z.-Y. Liu et al. / Cold Regions Science and Technology 86 (2013) 167–176
173
10
v /m · s-1
hlocal /W · m-2
Tlocal / ºC
Local ground surface temperature
Toe of slope
8 6 4 2 32 28 24 20 16 12 8 8
Local heat transfer coefficient
Right road shoulder
Left road shoulder Wind velocity of 0.3m height
6 4 2
Toe of slope
0
0
20
40
60
80
100
120
140
Distance of X direction / m Fig. 5. Wind speed of 0.3 m height, local heat transfer coefficient and temperature of different ground surface positions along the horizontal direction on March 20 of 20th year.
Qinghai to Tibet highway” (as shown in Fig. 3), which is drawn by Cold and Arid Regions Environmental and Engineering Research Institute of CAS and CCCC First Highway Consultants Co., Ltd. in 1999. The detailed computation parameters are the same as the predefined parameters of permafrost roadbed model and according to the 30 years meteorological data of Wudaoliang region (from 1971 to 2000), the mean annual air temperature is defined as −5.3 °C. The temperature distribution of soil with 0.5 m depth of different positions (natural permafrost, pavement, upwind and leeward road shoulder) and different times (from 21th year to 23th year since the highway engineering finished) are shown in Fig. 4, and it can be seen that the numerical results are in good agreement with the drilling measured data.
15 10
T / ºC
5 0 Upwind natural surface Upwind slope Pavement Leeward slope Leeward natural surface Enviroment
-5 -10 -15 0
30
60
90 120 150 180 210 240 270 300 330 360
Times / Day Fig. 6. The yearly surface temperature variation of different ground positions since the highway engineering has been established for 20 years.
4. Results and discussions 4.1. Numerical results of the permafrost roadbed Based on previous mathematical model and calculation method, the numerical simulation is carried out for complex coupled heat transfer within permafrost roadbed of Wudaoliang region. The wind speed at the height of 0.3 m (v), local heat transfer coefficient (hlocal) and temperature (Tlocal) of different positions along the horizontal direction on March 20 of 20th year (since the highway engineering finished) are shown in Fig. 5. It can be seen that owing to the shading effect of roadbed, the wind speed and local heat transfer coefficient decrease but the local surface temperature increases near the upwind slope; then the wind speed continue to ascend in the upwind slope and reaches its maximum value in the top of slope, and corresponding local heat transfer coefficient and surface temperature increase for the wind enhanced heat transfer; when the air flow reaches the leeward slope region, the wind speed and local heat transfer coefficient decrease and the local surface temperature increases firstly, then as the dropping circumfluence phenomenon appears in the air flow and strengthens the ground heat transfer, local heat transfer coefficient increases and surface temperature decreases rapidly; the wind speed of leeward natural ground is smaller than that of the upwind, however, due to the air flow backwash, the local heat transfer coefficient of leeward is larger than that of the upwind and thus the surface temperature of leeward is less than that of upwind.
Table 6 Parameters of the different ground surfaces' temperature fitting equations. Fitting parameters
Upwind natural surface
Upwind slope
Pavement
Leeward slope
Leeward natural surface
a b c σ2
2.54 11.65 0.68 0.98
4.22 13.48 0.64 1
6.21 15.12 0.66 1
4.34 14.01 0.66 0.99
2.47 11.43 0.66 0.99
174
Z.-Y. Liu et al. / Cold Regions Science and Technology 86 (2013) 167–176
6.5
And the temperature field is the foundation of interaction process research and roadbed stability analysis. Therefore, parametric analysis of the earth-atmosphere coupled system has been carried out in the following sections to find the system heat transfer characteristic.
6.0 Upwind natural surface Upwind slope Upwind road shoulder Pavement
5.5
Tave-T0
5.0
Leeward road shoulder Leeward slope Leeward natural surface
4.2. Effects of parameter on the coupled heat transfer
4.5 4.0 3.5 3.0 2.5 2.0 -6.0
-5.5
-5.0
-4.5
-4.0
-3.5
-3.0
Annual average enviroment temperature / ºC Fig. 7. Annual average temperature increase of various surfaces under different ambient temperature within the time range of 0 to 10 years.
Fig. 6 is the yearly surface temperature variation of different ground positions since the highway engineering has been established for 20 years. It shows that that the surface temperature of different ground positions all varies in sinusoid trend generally and is higher than the ambient temperature. Because there is no evaporation heat elimination from the ground surface, the temperature of the pavement keeps the highest value among all ground positions. And owing to snow covering, the temperature of natural ground surface has an abrupt step variation from the beginning of April to the end of October. Furthermore, based on the numerical results, the temperature of different ground surface positions, T(ts), have been fitted as the following time related sinusoidal function: T ðt s Þ ¼ a þ b sin
2π t þ cπ 2 30 12 s
ð16Þ
where ts is the computation time step (12 h) and a, b, and c are the fitting parameters of the above equations, the specific values of which are listed in the Table 6. Meanwhile, the numerical results show that the temperature distribution of the permafrost roadbed is influenced by many system factors.
(a) Average convective heat transfer coefficient
4.2.1. Annual average air temperature In the actual engineering constructions, the annual average temperature of ambient air changes its value by region variation. Fig. 7 shows the annual average temperature increase (Tave) of different ground surface positions under the annual average air temperature (T0) of –3.5 °C, –4.0 °C, –4.5 °C, –5.0 °C and –5.5 °C, within the time range of 0 to 10 years. It can be seen that the temperature of different ground surfaces positions increases as the ambient air temperature increase. However, the annual average temperature increase, Tave T0, changes little as the annual average air temperature increases and the difference between ground surfaces' annual average temperature increase and their average value, T ave , is less than 1%. 4.2.2. Annual average wind speed Fig. 8 shows the variation of local heat transfer coefficient and annual average temperature increase of different ground surface positions with different wind speed (vref keeps the constant value) and it can be seen that as the increase of wind speed, local heat transfer coefficient increases linearly as well but annual average temperature increase of different ground surface positions decreases. Among all the ground surface positions, the local heat transfer coefficient of upwind road shoulder exhibits the most serious change for the severe scouring of air flow. And it also shows that the annual average temperature increase of different ground surface positions changes non-linearly as the wind speed increases. This is because the decrease of ground surface temperature results in the decrease of radiation heat exchange between the ground surface and environment, thus the net energy absorption of ground surface varies non-linearly with the wind speed. 4.2.3. Roadbed height Fig. 9 shows the variation of local heat transfer coefficient and annual average temperature increase of different ground surface positions with the roadbed height (H) of 2 m, 3 m, 4 m and 5 m. It shows that due to the wind speed of different ground surface positions increasing as the roadbed height rising, the local heat transfer
(b) Annual average increase of temperature
35 Leeward road shoulder Leeward slope Leeward natural surface
Pavement Leeward road shoulder Leeward slope Leeward natural surface
Upwind natural surface Upwind slope Upwind road shoulder
7 6
25
Tave-T0
hlocal / W ·
m-2
30
8 Upwind natural surface Upwind slope Upwind road shoulder Pavement
20
5 4
15
10 3.0
3
3.5
4.0
4.5
vref / m · s-1
5.0
5.5
2
3.0
3.5
4.0
4.5
5.0
vref / m · s-1
Fig. 8. Annual mean surface heat convection thermal coefficient and temperature increase of various surfaces under different wind speed.
5.5
Z.-Y. Liu et al. / Cold Regions Science and Technology 86 (2013) 167–176
(a) Average convective heat transfer coefficient
(b) Annual average increase of temperature
32
7.0 Leeward road shoulder Leeward slope Leeward natural surface
Upwind natural surface Upwind slope Upwind road shoulder Pavement
30 28
6.5 6.0
26
Upwind natural surface Upwind slope Upwind road shoulder Pavement
5.5
24
Tave-T0
hlocal / W · m-2
175
22 20
5.0 4.5
Leeward road shoulder Leeward slope Leeward natural surface
4.0
18
3.5
16 14
3.0
12
2.5
10 2.0
2.5
3.0
3.5
4.0
4.5
2.0 2.0
5.0
2.5
3.0
H/m
3.5
4.0
4.5
5.0
H/m
Fig. 9. Annual mean surface heat convection thermal coefficient and temperature increase of various surfaces under different roadbed height.
coefficient increases and the annual average temperature increase decreases linearly as the roadbed height increases, especially for the pavement and road shoulder. The upwind road shoulder keeps the highest local heat transfer coefficient and its maximum value is up to 28.28 W/m 2 when the roadbed height reaches 5 m. Meanwhile, the annual average temperature increase of exhibits the most sharply decreasing. Under the same boundary conditions, the pavement's annual average temperature increase reduces 1.29 °C as the roadbed height rises from 2 m to 5 m. However, because the rising of roadbed height only affects the wind speed of natural ground surface that is near the slope, the annual average temperature increase of natural ground surface only shows a little change. 4.2.4. Route direction The route direction of the highway has a great influence on the coupled heat transfer process of the slope, which further affects the temperature field of permafrost roadbed. Accord to the spot investigation results, there are mainly three types of highway route directions in Wudaoliang region: EW, SW45° and SN directions. The annual average temperatures of different ground surface positions under the above three route directions and within the time range of 0 to 30 years are listed in Table 7. It can be seen that as the environmental temperature increases, the annual average temperatures of ground surfaces increase for all three route directions and their annual temperature increases are consistent with environmental annual temperature increases. Besides, it also shows that the route directions have direct influence on the annual average temperature of slope. The leeward slope (sunny side) of EW direction has the maximal annual average temperature among all three directions, and that of SN direction keeps the minimum value. For the upwind slope (shady side), the above rules are on the contrary.
Furthermore, there are obvious annual average temperature difference between upwind slope and leeward slope and such difference also reaches to the maximal under EW route direction. 5. Conclusions In the present work, the heat transfer processes of the comprehensive earth-atmosphere coupled system are simulated by the numerical method. In the computation, the surrounding air environment of the ground surface is brought into the research model and the corresponding unsteady and non-uniform influence factors of earth-atmosphere coupled system, such as the solar radiation, air temperature, wind speed and direction, surface vegetation and evaporation, are also considered in the numerical computation. Based on the above model, the spatial and temporal temperature distribution characteristics of permafrost roadbed surface under the typical working condition of Wudaoliang region are calculated and then the influences of system parameters, including the ambient temperature, wind speed, roadbed height and route trends, on the coupled heat transfer processes are investigated. The computation results show that: 1. The temperature of surface regions exhibit periodic sinusoid trend and always higher than ambient temperature. For different surface regions, the temperature of asphalt pavement keeps the maximum value and the natural surface has an abrupt step change in temperature at beginning of April and end of October; 2. The variation of ambient air temperature has little influence on annual average temperature of permafrost roadbed surface. The deviation between surfaces' annual average temperature increase and their mean value is less than 1%;
Table 7 Annual average surface temperature of different road directions. Time (year)
1 5 10 15 20 Annual mean amplification Annual mean amplification after 5 years
Natural surface
−2.79 −2.73 −2.61 −2.51 −2.40 0.021 0.022
Pavement
0.94 1.01 1.11 1.22 1.33 0.021 0.021
Upwind road shoulder
Upwind slope
Leeward road shoulder
Leeward slope
EW
SW45°
NS
EW
SW45°
NS
EW
SW45°
NS
EW
SW45°
NS
−1.70 −1.63 −1.52 −1.41 −1.31 0.021 0.021
−1.67 −1.60 −1.49 −1.38 −1.27 0.021 0.022
−1.57 −1.51 −1.40 −1.29 −1.18 0.021 0.022
−3.12 −3.09 −2.98 −2.87 −2.76 0.019 0.022
−2.52 −2.49 −2.38 −2.27 −2.16 0.019 0.022
−1.03 −1.00 −0.89 −0.78 −0.67 0.019 0.022
−1.29 −1.21 −1.10 −0.99 −0.88 0.022 0.022
−1.32 −1.24 −1.14 −1.03 −0.92 0.021 0.021
−1.43 −1.36 −1.25 −1.14 −1.03 0.021 0.022
0.34 0.41 0.51 0.62 0.73 0.021 0.021
0.01 0.08 0.18 0.28 0.40 0.021 0.021
−0.92 −0.88 −0.77 −0.66 −0.55 0.019 0.022
176
Z.-Y. Liu et al. / Cold Regions Science and Technology 86 (2013) 167–176
3. Local heat transfer coefficient increases as the wind speed and the roadbed height increases. However, the annual average temperature of various surfaces is on the contrary; 4. The annual mean temperature of leeward slope (sunny side) along EW direction is the highest one among all three route directions and such temperature decreased systematically from SW45° direction to SN direction. However, the annual temperature of upwind slope (shady side) is just contrary to the leeward slope. In the mean time, the annual average temperature difference between upwind slope and leeward slope also reaches the maximal value under EW route direction. Acknowledgments The study is supported by the Technical project of Ministry of Transport of the People's Republic of China (No. 2011-319-495-090). References Duffie, J.A., Beckman, W.A., 1996. Solar Engineering of Thermal Processes. John Wiley & Son Inc., New York. Lai, Y.M., Liu, S.Y., Wu, Z.W., Yu, W.B., 2002. Approximate analytical solution for temperature fields in cold regions circular tunnels. Cold Regions Science and Technology 34 (1), 43–49. Lai, Y.M., Wang, Q.S., Niu, F.J., Zhang, K.H., 2004. Three-dimensional nonlinear analysis for temperature characteristic of ventilated embankment in permafrost regions. Cold Regions Science and Technology 38, 165–184. Liu, Y.Z., Wu, Q.B., Zhang, J.M., Sheng, Y., 2002. Deformation of highway roadbed in permafrost of the Tibetan plateau. Journal of Glaciology and Geocryology 24 (1), 10–15. Liu, J., Zhang, W.W., Shao, J.T., 2009. CFD simulation of convective heat transfer coefficient of external surfaces of buildings. Journal of South China University of Technology 37 (8), 94–98. Mi, L., Lai, Y.M., Zhan, K.H., 2002. 3D nonlinear analyses for temperature field of ventilated embankment in cold regions. Journal of Glaciology and Geocryology 24 (6), 765–769.
Niu, F.J., Zhang, L.X., Qihao Yu, Q.H., Qun, X., 2002. Study on slope types and stability of typical slopes in permafrost regions of the Tibetan plateau. Journal of Glaciology and Geocryology 24 (5), 608–613. Niu, F.J., Liu, X.F., Ma, W., Wu, Q.B., Xu, J., 2008. Monitoring study on the boundary thermal conditions of duct-ventilated embankment in permafrost regions. Cold Regions Science and Technology 53, 305–316. Qiao, Y.L., Guo, S., Tang, Y.H., 2004. Characteristics of diffuse radiation on the QinghaiTibetan plateau. Acta Scientiarum Naturlium Universitatis Nankaisis 41 (3), 69–78. Sheng, Y., Liu, Y.Z., Zhang, J.M., Wu, J.M., et al., 2003. Analysis of the thawing of permafrost underlying Qinghai-Tibet highway. Journal of Glaciology and Geocryology 25 (1), 43–48. Su, B., Ning, X.J., 2004. The numerical study on the ventilated embankment in permafrost regions in QinghaiTibet railway. Cold Regions Science and Technology 38 (2–3), 229–238. Tao, W.Q., 2001. Heat Transfer. Xi'an Jiaotong University Press, Xi'an. Wang, T.X., 2000. Study on calculation principle and critical thickness of roadbed in permafrost area. Chang'an University, Xi'an. Wang, H.N., Wu, M.J., Wu, M.H., 2005. Nonlinear analysis of the influence of pavement types on embankment thermal regime in permafrost regions on the Tibetan plateau. Journal of Glaciology and Geocryology 27 (2), 169–175. Wang, D.P., Fu, Z., Fang, J.H., Li, H.Q., 2008. Influence of solar radiation on surface thermal regime of different pavement types and its permafrost underlying embankment on the Qinghai-Tibetan plateau. Journal of Highway and Transportation Research and Development 25 (3), 38–43. Wu, Q.B., Li, X., Li, W.J., 2001. The response model of permafrost along the QinhaiTibetan highway under climate change. Journal of Glaciology and Geocryology 23 (1), 1–6. Wu, Q.B., Lu, Z.J., Liu, Y.Z., 2005. Permafrost monitoring and its recent changes in Qinghai-Tibet Plateau. Advances in Climate Chang Research 1 (1), 26–28. Xu, X.Z., Wang, J.C., Zhang, L.X., 2001. Physics of Frozen Soils. Science Press, Beijing. Zhang, S.Q., Ding, Y.J., Lu, J., Liu, S.Y., 2004. Simulative study of water-heat process in the Tibetan plateau (I): soil moisture. Journal of Glaciology and Geocryology 26 (4), 384–388. Zhang, M.Y., Lai, Y.M., Liu, Z.Q., Gao, Z.H., 2005. Nonlinear analysis for the cooling effect of Qinghai-Tibetan railway embankment with different structure in permafrost regions. Cold Regions Science and Technology 42, 237–249. Zhong, Q., 1989. Solar radiation budget over the Qinghai-Xizhang Plateau region. Plateau Meteorology 8 (1), 1–12.
Contents lists available at SciVerse ScienceDirect
Cold Regions Science and Technology journal homepage: www.elsevier.com/locate/coldregions
Roadbed temperature study based on earth-atmosphere coupled system in permafrost regions of the Qinghai-Tibet Plateau Zhi-Yun Liu a,⁎, Jian-Bing Chen b, Long Jin b, Yu-Jie Zhang a, Chao Lei a a b
College of Geology Engineering & Geomatics, Chang'an University, Xi'an, Shaanxi 710054, China CCCC First Highway Consultants Co. Ltd, Xi'an, Shaanxi 710054, China
a r t i c l e
i n f o
Article history: Received 5 July 2012 Accepted 10 October 2012 Keywords: Permafrost roadbed Temperature field Earth-atmosphere coupled system Coupled heat transfer Parameter analysis
a b s t r a c t In this work, a two-dimensional unsteady model for analyzing the roadbed temperature in permafrost regions of the Qinghai-Tibet Plateau (QTP) is developed and the heat transfer processes of comprehensive earth-atmosphere coupled system are simulated by the numerical method. In the computation, the surrounding air environment of ground surface is brought into the research model and the corresponding unsteady and non-uniform influence factors of earth-atmosphere coupled system are also considered to obtain reasonable simulation results. Meanwhile, the model is validated with monitored data at a depth of 0.5 m in different roadbed positions. Based on the above model, the spatial and temporal temperature distributions of the roadbed surface in Wudaoliang region are determined. Then the influences of the system parameters, such as the ambient temperature, wind speed, roadbed height and route trends, on the local temperature distribution are investigated. The results show that, temperature of surface regions exhibits a periodic sinusoid trend and always higher than the ambient air temperature. For different surface regions, the temperature of asphalt pavement keeps the maximum value and the natural ground surface has an abrupt step change in temperature at beginning of April and end of October. The variation of ambient air temperature has little influence on annual average temperature of permafrost roadbed surface. The deviations between the surfaces' annual average temperature increase and their average values are less than 1%. Annual average temperature of various surfaces decreases as the wind speed and the roadbed height increase. Furthermore, the annual average temperature of leeward slope (also sunny side) along EW direction is the highest one among all three route directions and such temperature decreased systematically from SW45° direction to SN direction. However, the annual temperature of upwind slope (also shady side) is just contrary to the leeward slope. © 2012 Elsevier B.V. All rights reserved.
Contents 1. 2.
3. 4.
Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Modeling of the permafrost roadbed . . . . . . . . . . . . . . . . . . . . . . 2.1. The physical and mathematical model of earth-atmosphere coupled system 2.2. Boundary condition and other computation settings . . . . . . . . . . . 2.2.1. Velocity and direction of the wind . . . . . . . . . . . . . . . . 2.2.2. Environmental temperature . . . . . . . . . . . . . . . . . . . 2.2.3. Source term of the couple wall . . . . . . . . . . . . . . . . . 2.2.4. Special source term settings of the slope couple wall . . . . . . . 2.3. Numerical procedures . . . . . . . . . . . . . . . . . . . . . . . . . . Validation of the numerical procedure . . . . . . . . . . . . . . . . . . . . . Results and discussions . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4.1. Numerical results of the permafrost roadbed . . . . . . . . . . . . . . . 4.2. Effects of parameter on the coupled heat transfer . . . . . . . . . . . . . 4.2.1. Annual average air temperature . . . . . . . . . . . . . . . . . 4.2.2. Annual average wind speed . . . . . . . . . . . . . . . . . . . 4.2.3. Roadbed height . . . . . . . . . . . . . . . . . . . . . . . . 4.2.4. Route direction . . . . . . . . . . . . . . . . . . . . . . . .
⁎ Corresponding author. Tel.: +86 29 82339356. E-mail address: [email protected] (Z.-Y. Liu). 0165-232X/$ – see front matter © 2012 Elsevier B.V. All rights reserved. http://dx.doi.org/10.1016/j.coldregions.2012.10.005
. . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . .
168 168 168 169 170 170 170 171 171 171 173 173 174 174 174 174 175
168
Z.-Y. Liu et al. / Cold Regions Science and Technology 86 (2013) 167–176
5. Conclusions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Acknowledgments . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . References . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
1. Introduction Temperature distributions of the Qinghai-Tibet Highway (QTH) roadbed in permafrost regions are influenced by many factors, such as climatic conditions, physical parameters and construction scheme. So the temperature field of a highway roadbed would be different from that of the natural ground and such temperature difference shall further influence the stability of underlying permafrost (Niu et al., 2002; Wu et al., 2005). Meanwhile, the temperature field of the highway roadbed is also the prior condition to study the interaction process between the QTH and the permafrost. It not only directly influences the heat and mass transfer process within the roadbed, but also influences morphology and the freeze–thaw process of the underlying permafrost. Therefore, it is necessary to carry out numerical studies about the temperature field of permafrost roadbed (Sheng et al., 2003; Wang et al., 2005; Wu et al., 2001). The heat transfer process of permafrost roadbed includes two aspects: i) unsteady heat transfer process between roadbed and external environment; and ii) unsteady heat transfer process between roadbed and the underlying permafrost. Unlike the main conductive heat transfer process between roadbed and the underlying permafrost, the heat transfer process between roadbed and ambient environment consists of different and coupled heat transfer types, which includes: the solar radiation absorption of highway surface, the long-wave radiation emitting from the highway surface to surrounding environment and the forced-convection between the air and the highway surface. Furthermore, the different pavement structures make the change of absorption coefficient and evaporation capacity, and the different route trends influence the absorption coefficient of the slope surface. And also, some influencing factors (such as air temperature, wind velocity and direction, solar radiation density, absorption coefficient, evaporation capacity, etc.) change with time, which all aggravate the complexity of heat transfer process numerical computation work on the solving of temperature field of permafrost roadbed (Liu et al., 2002; Wang et al., 2008; Zhang et al., 2004). Due to the above reasons, it is very hard to obtain the temperature field of permafrost roadbed directly by the numerical method and this problem has became a challenging research subject in permafrost engineering research (Lai et al., 2002; Zhang et al., 2004). In recent years, many significant numerical and experimental research works had been conducted on the temperature field of permafrost roadbed. Niu et al. (2008) investigated the surface temperature of pipe-ventilated roadbed by placing temperature measuring points and summarized its annual variation rules, which provided the boundary conditions for solving temperature field of permafrost roadbed based on boundary layer theory. Mi et al. (2002) and Su and Ning (2004) simulated the temperature field of pipe-ventilated roadbed using the constant temperature boundary condition and the results showed that the duct-ventilated structure can effectively cool the roadbed. Zhang et al. (2005) and Lai et al. (2004) compared the temperature distributions of different roadbed structures and the results showed that the ripped-rock embankment would be the better choice for the high temperature permafrost region. However, the above research works on the temperature field of permafrost roadbed are on the premise of knowing the detailed information of the surface temperature distribution, which will need a long time field-observation or a suitable empirical formula. In this way, the influences of a large number of coefficients such as air temperature, wind direction and air speed, etc., are neglected. Consequently, in the numerical results there exists a large deviation from the actual situation and with rigid regional limitations. So in this work, to reveal the spatial and temporal temperature distribution characteristics of
175 176 176
permafrost roadbed surface and to optimize the design of highway in permafrost region, the complex heat transfer processes of the comprehensive earth-atmosphere coupled system have been simulated. In the simulation, the surrounding air environment of the ground surface is brought into the research model and the environmental influencing factors, such as the solar radiation, air temperature, wind speed and direction, surface vegetation and evaporation, are also considered. Based on the above model, the spatial and temporal temperature distributions of the roadbed surface in Wudaoliang region are determined. Then the influences of the system parameters, such as the ambient temperature, wind speed, roadbed height and route trends, on the local temperature distribution are investigated. 2. Modeling of the permafrost roadbed 2.1. The physical and mathematical model of earth-atmosphere coupled system Owing to direct heat exchanging processes between roadbed and underlying natural permafrost, the temperature distribution conditions of underlying natural permafrost have great influence on the temperature field of roadbed since the engineering construction is finished. To obtain a more reasonable calculation result, the temperature distribution of underlying natural permafrost that is on the time node of construction work finished should be imported into the computation module as the initial field. So in this work, there are two corresponding physical models for numerical calculation: i) the initial field model with natural permafrost and air environment; and ii) the permafrost roadbed model with comprehensive earth-atmosphere coupled system. The schematic of the above two models is shown in Fig. 1. As shown in the figure, the underlying natural permafrost is composed of silty clay, gravelly clay and fully weathered mudstone and the thickness of each layer is 2 m, 5 m and 23 m from the top down. The highway engineering is built by subgrade soil and asphalt pavement. The width and thickness of pavement structure is designed as 7 m and 0.3 m respectively. In order to decrease the influence of top environment, the height of the air layer is chosen as 20 m. The longitudinal length of the initial field mode is defined as 54 m. To avoid the influence of entrance effect, the length of bilateral natural permafrost of the roadbed is defined as 60 m. The width of bilateral road shoulder is defined as 1.5 m. Meanwhile, the height and slope ratio of the roadbed is defined as 3 m and 1:1.5 respectively, and the route direction highway is chosen as NS direction. The Wudaoliang region is located at the 4000 m height of the Qinghai-Tibet Plateau, so the reference pressure is set as 60,000 Pa in the numerical simulation. Furthermore, considering that the variation range of the ambient air temperature is within 20 centigrade, the air properties are defined as a fixed value in the numerical simulation. The detailed physical properties of air used in this work are listed as following: ρ=0.7441 kg/m3, Cp =1005 J/(kg⋅K), μ=1.751×10−5kg/(m⋅s). And because the rough ground surface can only be defined as the smooth surface in the calculation software, so the conductivity of air (λ) is amplified as 0.0586 W/(m·K) to obtain a reasonable results. The density of pavement structure, subgrade soil and different natural permafrost layers is also defined as a fixed value. But their specific heat capacity and conductivity are set as the temperature correlation function, except for the fixed definition of pavement structure. The specific value of the above mentioned physical properties are listed in Table 1 (Xu et al., 2001).
Z.-Y. Liu et al. / Cold Regions Science and Technology 86 (2013) 167–176
169
(b) Permafrost roadbed model
(a) Initial field model
Fig. 1. Schematics of the two computation model.
A two-dimensional unsteady model has been developed in the present work to predict the spatial and temporal temperature distribution characteristics of Qinghai-Tibet plateau permafrost roadbed surface. The heat transfer processes between the ground surface and ambient environmental are coupled in the boundary condition settings and the source terms of energy equations. The vector form of the continuity, momentum and energy equations are expressed as: Continuity equation: → div U ¼ 0
ð1Þ
Momentum equation: → ∂ðρui Þ ∂p þ divðη gradui Þ þ div ρui U ¼ − ∂t ∂xi
ð2Þ
2.2. Boundary condition and other computation settings In this paper, the surrounding air environment is brought into the research model and computation model is modified from the closed system (based on the boundary layer theory) to the opening system with air inlet and outlet. Therefore, more influencing factors are considered in the numerical simulation and the corresponding boundary conditions become more complex. The detail boundary condition settings are listed as the following: i) both the left and right sides of the air layer are defined as the velocity inlet boundary conditions. The wind enters the computation area vertically to the boundary and its speed, direction and temperature change annually during a year. Except for changing with time, the wind speed is also set as the function of the height from the ground. All the boundary condition parameter settings of the left and right sides are same, but only the wind direction is set as the contrast value to keep the quantity conservation of system. And the turbulence density and turbulence length scale can be calculated by: 1 I ¼ 0:16 ReDH 8
Energy equation: ! → ∂ðρT Þ k gradT þ ST þ div ρ U T ¼ div cp ∂t
ð3Þ
ð4Þ
l ¼ 0:07DH
ii) top of the air layer is defined as the wall boundary condition and its temperature varies annually during a year; iii) the bilateral side of the
Table 1 Physical property of the soil layers. Physical properties ρ (kg ⋅ m−) Cp (J ⋅ kg−1 ⋅ K−1)
λ (W ⋅ m−1 ⋅ K1)
T (K)
Pavement
Roadbed
Silty clay
Gravelly clay
Mudstone
233.15 253.15 263.15 268.15 271.15 272.15 272.65 272.95 273.05 273.15 313.15 233.15 272.75 272.91 273.07 273.23 313.15
2177.78 815.78 815.78 815.78 815.78 815.78 815.78 815.78 815.78 815.78 815.78 815.78 1.55 1.55 1.55 1.55 1.55 1.55
2183.6 8.38E + 02 8.43E + 02 8.52E + 02 8.70E + 02 1.14E + 03 2.06E + 03 4.12E + 03 8.24E + 03 1.65E + 04 1.02E + 03 1.02E + 03 1.5 1.5 1.49 1.48 1.4 1.4
2080 1.05E + 03 1.06E + 03 1.20E + 03 1.44E + 03 2.88E + 03 6.01E + 03 1.20E + 04 2.40E + 04 4.94E + 04 1.17E + 03 1.17E + 03 1.8 1.8 1.78 1.73 1.5 1.5
2340 8.78E + 02 8.97E + 02 9.83E + 02 1.07E + 03 1.92E + 03 3.85E + 03 7.69E + 03 1.50E + 04 2.95E + 04 1.22E + 03 1.22E + 03 2.2 2.2 2.13 1.94 1.0 1.0
2310 2.64E + 06 2.70E + 06 2.80E + 06 3.00E + 06 4.00E + 06 8.00E + 06 1.70E + 07 3.30E + 07 6.56E + 07 2.97E + 06 2.97E + 06 2.50 2.50 2.47 2.39 2.0 2.0
170
Z.-Y. Liu et al. / Cold Regions Science and Technology 86 (2013) 167–176
Table 2 Average wind velocity of different months in the Wudaoliang region. Month
Jan.
Feb.
Mar.
Apr.
May
June
July
Aug.
Sep.
Oct.
Nov.
Dec.
v (m ⋅ s−1)
5.69
5.83
5.92
4.79
4.37
4.11
3.72
3.37
3.38
3.58
4.29
5.23
permafrost layers are all defined as the adiabatic boundary condition; iv) the bottom fully weathered mudstone layer is defined as the heat flux boundary condition with the specific value of 0.06 W/m2 (Zhang et al., 2005); v) the natural permafrost and road engineering surface are defined as the 0.1 m thick couple heat transfer wall. The above mentioned complex heat transfer processes between the ground surface and ambient environment are coupled in the source term computation of couple wall. In the following sections, the specific value of air temperature, wind speed and direction of different times as well as the computation method of couple wall's source terms will be given. 2.2.1. Velocity and direction of the wind According to the research, the wind velocity of atmospheric environment is the function of the height from the ground (Liu et al., 2009): v ¼ vref ðh=href Þ
ð5Þ
where, vref is the wind speed of reference height; href is the height reference spot and chosen as 3 m in this work. Meanwhile, the wind speed and direction also change periodically for years. The average wind speed of different months in the Wudaoliang region is listed in Table 2. To simplify the numerical computation, the wind direction is selected as the maximum wind frequency during a certain time. And based on the meteorological data, the wind direction is defined as from left to right in the time range of Feb. 20 to Aug. 20 and from right to left in the other time of year. 2.2.2. Environmental temperature The ambient air temperature changes periodically during a year and its calculation function is set as the sinusoid in this paper: 2πt 2π þ 3600 24 30 12 3
6
4
3
2
þ 0:4266d −0:1875d−1:455d þ 135:8
ð8Þ
where d is the period order of ten days during a year (since Aug. 20) and its value range from 1 to 36 (assuming 30 days a month). Because the local ground surface temperature is different for the surface position and changes with times, so the value of Qrad also changes with time and local surface temperature: 4 4 Q rad ¼ εσ T local −T sky
ð9Þ
where ε is the emissivity of ground surface and defined as 0.9; Tlocal is the temperature of local ground surface region and can be derived through the numerical calculation for every time step; Tsky is the sky temperature and can be calculated as (Duffie and Beckman, 1996): 1:5
T sky ¼ 0:0552T air
ð10Þ
ð6Þ Qeva can be calculated as:
In this equation, t is the influence time since the construction work of highway is finished (Aug. 20); T0 is the mean annual temperature of Wudaoliang region and its specific value is defined as −5.3 °C based on the local meteorological data; g(t) is the annual temperature ascending speed and defined as 0.022 °C/year (Wu et al., 2005). Because it is assumed that all the parameters keep a constant value or vary periodically in the initial field model, so the term of g(t) is not included in the air temperature computation equation in such model. 2.2.3. Source term of the couple wall The source term of couple wall is the net energy absorbed by ground surface and it can be obtained by analyzing the heat transfer process between roadbed and ambient environment: Q sou ¼ αQ sol −Q rad −Q eva
7
Q sol ¼ 3:214e−7d −4:111e−5d þ 0:002025d −0:04629d
0:14
T air ðt Þ ¼ T 0 þ g ðt Þ þ 11:3 sin
road shoulder are both 0.72 and pavement is 0.9. However, the absorption coefficient of the natural ground surface varies with time for the ground vegetations. Considering that Wudaoliang belongs to the typical desert steppe region, the absorption coefficient of the natural ground surface is simplified as changing with seasons and its specific value is listed in Table 3. The solar radiation intensity casting on the ground surface includes the direct normal solar radiation and diffuse solar radiation and its value changes periodically during a year in generally. For lack of comprehensive meteorological data, the solar radiation intensity data used in the numerical computation are based on the whole Qinghai-Tibet Plateau solar radiation data (Zhong, 1989) and obtained by the following time correlation function:
ð7Þ
where α is the absorption coefficient of the ground surface; Qsol is the total solar radiation projected on the ground surface; Qrad is the long-wave radiation heat loss from ground surface to environment; Qeva is the heat loss of surface evaporation. Because the specific values of the above terms in Eq. (7) fluctuate for different surface types and some of them change with time, so the source terms of the couple wall should be determined according to the surface type. The detailed computation methods of every term in Eq. (7) are accounted below. The absorption coefficients of the slope, road shoulder and pavement surface are defined as the fixed values, where the slope and
Q eva ¼ UG
ð11Þ
where G is the latent heat of vaporization of water and defined a fixed value of 2500 kJ/kg; U is the month surface evaporation quantity and can be calculated as (Wang, 2000): U ¼ 0:3082⋅uw −0:849
ð12Þ
and uw is the month water surface evaporation of the Wudaoliang region, which can be obtained by inquiring the hydrology data (see Table 4). Meanwhile, considering the climatic condition of the Wudaoliang region, Qeva has a certain value only from April to November in this paper for the ground surface is exposed to the ambient environment and reduces to Table 3 Absorbing coefficient of the natural ground surface under different seasons. Time range
Spring
Summer
Autumn
Winter
Alpine steppe
0.644
0.719
0.691
0.668
Table 4 Water surface evaporation of different months. Month
Apr.
May
June
July
Aug.
Sep.
Oct.
Nov.
uw (mm)
114
156
158
167
161
112
79
68
Z.-Y. Liu et al. / Cold Regions Science and Technology 86 (2013) 167–176 Table 5 Slope coefficient of different months and different route directions. Month
Jan. Feb. Mar. Apr. May June July Aug. Sep. Oct. Nov. Dec.
NS
SW45°
EW
Upwind
Leeward
Upwind
Leeward
Upwind
Leeward
0.958 0.946 0.934 0.931 0.928 0.922 0.928 0.931 0.934 0.946 0.958 0.97
0.958 0.946 0.934 0.931 0.928 0.922 0.928 0.931 0.934 0.946 0.958 0.97
0.58 0.658 0.76 0.82 0.886 0.904 0.886 0.82 0.76 0.658 0.58 0.554
1.3 1.18 1.066 1.0 0.94 0.916 0.94 1.0 1.066 1.18 1.3 1.366
0.418 0.52 0.67 0.79 0.88 0.898 0.88 0.79 0.67 0.52 0.418 0.4
1.414 1.24 1.126 1.03 0.992 0.898 0.992 1.03 1.126 1.24 1.414 1.498
2.2.4. Special source term settings of the slope couple wall Compared to the surfaces of the natural ground, road shoulder and pavement, the slope surface has different solar radiation absorption quantities for the different route directions and its angle factor does not keep to 1.0 for the slope inclination. Therefore, to reflect the effects of the route direction and angle factor, the calculation equation of the slope couple wall's source term has been modified as: ð13Þ
where k and X is the slope coefficient and angle factor respectively. Because the direct solar radiation accounts for 60% and diffuse solar radiation is 40% in the total solar radiation of Qinghai-Tibet Plateau (Qiao et al., 2004), so the integral slope coefficient used in present work is defined as: κ ¼ 0:6κ dir þ 0:4κ dif
ð14Þ
where κdif is the slope coefficient of diffuse solar radiation and defined as 1.0 for its non-directional characteristic; κdir is the slope coefficient of direct solar radiation and also different from the sunny slope and shady slope as well (Wang, 2000):
κ dir−shady ¼ ð cot α s − cot βj cosðδ−θÞjÞ sin α s κ dir−sun ¼ ð cot α s − cot βj cosðδ−θÞjÞ sin α s
above equations, it can be inferred that the slope coefficient not only changes with time but also varies with the route direction. To simplify the numerical computation, the slope coefficient is assumed as the month of piecewise function and the detailed slope coefficients of different months (from June to Dec.) and different route directions (NS, SW45° and EW) used in this work are listed in the Table 5. The angle factor, X, can be obtained through the radiation heat transfer related theory formulas (Tao, 2001). Under the physical model used in this work, the specific value of the angle coefficient is calculated as 0.852. 2.3. Numerical procedures
0 in other times for the snow covering. And the detailed time of snow melting and freezing is judged by whether the ground surface temperature is larger than 0 °C or not. Furthermore, owing to the structure of the asphalt pavement, the heat loss carried away by surface vaporization is not included in its source term computation.
Q sou ¼ αkQ sol −XQ rad −Q eva
ð15Þ
where αs is subgrade slope angle; β is sun elevation angle; θ is sun azimuth angle; and δ is the angle between east and north. Under the
The geometry and quadrilateral mesh of the permafrost roadbed are created by GAMBIT 2.3.16 software, as shown in Fig. 2. To get better numerical results, local grid refinement is applied on the roadbed slope, road shoulder and pavement and the final ratio between minimal grid scale and maximal grid scale is 1/15. The governing equations for the forced convection and heat conduction are solved in the FLUENT 6.3 software package. The 2-D, unsteady state and turbulence model with implicit solver is adopted in the simulation. In the boundary condition settings, except for the insulation boundary condition of the bilateral side of permafrost layers, other complex boundary conditions and source term of the couple wall are imported into the FLUENT by UDF program. The properties of working fluid are chosen as the air at 273.15 K and the physical properties of different soils are taken as the step function of temperature. Meanwhile, in order to embody the influence of periodically changing environment on the coupled heat exchange of whole computation model, ambient temperature, wind speed, solar radiation and evaporation capacity are set as the annual periodical change function. The SIMPLE algorithm has been used for pressure velocity coupling. Standard scheme is adopted for the discretization of the pressure term and the second order upwind scheme is applied to discretize other terms. A convergence criterion of 10 −6 is imposed on the residuals of the momentum, energy and 10−5 on continuity the equations. The time step is defined as 12 h and the inner iterations of each time step are chosen as 200. Besides, the calculation time of initial field model is defined as 100 years and then its calculation results are imported into the permafrost roadbed model as the initial field of non-steady state computation. 3. Validation of the numerical procedure In order to validate the present numerical procedure, the comparative studies have been carried out for the temperature monitoring data and numerical calculation result of No.7 Tangbei permafrost roadbed in Wudaoliang. In the numerical computation, the prototype of the computation model is based on the geological section map of the Wudaoliang region in “Frozen Soil Engineering Geology Section for
20
Y /m
10 0 -10 -20 -30 -30
-20
-10
0
10
20
30
40
171
50
60
70
80
90
X /m Fig. 2. Computation grid of the permafrost roadbed.
100
172
Z.-Y. Liu et al. / Cold Regions Science and Technology 86 (2013) 167–176
4.0 0.0
Y /m
-4.0 -8.0
-12.0 -16.0
30
20
10
0.0
10
20
30
X /m Fig. 3. Geological section map of highway in Wudaoliang region.
(a) Natural ground surface
(b) Pavement
6
15 Numerical results
Experimental results
2 0 -2 -4 -6 -8
-10 -12
9 6 3 0 -3 -6 -9
-12 2008/2/20
2008/8/21
2009/2/20
2009/8/20
2010/3/18
-15
2010/9/21
Date / yy-mm-dd
2008/8/21
2009/2/20
2009/8/20
2010/3/5
2010/9/5
--
(d) Leeward road shoulder Temperature of upwind road shoulder / ºC
12 Experimental results
2008/2/18
Date / yy-mm-dd
(c) Upwind road shoulder Temperature of leeward road shoulder / ºC
Numerical results
12
4
Temperature of pavement / ºC
Temperature of natural surface / ºC
Experimental results
Numerical results
9 6 3 0 -3 -6 -9
15 12
Experimental results
Numerical results
9 6 3 0 -3 -6 -9
-12
-12
-15
2008/2/18 2008/8/21 2009/2/20 2009/8/20 2010/3/5 2010/9/5
Date / yy-mm-dd
--
2008/2/18 2008/8/21 2009/2/20 2009/8/20 2010/3/5 2010/9/5
Date / yy-mm-dd
Fig. 4. The temperature distributions of 0.5 m depth soil under different positions and times.
--
Z.-Y. Liu et al. / Cold Regions Science and Technology 86 (2013) 167–176
173
10
v /m · s-1
hlocal /W · m-2
Tlocal / ºC
Local ground surface temperature
Toe of slope
8 6 4 2 32 28 24 20 16 12 8 8
Local heat transfer coefficient
Right road shoulder
Left road shoulder Wind velocity of 0.3m height
6 4 2
Toe of slope
0
0
20
40
60
80
100
120
140
Distance of X direction / m Fig. 5. Wind speed of 0.3 m height, local heat transfer coefficient and temperature of different ground surface positions along the horizontal direction on March 20 of 20th year.
Qinghai to Tibet highway” (as shown in Fig. 3), which is drawn by Cold and Arid Regions Environmental and Engineering Research Institute of CAS and CCCC First Highway Consultants Co., Ltd. in 1999. The detailed computation parameters are the same as the predefined parameters of permafrost roadbed model and according to the 30 years meteorological data of Wudaoliang region (from 1971 to 2000), the mean annual air temperature is defined as −5.3 °C. The temperature distribution of soil with 0.5 m depth of different positions (natural permafrost, pavement, upwind and leeward road shoulder) and different times (from 21th year to 23th year since the highway engineering finished) are shown in Fig. 4, and it can be seen that the numerical results are in good agreement with the drilling measured data.
15 10
T / ºC
5 0 Upwind natural surface Upwind slope Pavement Leeward slope Leeward natural surface Enviroment
-5 -10 -15 0
30
60
90 120 150 180 210 240 270 300 330 360
Times / Day Fig. 6. The yearly surface temperature variation of different ground positions since the highway engineering has been established for 20 years.
4. Results and discussions 4.1. Numerical results of the permafrost roadbed Based on previous mathematical model and calculation method, the numerical simulation is carried out for complex coupled heat transfer within permafrost roadbed of Wudaoliang region. The wind speed at the height of 0.3 m (v), local heat transfer coefficient (hlocal) and temperature (Tlocal) of different positions along the horizontal direction on March 20 of 20th year (since the highway engineering finished) are shown in Fig. 5. It can be seen that owing to the shading effect of roadbed, the wind speed and local heat transfer coefficient decrease but the local surface temperature increases near the upwind slope; then the wind speed continue to ascend in the upwind slope and reaches its maximum value in the top of slope, and corresponding local heat transfer coefficient and surface temperature increase for the wind enhanced heat transfer; when the air flow reaches the leeward slope region, the wind speed and local heat transfer coefficient decrease and the local surface temperature increases firstly, then as the dropping circumfluence phenomenon appears in the air flow and strengthens the ground heat transfer, local heat transfer coefficient increases and surface temperature decreases rapidly; the wind speed of leeward natural ground is smaller than that of the upwind, however, due to the air flow backwash, the local heat transfer coefficient of leeward is larger than that of the upwind and thus the surface temperature of leeward is less than that of upwind.
Table 6 Parameters of the different ground surfaces' temperature fitting equations. Fitting parameters
Upwind natural surface
Upwind slope
Pavement
Leeward slope
Leeward natural surface
a b c σ2
2.54 11.65 0.68 0.98
4.22 13.48 0.64 1
6.21 15.12 0.66 1
4.34 14.01 0.66 0.99
2.47 11.43 0.66 0.99
174
Z.-Y. Liu et al. / Cold Regions Science and Technology 86 (2013) 167–176
6.5
And the temperature field is the foundation of interaction process research and roadbed stability analysis. Therefore, parametric analysis of the earth-atmosphere coupled system has been carried out in the following sections to find the system heat transfer characteristic.
6.0 Upwind natural surface Upwind slope Upwind road shoulder Pavement
5.5
Tave-T0
5.0
Leeward road shoulder Leeward slope Leeward natural surface
4.2. Effects of parameter on the coupled heat transfer
4.5 4.0 3.5 3.0 2.5 2.0 -6.0
-5.5
-5.0
-4.5
-4.0
-3.5
-3.0
Annual average enviroment temperature / ºC Fig. 7. Annual average temperature increase of various surfaces under different ambient temperature within the time range of 0 to 10 years.
Fig. 6 is the yearly surface temperature variation of different ground positions since the highway engineering has been established for 20 years. It shows that that the surface temperature of different ground positions all varies in sinusoid trend generally and is higher than the ambient temperature. Because there is no evaporation heat elimination from the ground surface, the temperature of the pavement keeps the highest value among all ground positions. And owing to snow covering, the temperature of natural ground surface has an abrupt step variation from the beginning of April to the end of October. Furthermore, based on the numerical results, the temperature of different ground surface positions, T(ts), have been fitted as the following time related sinusoidal function: T ðt s Þ ¼ a þ b sin
2π t þ cπ 2 30 12 s
ð16Þ
where ts is the computation time step (12 h) and a, b, and c are the fitting parameters of the above equations, the specific values of which are listed in the Table 6. Meanwhile, the numerical results show that the temperature distribution of the permafrost roadbed is influenced by many system factors.
(a) Average convective heat transfer coefficient
4.2.1. Annual average air temperature In the actual engineering constructions, the annual average temperature of ambient air changes its value by region variation. Fig. 7 shows the annual average temperature increase (Tave) of different ground surface positions under the annual average air temperature (T0) of –3.5 °C, –4.0 °C, –4.5 °C, –5.0 °C and –5.5 °C, within the time range of 0 to 10 years. It can be seen that the temperature of different ground surfaces positions increases as the ambient air temperature increase. However, the annual average temperature increase, Tave T0, changes little as the annual average air temperature increases and the difference between ground surfaces' annual average temperature increase and their average value, T ave , is less than 1%. 4.2.2. Annual average wind speed Fig. 8 shows the variation of local heat transfer coefficient and annual average temperature increase of different ground surface positions with different wind speed (vref keeps the constant value) and it can be seen that as the increase of wind speed, local heat transfer coefficient increases linearly as well but annual average temperature increase of different ground surface positions decreases. Among all the ground surface positions, the local heat transfer coefficient of upwind road shoulder exhibits the most serious change for the severe scouring of air flow. And it also shows that the annual average temperature increase of different ground surface positions changes non-linearly as the wind speed increases. This is because the decrease of ground surface temperature results in the decrease of radiation heat exchange between the ground surface and environment, thus the net energy absorption of ground surface varies non-linearly with the wind speed. 4.2.3. Roadbed height Fig. 9 shows the variation of local heat transfer coefficient and annual average temperature increase of different ground surface positions with the roadbed height (H) of 2 m, 3 m, 4 m and 5 m. It shows that due to the wind speed of different ground surface positions increasing as the roadbed height rising, the local heat transfer
(b) Annual average increase of temperature
35 Leeward road shoulder Leeward slope Leeward natural surface
Pavement Leeward road shoulder Leeward slope Leeward natural surface
Upwind natural surface Upwind slope Upwind road shoulder
7 6
25
Tave-T0
hlocal / W ·
m-2
30
8 Upwind natural surface Upwind slope Upwind road shoulder Pavement
20
5 4
15
10 3.0
3
3.5
4.0
4.5
vref / m · s-1
5.0
5.5
2
3.0
3.5
4.0
4.5
5.0
vref / m · s-1
Fig. 8. Annual mean surface heat convection thermal coefficient and temperature increase of various surfaces under different wind speed.
5.5
Z.-Y. Liu et al. / Cold Regions Science and Technology 86 (2013) 167–176
(a) Average convective heat transfer coefficient
(b) Annual average increase of temperature
32
7.0 Leeward road shoulder Leeward slope Leeward natural surface
Upwind natural surface Upwind slope Upwind road shoulder Pavement
30 28
6.5 6.0
26
Upwind natural surface Upwind slope Upwind road shoulder Pavement
5.5
24
Tave-T0
hlocal / W · m-2
175
22 20
5.0 4.5
Leeward road shoulder Leeward slope Leeward natural surface
4.0
18
3.5
16 14
3.0
12
2.5
10 2.0
2.5
3.0
3.5
4.0
4.5
2.0 2.0
5.0
2.5
3.0
H/m
3.5
4.0
4.5
5.0
H/m
Fig. 9. Annual mean surface heat convection thermal coefficient and temperature increase of various surfaces under different roadbed height.
coefficient increases and the annual average temperature increase decreases linearly as the roadbed height increases, especially for the pavement and road shoulder. The upwind road shoulder keeps the highest local heat transfer coefficient and its maximum value is up to 28.28 W/m 2 when the roadbed height reaches 5 m. Meanwhile, the annual average temperature increase of exhibits the most sharply decreasing. Under the same boundary conditions, the pavement's annual average temperature increase reduces 1.29 °C as the roadbed height rises from 2 m to 5 m. However, because the rising of roadbed height only affects the wind speed of natural ground surface that is near the slope, the annual average temperature increase of natural ground surface only shows a little change. 4.2.4. Route direction The route direction of the highway has a great influence on the coupled heat transfer process of the slope, which further affects the temperature field of permafrost roadbed. Accord to the spot investigation results, there are mainly three types of highway route directions in Wudaoliang region: EW, SW45° and SN directions. The annual average temperatures of different ground surface positions under the above three route directions and within the time range of 0 to 30 years are listed in Table 7. It can be seen that as the environmental temperature increases, the annual average temperatures of ground surfaces increase for all three route directions and their annual temperature increases are consistent with environmental annual temperature increases. Besides, it also shows that the route directions have direct influence on the annual average temperature of slope. The leeward slope (sunny side) of EW direction has the maximal annual average temperature among all three directions, and that of SN direction keeps the minimum value. For the upwind slope (shady side), the above rules are on the contrary.
Furthermore, there are obvious annual average temperature difference between upwind slope and leeward slope and such difference also reaches to the maximal under EW route direction. 5. Conclusions In the present work, the heat transfer processes of the comprehensive earth-atmosphere coupled system are simulated by the numerical method. In the computation, the surrounding air environment of the ground surface is brought into the research model and the corresponding unsteady and non-uniform influence factors of earth-atmosphere coupled system, such as the solar radiation, air temperature, wind speed and direction, surface vegetation and evaporation, are also considered in the numerical computation. Based on the above model, the spatial and temporal temperature distribution characteristics of permafrost roadbed surface under the typical working condition of Wudaoliang region are calculated and then the influences of system parameters, including the ambient temperature, wind speed, roadbed height and route trends, on the coupled heat transfer processes are investigated. The computation results show that: 1. The temperature of surface regions exhibit periodic sinusoid trend and always higher than ambient temperature. For different surface regions, the temperature of asphalt pavement keeps the maximum value and the natural surface has an abrupt step change in temperature at beginning of April and end of October; 2. The variation of ambient air temperature has little influence on annual average temperature of permafrost roadbed surface. The deviation between surfaces' annual average temperature increase and their mean value is less than 1%;
Table 7 Annual average surface temperature of different road directions. Time (year)
1 5 10 15 20 Annual mean amplification Annual mean amplification after 5 years
Natural surface
−2.79 −2.73 −2.61 −2.51 −2.40 0.021 0.022
Pavement
0.94 1.01 1.11 1.22 1.33 0.021 0.021
Upwind road shoulder
Upwind slope
Leeward road shoulder
Leeward slope
EW
SW45°
NS
EW
SW45°
NS
EW
SW45°
NS
EW
SW45°
NS
−1.70 −1.63 −1.52 −1.41 −1.31 0.021 0.021
−1.67 −1.60 −1.49 −1.38 −1.27 0.021 0.022
−1.57 −1.51 −1.40 −1.29 −1.18 0.021 0.022
−3.12 −3.09 −2.98 −2.87 −2.76 0.019 0.022
−2.52 −2.49 −2.38 −2.27 −2.16 0.019 0.022
−1.03 −1.00 −0.89 −0.78 −0.67 0.019 0.022
−1.29 −1.21 −1.10 −0.99 −0.88 0.022 0.022
−1.32 −1.24 −1.14 −1.03 −0.92 0.021 0.021
−1.43 −1.36 −1.25 −1.14 −1.03 0.021 0.022
0.34 0.41 0.51 0.62 0.73 0.021 0.021
0.01 0.08 0.18 0.28 0.40 0.021 0.021
−0.92 −0.88 −0.77 −0.66 −0.55 0.019 0.022
176
Z.-Y. Liu et al. / Cold Regions Science and Technology 86 (2013) 167–176
3. Local heat transfer coefficient increases as the wind speed and the roadbed height increases. However, the annual average temperature of various surfaces is on the contrary; 4. The annual mean temperature of leeward slope (sunny side) along EW direction is the highest one among all three route directions and such temperature decreased systematically from SW45° direction to SN direction. However, the annual temperature of upwind slope (shady side) is just contrary to the leeward slope. In the mean time, the annual average temperature difference between upwind slope and leeward slope also reaches the maximal value under EW route direction. Acknowledgments The study is supported by the Technical project of Ministry of Transport of the People's Republic of China (No. 2011-319-495-090). References Duffie, J.A., Beckman, W.A., 1996. Solar Engineering of Thermal Processes. John Wiley & Son Inc., New York. Lai, Y.M., Liu, S.Y., Wu, Z.W., Yu, W.B., 2002. Approximate analytical solution for temperature fields in cold regions circular tunnels. Cold Regions Science and Technology 34 (1), 43–49. Lai, Y.M., Wang, Q.S., Niu, F.J., Zhang, K.H., 2004. Three-dimensional nonlinear analysis for temperature characteristic of ventilated embankment in permafrost regions. Cold Regions Science and Technology 38, 165–184. Liu, Y.Z., Wu, Q.B., Zhang, J.M., Sheng, Y., 2002. Deformation of highway roadbed in permafrost of the Tibetan plateau. Journal of Glaciology and Geocryology 24 (1), 10–15. Liu, J., Zhang, W.W., Shao, J.T., 2009. CFD simulation of convective heat transfer coefficient of external surfaces of buildings. Journal of South China University of Technology 37 (8), 94–98. Mi, L., Lai, Y.M., Zhan, K.H., 2002. 3D nonlinear analyses for temperature field of ventilated embankment in cold regions. Journal of Glaciology and Geocryology 24 (6), 765–769.
Niu, F.J., Zhang, L.X., Qihao Yu, Q.H., Qun, X., 2002. Study on slope types and stability of typical slopes in permafrost regions of the Tibetan plateau. Journal of Glaciology and Geocryology 24 (5), 608–613. Niu, F.J., Liu, X.F., Ma, W., Wu, Q.B., Xu, J., 2008. Monitoring study on the boundary thermal conditions of duct-ventilated embankment in permafrost regions. Cold Regions Science and Technology 53, 305–316. Qiao, Y.L., Guo, S., Tang, Y.H., 2004. Characteristics of diffuse radiation on the QinghaiTibetan plateau. Acta Scientiarum Naturlium Universitatis Nankaisis 41 (3), 69–78. Sheng, Y., Liu, Y.Z., Zhang, J.M., Wu, J.M., et al., 2003. Analysis of the thawing of permafrost underlying Qinghai-Tibet highway. Journal of Glaciology and Geocryology 25 (1), 43–48. Su, B., Ning, X.J., 2004. The numerical study on the ventilated embankment in permafrost regions in QinghaiTibet railway. Cold Regions Science and Technology 38 (2–3), 229–238. Tao, W.Q., 2001. Heat Transfer. Xi'an Jiaotong University Press, Xi'an. Wang, T.X., 2000. Study on calculation principle and critical thickness of roadbed in permafrost area. Chang'an University, Xi'an. Wang, H.N., Wu, M.J., Wu, M.H., 2005. Nonlinear analysis of the influence of pavement types on embankment thermal regime in permafrost regions on the Tibetan plateau. Journal of Glaciology and Geocryology 27 (2), 169–175. Wang, D.P., Fu, Z., Fang, J.H., Li, H.Q., 2008. Influence of solar radiation on surface thermal regime of different pavement types and its permafrost underlying embankment on the Qinghai-Tibetan plateau. Journal of Highway and Transportation Research and Development 25 (3), 38–43. Wu, Q.B., Li, X., Li, W.J., 2001. The response model of permafrost along the QinhaiTibetan highway under climate change. Journal of Glaciology and Geocryology 23 (1), 1–6. Wu, Q.B., Lu, Z.J., Liu, Y.Z., 2005. Permafrost monitoring and its recent changes in Qinghai-Tibet Plateau. Advances in Climate Chang Research 1 (1), 26–28. Xu, X.Z., Wang, J.C., Zhang, L.X., 2001. Physics of Frozen Soils. Science Press, Beijing. Zhang, S.Q., Ding, Y.J., Lu, J., Liu, S.Y., 2004. Simulative study of water-heat process in the Tibetan plateau (I): soil moisture. Journal of Glaciology and Geocryology 26 (4), 384–388. Zhang, M.Y., Lai, Y.M., Liu, Z.Q., Gao, Z.H., 2005. Nonlinear analysis for the cooling effect of Qinghai-Tibetan railway embankment with different structure in permafrost regions. Cold Regions Science and Technology 42, 237–249. Zhong, Q., 1989. Solar radiation budget over the Qinghai-Xizhang Plateau region. Plateau Meteorology 8 (1), 1–12.