我国冻土面积位居世界第三位，常年冻土和季节性冻土面积分别约占陆地总面积的21.5%和53.5%[1]。季节性冻土的冻胀与融沉变形会对建于其上的工程设施造成危害。特别是近年来，我国高速公路、高速铁路不断兴建，这类设施对路基变形的控制要求较高。因此，冻胀变形成为工程界关注的一个重要问题。我国西北、内蒙以及东北地区修建的铁路大都受到路基冻胀的危害。特别是西北地区的兰新线等每年入冬及开春季节都因路基冻胀或融沉而不得不限速运行。
土体冻胀的机理十分复杂。20世纪30年代之前，人们认为孔隙水结冰膨胀是冻胀发生的主要原因[2]，后来的研究表明未冻水向冻结区迁移形成分凝冰才是冻胀发生的更重要原因[3]。Taber[2]使用凝固后体积收缩的苯代替水作为孔隙液体，同样观测到冻胀现象，这有力地证明了上述论断。
在气温足够低时，地基土会从上向下逐步冻结。在这一过程中，深处未冻土层内的水分会在温度梯度、压力梯度等因素作用下向冻结区域迁移，生成冰透镜体。对冰透镜体的生长过程目前主要有两种理论，即早期的第一冻胀理论(毛细管模型)和后来发展的第二冻胀理论(冻结缘模型)[3]。
毛细管模型最早由Taber提出[2]，假设土体中的孔隙为毛细管，通过克拉佩龙方程[4]描述土体中温度T、冰压力Pi以及水压力Pw之间的平衡关系。
温度下降时，原本的压力平衡会被打破，底部储水层水分在压力作用下向上迁移，冻结为冰透镜体的一部分。该模型认为当温度足够低时，冰透镜体会在压力作用下侵入孔隙，阻断水分迁移，导致冰透镜体无法继续生长。通过将克拉佩龙方程与杨-拉普拉斯方程相结合，可以得到冻胀停止的临界温度。
但是，毛细管模型预测的最大冻胀压力、水分迁移速率与实际情况误差较大，且无法解释冰透镜体为何会形成层状结构，也就无法预测分层冰的形成。因此，Miller[5]和Harlan[6]等提出了冻结缘模型。冻结缘模型最重要的假设是：温度低于临界温度后，冻胀依然可能继续发生。所谓的冻结缘，是指最底部冰透镜体以下，到冻结锋面(即温度为水分冻结温度的面)这一段含有冰水混合物的土体区域。在冻结缘中，水分仍然可以流动，从而使水分向冰透镜体不断迁移，导致冰透镜体不断变厚，冻胀量继续增大。
基于上述两类理论，一些数学模型被用于对冻胀情况下的水分运移规律进行模拟分析。20世纪70年代，Harlan[6]忽略土体中的应力场变化，建立了渗透系数与含冰量相关的水热耦合模型。Miller [7]假设冰透镜体与底部的土体刚接，提出了刚性冰模型。Kornad等[8]从势能的角度分析，从冰的吸水速度入手，建立了分凝势模型。其中，刚性冰模型与分凝势模型都以冻结缘理论为基础。刚性冰模型假设土体中一点的有效应力为0时，该处的冰将土体挤开，就会有新的冰透镜体生成。冰透镜体生长的速率，则与下方水分流向冰透镜体的速率相同。而分凝势模型则从化学中的势能理论出发，根据冰透镜体处水分分凝成冰的势能大小，计算冰透镜体的生长速度。
现有针对冻胀进行模拟分析的模型主要考虑的都是未冻区的液态水向冻结区的迁移，对气态水的迁移考虑较少。近年来人们观察到，在干旱地区地下水位较深的路基中，其不透水道面下方也会有较多水分的聚集，由此推断气态水的迁移不可忽略。李强、姚仰平等[9]将这后一种的现象形象地称为“锅盖效应”，工程界也有人称其为“覆盖效应”[10―11]。盛岱超课题组[10,12―13]也对此进行了研究，指出“锅盖”温度在冰点以下时路基中水分的向上迁移量才较显著。本课题组此前也对覆盖效应进行了建模分析[11]，得到相近的认识。
受冻胀困扰较大的寒区铁路路基，其“覆盖效应”的表现与上有道面覆盖的公路路基并不相同。在浅层土体冻结之前，一般铁路路基表层的渗透性并非很低，表现出弱覆盖或不覆盖的特性。当浅层因雨季降雨等而有一定残留水分时，冬季随冻结发展，路基表层蒸发量逐渐降低，由不覆盖逐渐变为弱覆盖，甚至强覆盖。这种顶部透气透水条件随时间变化从而影响路基内水分迁移状况的现象，我们称之为“时变覆盖效应”。
本文针对“时变覆盖效应”这一现象，考虑路基表层蒸发条件随时间变化，以Miller[7]基于第二冻胀理论建立的刚性冰模型为基础，对入冬时节路基土体内的水分迁移、冰透镜体生成与生长进行模拟分析。为此，结合Philip等 [14]提出的基于质能平衡的非等温水热气耦合运移模型(PDV模型)，建立了相应的数学模型，实现了对时变覆盖效应的数值模拟，通过与试验对比检验了其正确性，进而通过对多种不同条件下的模拟计算对有关影响因素进行了讨论。
1 数值模型
本节基于已有研究给出多孔介质内水分及热量迁移的控制方程，以及冰透镜体生成、生长的计算方法。
1.1 质量与热量守恒方程
路基内土体组分如图1所示，包括土颗粒、水蒸气、冰以及液态水共四相。
图1 路基土体各组分示意图
Fig.1 Schematic diagram of different components in subgrade soil
1.1.1 质量守恒方程
基于Philip等[14]提出的考虑质能平衡的非等温水热气耦合运移模型(PDV模型)，结合表征液态水迁移的Richards方程和表征气态水迁移的Fick定律，取深度方向为坐标轴正向，考虑基质吸力梯度、重力梯度和温度梯度的影响，非饱和区域水气运移的控制方程可写为：
式中：θu为未冻结水体积分数；ρi为冰的密度；ρw为水的密度；θi为冰的体积分数；h为基质吸力对应的水头；KLh为水力梯度导致的等温水力传导率；Kvh为等温气体传导率；KLT为温度梯度导致的非等温水力传导率；KvT为非等温气体传导率；T为绝对温度。
基质吸力对应水头h与含水量之间通过土水特征曲线建立联系。针对本文模拟的初始含水量体积分数0.1以上的土体，忽略干湿循环过程中的滞回现象，采用van Genuchten模型，有关系式：
式中，α、n(m=1-1/n)为需要根据试验确定的参数；Sel为有效饱和度，其定义为：
式中：θs为饱和含水率；θr为残余含水率。
在van Genuchten模型中，等温水力渗透系数可以表示为：
考虑到冻结缘内存在冰水混合区，该区内冰会阻碍水分迁移，需要根据含冰量对等温水力渗透系数进行修正。参考Taylor等[15]的研究内容，修正后的等温水力渗透系数为：
式中：ks为饱和渗透系数；S为有效饱和度。
Satio等[16]给出其余三个系数表达式如下：
式中：Hc为与基质吸力对应的水头；Hr为相对湿度，
为与土性相关的放大系数，一般取7；γ为表面张力系数；γ0为25℃时的表面张力系数。表面张力系数γ与温度关系可以表达为[17]：
Mw为水的摩尔质量；R为气体常量；Dv是水蒸气扩散系数，可表示为：
其中：ρvs为水蒸气密度，根据修正Kelvin公式可表示为[18]：
放大因子ηe表达式为：
式中，fc为土中黏粒的质量分数。当土中不含黏粒时，式(12)最后一项取为0。
2.1.2 热量守恒方程
非饱和土体的热量守恒方程可表示为：
式中：L为冰水相变潜热，取值为334.5 kJ/kg；C(θu)为体积热容，由各相体积热容根据比重加权平均给出[19]。λ(θu)为非饱和土的导热系数，可表示为[20]：
其中，
对粉土λdry=0.5×10-1.5φ，φ为孔隙率。kr经验公式如下：
式中，Tf为土壤中的水分冻结温度。
对土未冻水含量与温度的关系，目前还存在争论。争论的核心是，未冻水含量与初始含水量是否相关。早期研究者通过式(16)所示的克拉佩龙方程建立土中未冻水含量与温度间的关系，并通过试验予以验证[21-24]。
但是，克拉佩龙方程未考虑孔隙中冰点的降低，也未考虑土颗粒表面与未冻水的化学结合。此后随着未冻水含量监测技术的发展，TDR方法与NMR方法的观测结果给出了不同的结论。Tice等[25]使用NMR方法的观测结果表明，初始含水量越高，未冻水含量越高。Watanabe与Wake[26]的试验结果也佐证了这观点。尽管有学者指出，TDR与NMR方法存在相变平衡时间过晚的问题[27]，但未冻水含量与初始含水量相关正成为主流认知。
本文选用McKenzie等[28]2007年建立的半经验公式确定未冻水含量与总含水量的关系：
其中，W为根据观测数据拟合的参数。
1.2 冻胀生长模型
现有的冻胀预测模型主要有水热耦合模型、刚性冰模型和分凝势模型三种。其中，水热耦合模型计算较为简便，但模型过于简化，未考虑水分结冰膨胀对土体的影响，与实际不尽相符，计算精度不高[29]。分凝势模型计算较为准确，但需要用到的参数较多。而实际工程中，土中水分的化学势等参数很难确定，因此难以实用。刚性冰模型考虑了土中冰透镜体形成与发展的过程，与实际较接近[7]，同时模型所需的孔隙率、渗透系数、导热系数参数较易获取，因此本文采用刚性冰模型。
基于冻胀模拟的刚性冰模型，根据土体中有效应力的大小判定冰透镜体的生成。当冻结缘内一点处有效应力为0时，该点处的冰将土挤开，便会有新的冰透镜体形成。根据Bishop[30]1959年给出的非饱和土有效应力公式，考虑冰所产生的应力，可将土中有效应力表示为：
式中：φ表示土的孔隙率；Pw表示水分所承受的分压；Pi表示冰透镜体所承受的分压。结合式(16)克拉佩龙方程，消去式中的冰压项，可得到有效应力表达式：
式中：σ为总应力，是上覆压力与土体自重应力之和；T0为水的冻结温度。
一旦冰透镜体生成，就将阻碍水分迁移。在冻胀发展过程中，水分迁移的边界面将随着冰透镜体的生成发生变化。而通过冰透镜体下方的质量与热量守恒方程，可以解出水分流向最下方冰透镜体的速率，该速率即为冰透镜体的生长速率。
1.3 蒸发速率
在对时变覆盖效应模拟时，顶部蒸发速率对土内水分迁移状况有较大影响。土体不饱和时，常用与风速呈线性相关的半经验蒸发公式[16]：
式中：u为风速；es为空气湿度；ed为土体顶部湿度；a和b为试验拟合参数。
1.4 数值计算方法
在实际的铁路路基工程中，水气除竖向运移外，还会向两侧坡面进行运移，缓解路基顶部水气阻滞效应。当降雨发生时，路基两侧坡面上水分又会向路基内部迁移，加剧路基土体水气阻滞效应。同时，路基两侧空气与中间的热量传输也会对土体内的温度分布产生影响。因铁轨下方土体距两侧坡面距离较远，受两侧水热变化影响较小，本文将铁轨下方土体近似为一维土柱，使用一维模型进行模拟。
待求解的偏微分方程式(1)、式(13)连同冰水含量关系式(17)写为：
式(21)与平衡方程及有效应力方程式(18)联立求解可实现土体冻结过程中水气迁移及冻胀过程的模拟。有限差分表示待求函数的空间偏导数，从而将偏微分方程转化为关于时间的常微分方程[31]，进而通过二阶总变差最小(TVD)方法对其求解。
将关于时间的常微分方程记为：
则二阶TVD方法公式为：
2 模型验证
这里先对基本的水气迁移试验和冻结试验进行模拟计算，以验证模型的合理性。这里的两个试验均来自于文献[1]。
2.1 冻结融化试验
冻结试验中，将粉土做成高度10 cm的圆柱体，然后将土体冻结到-1℃。融化过程中，保持顶部温度为-1℃，底部温度为1℃。因为试样周围设置保温层，可以简化为一维模型。由于模拟过程为冻结试样融化的过程，因此不考虑冻胀的影响，但需要考虑冰对水分迁移的阻滞作用。
2.1.1 参数取值
粉土的参数取值如表1所示，数据取自文献[11]。
2.1.2 模拟结果
将总含水量定义为θ=θu+=(ρi/ρw)θi ，即将冰折算为等质量的水。模拟得到的总含水量与实测总含水量对比如图2所示。
图2显示，模拟结果与实测结果较为接近，都表现出在试样底部，因为水分向下迁移，存在含水量峰值的特点。
表1 粉土经验参数表
Table 1 Empirical parameters of silt
符号 取值 参数说明θr 0.05 残余水分体积分数φ 0.42 孔隙率θs 0.42 饱和水分体积分数m 0.26 经验参数Cs 724 J/(kg·K)多孔介质热容Cw 4200 J/(kg·K))水的热容Ci 2100 J/(kg·K)冰的热容α 1.6 m-1 经验参数W 0.47 经验参数Tf-0.47 ℃ 粉土冻结温度ks1×10-7 m/s 饱和土渗透系数λi 2.2 W/(m·K)冰的导热系数λw 0.54 W/(m·K)水的导热系数λs 2.52 W/(m·K)粉土的导热系数Mw18 g/mol 水的摩尔质量
图2 冻结试验模拟结果对比图
Fig.2 Comparison of simulation results and measured results of frozen experiment
2.2 冻胀试验
试验使用的土样为兰州黄土，其参数取值与表1相同。试样为直径15 cm、高12 cm的土柱，初始温度为1℃。冻结过程中，试样顶部温度呈直线下降，在50 h内从-1℃降到-4℃。底部使用恒定水位的马廖特瓶与土柱底面连接，以提供水源。
模拟得到的冻胀量对比如图3所示。
图3 兰州黄土冻胀试验模拟结果对比图
Fig.3 Comparison of simulation results and measured results of Lanzhou loess frost heave experiment
图3显示，在降温速度恒定的情况下，冻胀过程中，初始冻胀速度较快，此后冻胀速度逐渐降低。模拟值与实测值较为接近，但与实测值相比，模拟的冻胀速度稍快一些。
3 “时变覆盖效应”模拟分析
3.1 模型参数
模型选取为3 m高的一维粉土柱，底部含水量体积分数恒定为0.2。顶部分别选取为蒸发边界条件和覆盖边界条件，以分别模拟时变覆盖效应及覆盖效应 。初始含水量体积分数设定为全域0.2。
粉土的经验参数取值与表1相同，土体顶部温度边界条件如下：
式中：Ttop单位为℃。式(24)反映顶部日均温度随季节按余弦曲线周期变化、每日温度随时间周期变化。模拟时间为92天时，日最低气温降到0℃以下。土体底部温度恒定为15℃，土体初始温度选取为全域15℃。
蒸发公式中，土体顶部空气湿度接近正弦变化，基于本课题组于兰新铁路武威段实测的空气湿度数据，选取拟合公式如下：
式中，t为每日的时刻，单位为小时。
蒸发公式中的经验参数选取为a=0.40，b=0.60，对应风速单位为m/s，蒸发速率单位为mm/d。模拟风速基于本课题组于兰新铁路武威段实测的风速数据，选取为1.5 m/s。
模拟时段从夏季气温最高时开始至冬季温度最低时。
3.2 模拟结果
模拟得到的冻胀量以及土体内的冰透镜体分布如图4所示。
图4 覆盖/不覆盖土体冰透镜体生成图
Fig.4 Ice lens in covered/uncovered soil
从图4可以看出，全时覆盖情况下土体冻结锋面更高，但冻胀量更大。这是因为顶部覆盖的情况下，顶部的初始含水量就比较高，导致冻胀量更大。但二者并不存在数量级上的差距，不覆盖情况下的冻胀仍然可能造成工程危害。
模拟得到的非覆盖土体顶部水分日蒸发量变化如图5所示，非覆盖土体底部水分迁移量如图6所示。
结合图5和图6可以看出，在时变覆盖条件下，夏季蒸发量较高，这导致冻胀发生前，顶部水分含量较低。随着气温降低，土体内水分在温度作用下向上迁移，顶部含水量升高，蒸发量反而有所提高。而当冻胀现象开始发生后，冰透镜体阻碍了顶部水分迁移，蒸发量迅速下降到接近0毫米每日。与此同时，底部水分迁移量迅速增加，表现为与完全覆盖类似的特点。
图5 非覆盖土体顶部水分日蒸发量变化图
Fig.5 Evaporation rate at the top of uncovered soil
图6 非覆盖土体底部水分迁移速率变化图
Fig.6 Water migration rate at the bottom of uncovered soil
4 影响参数分析
这里通过计算对比分析土的性质、温差大小、初始含水量大小对“时变覆盖效应”的影响。
4.1 土体性质的影响
为分析不同性质土体时覆盖效应的强弱，分别选取砂土、粉土、粘土三种土，对“时变覆盖效应”进行模拟分析。
参考van Genuchten模型给出的典型砂土、粉土和粘土的参数[32]，表1有变化的参数列于表2，模拟得到的冻胀量对比如表3所示。
表2 砂土、粉土、粘土参数表
Table 2 Parameters of sand, silt and clay
符号 砂土 粉土 粘土n 2.68 1.37 1.09 α/m-1 14.5 1.6 0.8 B 0.61 0.47 0.56 ks/(m/s) 1×10-5 1×10-7 1×10-9
表3 砂土、粉土、粘土时变覆盖效应冻胀量表
Table 3 Frost heave of sand, silt and clay under‘time-varying canopy effect'
土体类型 砂土 粉土 粘土冻胀量/mm 0 110 0
可以看出，在砂土以及黏土中，并没有发生冻胀现象。这与实验(实际情况)中，初始含水量远低于饱和的情况下，砂土与粘土不发生冻胀的情况相符[32]。
由于粉土中发生了冻胀现象，粉土土柱总高度高于砂土和粘土。为便于比较，去掉粉土中的冰透镜体，将剩余土柱的含水量体积分数与砂土、黏土中的含水量体积分数进行对比。土体中的冰折算为等质量的水的体积分数。从而得到砂土、粉土、粘土最终含水量体积分数对比图如图7所示。
图7 不同土性含水量体积分数对比图
Fig.7 Comparison of water volumetric content for different soil types
由图7可见，砂土中顶部含水量降低，黏土中含水量基本不变，而粉土则在顶端发生较为明显的水分聚集现象。
模拟得到的砂土、粉土、粘土顶部水分日蒸发量以及底部水分迁移速率随时间的变化分别如图8、图9所示。
图8 不同土体顶部水分日蒸发量变化图
Fig.8 Comparison of evaporation rate for different soil types
从图8可以看出，砂土与粘土的蒸发量都随时间逐渐降低，当入冬之后，蒸发量都接近于0。从图9可以看出，砂土中水分以向下迁移为主，而粘土中水分迁移困难，含水量变化较小，底部水分迁移量很少。
结合表3与图7～图9，可以看出，砂土中水分以向下迁移为主，顶部含水量较低，不会发生冻胀现象；粘土中水分迁移较为困难，结冰后缺少水分补给，同样难以形成冰透镜体发生冻胀。而粉土中，在温度梯度作用下水分会发生向上迁移，同时结冰后水分补给也较为容易，最易发生冻胀灾害。
图9 不同土性底部水分迁移量变化图
Fig.9 Comparison of water migration for different soil types
4.2 温差大小的影响
为分析不同温差大小对时变覆盖效应的影响，保持顶部温度不变，分别设底部温度保持为10℃、15℃、20℃，进行模拟分析。其余模拟参数与第4节相同。模拟得到的冻胀量对比如表4所示，可以看出，底部温度越高，温差越大，冻胀量也越大。
表4 不同温差大小时变覆盖效应冻胀量表
Table 4 Frost heave of soil under different temperature differences
底部温度/(℃)10 15 20冻胀量/mm 83 110 122
模拟得到的顶部水分日蒸发量变化图如图10所示；底部水分迁移速率变化情况如图11所示。
图10 不同温度顶部水分日蒸发量变化图
Fig.10 Comparison of evaporation rate under different temperature differences
从图11可见，底部温度越高，温差越大，底部水分迁移量越大。温度梯度的加大，也是导致冻胀量变大的主要原因之一。同时，因为水分迁移量大，顶部水分积聚较多，蒸发量较大。但是当冬季冻胀开始发生后，蒸发量还是都降低到0附近。
图11 不同温度底部水分迁移量变化图
Fig.11 Comparison of water migration under different temperature differences
4.3 初始含水量的影响
分别选取初始含水量对应的体积分数为0.1、0.2、0.3、0.4，其余参数与第4节相同。模拟得到的冻胀量对比如表5所示。从表5可以看出，初始含水量越高，冻胀量越大。
表5 不同初始含水量时变覆盖效应冻胀量表
Table 5 Frost heave of soil with different initial water content
初始含水量 0.1 0.2 0.3 0.4冻胀量/mm 42 110 153 221
模拟得到的顶部水分日蒸发量变化如图12所示；底部水分迁移速率变化如图13所示。
从图12可见，初始含水量越大，蒸发量越大。同时，由图13可见，在冻胀发生前，初始含水量较高或者较低都会导致底部水分迁移量较少，初始含水量适中时水分向上迁移量最大。因为初始含水量较低时，土体中可供迁移的水分较少；而初始含水量较高时，重力降水作用加强，水分迁移量同样较少。这与文献[33]针对完全覆盖情况下的分析结论一致[33]。但当冻胀发生后，初始含水量越高，底部水分迁移量越大。这是因为液态水不断流向冰透镜体冻结，较高的含水量并不会阻碍水分迁移。
图12 不同初始含水量顶部水分日蒸发量变化图
Fig.12 Comparison of evaporation rate with different initial water contents
图13 不同初始含水量底部水分迁移量变化图
Fig.13 Comparison of water migration with different initial water content
5 结论
路基冻胀是铁路、公路运营及维护中很值得重视的问题之一。本文就寒区路基内水分迁移及冻胀进行分析，针对浅层冻结状况变化使顶部透水透气性发生改变从而影响路基内水分迁移集聚的现象，提出“时变覆盖效应”概念，，并提出相应的计算分析模型，实现了时变覆盖条件下路基的冻胀分析。经计算分析得到的主要认识如下：
(1)在铁路路基中，冬季顶部水分冻结形成冰透镜体，成为覆盖层，会阻碍水分蒸发，导致路基顶部水分集聚，进而加大土体冻胀。这种顶部覆盖条件随时间变化的现象，称为“时变覆盖效应”。
(2)与“覆盖效应”相比，“时变覆盖效应”冻胀量相对较小，但仍会造成一定的工程危害。
(3)与砂土和粘土相比，粉土中最易因“时变覆盖效应”发生冻胀灾害。
(4)土体上下温差越大，底部水分迁移量越大，“时变覆盖效应”下的冻胀量也越大。
(5)初始含水量越高，冻结开始后的底部水分迁移量越大，“时变覆盖效应”下的冻胀量也越大。
“时变覆盖效应”的模拟为一复杂的多场耦合问题，且需考虑冻胀过程中的冰透镜体生成与生长。本文建立模型对这一现象进行了模拟，揭示了铁路路基等不覆盖土体顶部蒸发边界条件变化的现象，及其对土体内水分迁移的影响。本文尚未考虑上覆压力等影响因素对“时变覆盖效应”的影响，同时模型为一维，尚不能精确反映实际铁路路基的不均匀冻胀现象。
[1]徐斅祖, 王家澄, 张立新.冻土物理学[M].北京: 科学出版社.2001 Xu Xiaozu, Wang Jiacheng, Zhang Lixin.Frozen soil physics [M].Beijing: The Science Publishing Company.2001.(in Chinese)
[2]Taber S.The mechanics of frost heaving [J].Journal of Geology, 1930, 38(4): 303―317.
[3]Peppin S L, Style R W.The physics of frost heave and ice-lens growth [J].Vadose Zone Journal, 2013, 12(1):1―12.
[4]Wark K J.Generalized thermodynamic relationships.Thermodynamics, 5th ed [M].New York: McGraw-Hill,Inc., 1988.
[5]Miller R D.Freezing and heaving of saturated and unsaturated soils [J].Highway Research Record, 1972,393: 1―11.
[6]Harlan R L.Analysis of coupled heat–fluid transport in partially frozen soil [J].Water Resource Research, 1973,9: 1314―1323.
[7]Miller.Frost heaving in non-colloidal soils [C].Proceedings of the 3rd International Conference on Permafrost,.Edmonton, AB, Canada.10–13 July, 1978.Natl.Res.Counc.of Canada, Ottawa, 1978, 1: 708―713.
[8]Konrad J, Duquennoi C.A Model for water transport and ice lensing in freezing soils [J].Water Resources Research, 1993, 29(9): 3109―3124.
[9]李强, 姚仰平, 韩黎明, 等.土体的“锅盖效应”[J].工业建筑, 2014, 44(2): 69―71.Li Qiang, Yao Yangping, Han Liming, et al.Pot-cover effect of soil [J].Industrial Construction, 2014, 44(2):69―71.(in Chinese)
[10]滕继东, 贺佐跃, 张升, 等.非饱和土水气迁移与相变: 两类“锅盖效应”的发生机理及数值再现[J].岩土工程学报, 2016, 38(10): 1813―1821.Teng Jidong, He Zuoyue, Zhang Sheng, et al.Moisture transfer and phase change in unsaturated soils: physical mechanism and numerical model for two types of“canopy effect” [J].Chinese Journal of Geotechnical Engineering.2016, 38(10): 1813―1821.(in Chinese)
[11]宋二祥, 罗爽, 孔郁斐, 等.路基土体“锅盖效应”的数值模拟分析[J].岩土力学, 2017, 38(6): 1781―1788.Song Erxiang, Luo Shuang, Kong Yufei, et al.Simulation and analysis of pot-cover effect on moisture transport in subgrade soil [J].Rock and Soil Mechanics, 2017, 38(6):1781―1788.(in Chinese)
[12]Zhang S, Teng J, He Z, et al.Importance of vapor flow in unsaturated freezing soil: a numerical study [J].Cold Regions Science & Technology, 2016, 126(6): 1―9.
[13]Sheng D, Zhang S, Yu Z, et al.Assessing frost susceptibility of soils using PC Heave [J].Cold Regions Science & Technology, 2013, 95(11): 27―38.
[14]Philip J R, De Vries D A.Moisture movement in porous materials under temperature gradient [J].Trans am Geophys Union, 1957, 38(2): 222―232.
[15]Taylor G S, Luthin J N.A model for coupled heat and moisture transfer during soil freezing [J].Canadian Geotechnical Journal, 1978, 15(4): 548―555.
[16]Saito H, Simunek J, Scanlon B R, et al.Numerical analysis of coupled water, vapor and heat transport in the vadose zone using HYDRUS [J].Vadose Zone Journal,2006, 5(2): 784―800.
[17]Simunek, Saito, Sakai, et al.The hydrus-1D software package for simulating the one-dimensional movement of water, heat, and multiple solutes in variably-saturated media [EB].Departmeat of Environmental Sciences,University of California Riverside, 2005.
[18]Gaylon S Campbell.Soil physics with BASIC: transport models for soil-plant systems [M].Amsterdam: Vol.14Elsevier, 1985.
[19]夏锦红, 陈之祥, 夏元友, 等.不同负温条件下冻土导热系数的理论模型和试验验证[J].工程力学, 2018,35(5): 109―117.Xia Jinhong, Chen Zhixiang, Xia Yuanyou, et al.Theoretical model and experimental verification on thermal conductivity of frozen soil under different negative temperature conditions [J].Engineering Mechanics, 2018, 35(5): 109―117.(in Chinese)
[20]原喜忠, 李宁, 赵秀云, 等.非饱和(冻)土导热系数预估模型研究[J].岩土力学, 2010, 31(9): 2689―2694.Yuan Xizhong, Li Ning, Zhao Xiuyun, et al.Study of thermal conductivity model for unsaturated unfrozen and frozen soils [J].Rock and Soil Mechanics, 2010, 31(9):2689―2694.(in Chinese)
[21]Williams.Suction and its effects in unfrozen water of frozen soils.In: Permafrost, proceedings of an international conference [C].Washington, DC: National Academy of Sciences, 1966: 225―229.
[22]Williams P J.Properties And Behavior Of Freezing Soils[M].Norwegian: Norwegian Geotechnical Institute,1967.
[23]Anderson, Tice A R.Predicting unfrozen water contents in frozen soils from surface area measurements [J].Highway Research Record, 1972, 393(2): 12―18.
[24]Fisher E A.The freezing of water in capillary systems[J].Journal of Physical Chemistry, 2002, 28(4): 360―367.
[25]Tice A R.Determination of unfrozen water in frozen soil by pulsed nuclear magnetic resonance [C].Permafrost,Proceedings of the Third International Conference, 1978,1: 150―155.
[26]Watanabe K, Wake T.Measurement of unfrozen water content and relative permittivity of frozen unsaturated soil using NMR and TDR [J].Cold Regions Science &Technology, 2009, 59(1): 34―41.
[27]Kurylyk B L, Watanabe K.The mathematical representation of freezing and thawing processes in variably-saturated, non-deformable soils [J].Advances in Water Resources, 2013, 60(60): 160―177.
[28]Mckenzie J M, Voss C I, Siegel D I.Groundwater flow with energy transport and water–ice phase change:Numerical simulations, benchmarks, and application to freezing in peat bogs [J].Advances in Water Resources,2007, 30(4): 966―983.
[29]武建军, 韩天一.饱和正冻土水-热-力耦合作用的数值研究[J].工程力学, 2009, 26(4): 246―251.Wu Jianjun, Han Tianyi.Numerical research on the coupled process of the moisture-heat-stress fields in saturated soil during freezing [J].Engineering Mechanics, 2009, 26(4): 246―251.(in Chinese)
[30]Bishop A W.The principle of effective stress [J].Teknisk Ukeblad, 1959, 106(39): 859―863.
[31]毛卫南, 刘建坤.不同离散化方法在正冻土水热耦合模型中的应用[J].工程力学, 2013, 30(10): 128―132.Mao Weinan, Liu Jiankun.Different discretization method using in coupled water and heat transport mode for soil under freezing conditions [J].Engineering Mechanics, 2013, 30(10): 128―132.(in Chinese)
[32]Simunek J, Huang K, van Genuchten M.The HYDRUS-ET software package for simulating the one-dimentional movement of water, heat and multiple solutes in variably-saturated media, version 1.1 [M].Bratislava: Inst.Hydrology Slovak Acad.Sci, 1997.
[33]姚仰平, 王琳.影响锅盖效应因素的研究[J].岩土工程学报, 2018, 40(8): 1373―1382.Yao Yangping, Wang Lin.Research on the influence factors on the “Pot-cover effect” [J].Chinese Journal of Geotechnical Engineering, 2018, 40(8): 1373―1382.(in Chinese)