陈玉燕，马学宁，张 正
（1.兰州交通大学土木工程学院，甘肃兰州 730070；
2.中机中联工程有限公司，重庆 400039；
3.甘肃省交通规划勘察设计院股份有限公司，甘肃兰州 730010）
冻土是一种含冰的特殊土体，是由固体土颗粒、冰、孔隙水及孔隙气体组成的多相体系，其温度一般低于0℃。冻土广泛分布于我国的西部及北部地区，约占我国国土面积的75%，多年冻土面积约为2.068×106km2，占我国国土面积的21.5%，占世界冻土面积的10%［1］，在世界屋脊——青藏高原地区分布着我国绝大多数多年冻土。随着我国一带一路战略的持续推进和中西部地区的进一步发展，将有大量的工程建设在冻土地区开展，如川藏铁路、兰张三四线、敦格铁路等的建设，这些铁路所经过的地区大部分为冻土地区。
无论是多年冻土地区的道路工程还是季节冻土地区的道路工程，都会由于土体的冻胀融沉产生一系列的病害，如由于冻胀融沉作用引起的支挡结构的破坏、公路路面的开裂、铁路轨面不平顺等，这些都是冻土地区工程建设所必须要解决的问题，如果这些问题处理不当将会造成巨大经济损失［2］。冻土的性质除了受颗粒大小、矿物组成、含水量及压实度等影响之外还主要受到含冰量的影响，而冻土的含冰量又与温度与含水量密切相关，路基中的温度场和水分场的分布、状态及其变化规律是影响路基冻胀与否的主要原因［1］。在冻土地区的道路工程中，预测路基中不同深度处的温度场及水分场的状态及其变化规律是一个十分重要且又相当复杂的问题，因为其影响因素并不单一，而是各种因素相互作用的结果，也是判断道路工程是否会发生冻害的重要条件，因此，对冻土温度场及水分场的状态及其变化规律进行深入研究显得尤为重要［1-3］。
随着冻土地区工程项目的逐渐增多，为了减少冻土地区路基工程的病害，国内外学者从模型试验、数值模拟及理论分析等方面对冻土的冻胀融沉机理及冻土中温度场及水分场的的变化等进行了大量研究［4-14］，并取得了一定的成果，提出了许多预防和处理路基冻害的措施，但目前的相关研究大多都只是进行了单纯的数值模拟或者模型试验，并没有将两者结合起来，通过相互验证、相互对比来研究冻土地区路基温度场的变化规律。文中以兰州地区典型黄土为研究对象，考虑在有限区域内地面温度无明显差异，假设相邻区域间无热能交换，将实际自然条件下的冻结条件在模型试验中简化为一维冻结条件，采用数值模拟结合室内试验的方法研究了不同环境温度和土体含水率对地基温度场的影响。
1.1 试验装置及材料
本次模型试验所用到的装置是一款自行设计的一维封闭冻胀试验系统，该系统由试验模型箱、温度控制系统、测量系统等组成。试验装置示意图如图1所示。
模型长、宽、高分别为80 cm、60 cm、80 cm。温度控制装置使用的是XT5704LT-BF28低温恒温循环器，制冷循环液使用的是无水乙醇，温度输出范围为-30℃～＋50℃，精度为0.2℃，温度传感器测量精度为0.1℃，工作范围为-50℃～＋200℃，数据采集使用LCXAF巡检仪及计算机。在土样内部左右两侧各布置一组温度传感器，每组共有5个，由上至下每隔12 cm布设一个。试验所用土体为兰州地区的典型黄土，其基本物理性质指标为：比重GS=2.71；
塑限WP=15.21%；
液限WL=28.38%；
塑性指数IP=13.17；
最优含水率WOPC=13.6%，最大干密度ρdmax=1.89 g/cm3。
1.2 试验方案和步骤
图1 模型试验系统示意图Fig.1 Schematic diagram of model test system
表1 冻结试验方案Table 1 Freezing test scheme
试验所模拟的是封闭、无外荷载条件下的自上而下的单向冻结，试验方案见表1。试验步骤为：（1）按设计要求制作模型箱，并在其绝热边界填充保温棉，土体导热系数为1.0～1.6 W/（m·K），保温棉导热系数为0.03～0.035 W/（m·K），土体导热系数约为保温棉导热系数的40倍，故认为试样除冷端外其他边界为绝热边界，模型试验符合一维冻结条件。为模拟自然边界条件，在模型箱底部铺设5 cm中砂；
（2）将土体分层填入模型箱，并击实至设计压实度0.85，分层填筑过程中对土体进行补水，使其含水率为设计含水率（14%、16%、18%），并将温度传感器埋设在相应位置，并记录其初值，填筑完成后将试验箱静置72 h，使土中的水分分布均匀；
（3）将冷端温度由环境温度（14℃）调至试验设计的冻结温度并保持恒定，当土体温度达到设计值时停止冻结，将温度调为设计的融化温度使其融化，测量温度取同一深度两组测量值的平均值。
因本研究中所涉及的冻土路基中的温度场和水分场是动态变化的，且两者之间的变化是相互作用、相互影响的，是复杂的耦合关系，COMSOL Multiphysics是一款基于有限元的具有强大的多物理场耦合能力和强非线性微分方程组求解能力的数值仿真软件，可较好的模拟冻土中水分场和温度场的耦合作用［15-16］。文中利用COMSOL Multiphysics中的理查兹方程模块和多孔介质模块对冻土中温度场和水分场的变化进行模拟的。
土体在冻融循环过程中伴随着水分相变的导热，基于含相变潜热项建立的热传导方程为：
式中：C为土体的体积比热容，J·m-3·K-1；
λ为土体的导热系数，W·m-1·K-1；
L为水的相变潜热，J·kg-1，取值为3.35×105J/kg；
t为时间，s；
T为温度，℃；
∇为微分算子，对于二维问题为［∂/∂x，∂/∂z］；
θi为冻土中冰的体积含量，%；
ρi为冰的密度，kg·m-3。
冻融状态下路基内始终存在未冻水，而土体中未冻水的迁移遵循达西定律［15-17］。根据有相变项的Richard方程，非饱和冻土中未冻水迁移的微分方程为：
式中：θu为冻土中未冻水的体积含量，%；
Kg（θu）为非饱和冻土重力方向的导水系数，单位为m·s-1，表示单位水力梯度下的流量；
D（θu）为非饱和冻土的水分扩散系数，单位为m2·s-1，其主要由于含水量的差异引起。
数值计算时网格划分采用超细化网格，计算步长为0.05 h，相对容差为0.000 1，模型四周和底部设置为热绝缘边界。为方便与试验结果进行比较，文中针对冷端温度为-10、-15、-20℃工况下进行数值模拟。
3.1 模型试验与数值模拟可靠性验证及温度沿深度的分布
通过对比模型试验与数值模拟的结果，对模型试验与数值模拟的可靠性进行相互验证。文中以含水率为16%、冻结时间为250 h、冻结温度分别为-10℃、-20℃进行验证，对比分析数值模拟和试验过程中土体温度沿深度的的分布规律，对比结果如图2所示，可以看出，数值计算结果与试验结果吻合度较高，说明所建数值模型和所取计算参数是可靠。
从图2（a）可见，室内试验和数值模拟所得到的曲线在深度为12～24 cm和36～64 cm处基本平行，在深度为24～36 cm处两曲线相交；
在12～24 cm处曲线斜率较小，即温度梯度（同一时刻相邻两测点温差与距离之比）较大，在同一深度处数值模拟所得温度略小于试验所得温度；
在36～64 m处曲线斜率较大，即温度梯度较小，试验曲线越靠近-0.7℃斜率越小，模拟曲线越靠近-0.5℃斜率越小，且在同一深度处试验温度小于数值模拟温度，由此可知，本试验中水的相变发生在-0.7℃左右，即试验中所涉及的水的冻结温度可能在-0.7℃左右，并非数值模拟中所设置的-0.5℃，所以才出现上述现象。观察图2（b）发现，在冻结温度为-20℃，冻结时间为250 h时，整个监测范围内土体的温度均在-13～-4℃之间，且模拟曲线与试验曲线基本平行，温度梯度随着深度的增加逐渐减小，变换规律与图2（a）中还未发生相变时的变化规律相符。对比图2（a）、（b）发现，图2（b）中两曲线间距大于图2（a）中两曲线间距（这是由于冻结温度不同造成），且图2（b）中同一时刻同一深度处试验温度在各个深度处均小于模拟温度，这与图2（a）中再降温阶段的模拟温度小于试验温度相反，出现这种现象的原因是在冻结温度为-10℃时试验所得曲线与数值模拟所得曲线间距较小，且试验所得平稳变化阶段持续时间略大于数值模拟所得平稳变化阶段时间，导致两曲线相交，使后期数值模拟所得温度小于试验所得温度。
图2 土体温度沿深度的分布对比图Fig.2 Distribution diagram of soil temperature along depth
3.2 土体温度变化规律
图3为室内试验和数值模拟所得到的冻结温度为-15℃、土体含水率为16%时各测温点在一个冻融循环周期内温度变化曲线。可以看出，各测点模拟曲线变化趋势与实测曲线变化规律基本一致，此次试验过程中未出现过冷阶段（过冷阶段：当土体温度下降后土中水处于负温但无冰晶存在的阶段称为过冷阶段），造成这种现象的原因可能是试验土样含水率太小且冻结速度过快所致。冻结过程中各测温点温度变化大致可分3个阶段，即降温阶段、平稳变化阶段和再降温阶段。
图3 各测温点温度随时间变化曲线（单位：cm）Fig.3 Temperature curve of each temperature measurement point with time（Unit：cm）
从0～100 h为降温阶段，这一阶段土体未冻结，土中水分状态未发生变化，土体还处于融土状态。降温阶段土体温度梯度最大，即这一阶段土体温度变化速率最大，且距冷端越近温度梯度越大，这一现象模拟曲线比试验曲线表现的更加明显。
因冻结过程中未出现过冷阶段，则第二阶段为平稳变化阶段，这一阶段持续的时间因深度不同而有所差异。平稳变化段土体中的液态水开始变为固态冰，水转化为冰的过程中释放出大量潜热，与冷端输入的冷能相互抵消，达到一个相对平衡的状态，因文中所涉及的水并非纯水，其冻结温度小于0℃，所以图中平稳变化阶段的曲线在在0℃之下，符合实际情况。平稳变化阶段持续时间的长短与到冷端的距离有明显关系，距冷端越远的点进入第二阶段时间越迟且持续时间越长，距冷端越近的点进入第二阶段越早且持续时间越短。
当冷端输入的冷能与相变释放的潜热完全抵消后土体温度变化进入第三阶段。第三阶段为再降温阶段，这一阶段中土体中大量水都转换为冰，再降温阶段的降温速率大于第二阶段小于第一阶段，这一阶段各测点温度曲线近似平行，且各曲线间距基本相等，所以在这一阶段同一时刻相邻各测点温差基本相等，约为1.2℃。
融化阶段各测温点温度变化曲线与冻结阶段有一定的对称性；
融化阶段初期各测点温度迅速升高至-0.5℃左右后进入平稳变化阶段；
距离上表面越近温度变化越快并越早进入平稳变化阶段；
与冻结阶段相反，融化阶段出现平稳变化阶段的原因是此时土中水分由固态冰逐渐转化为液态水，此相变吸热，在相变期间传输进土体的热量全部被冰吸收，所以土体温度保持不变。
3.3 不同冷端温度对土体中温度场的影响
图4为土体含水率为16%时不同冻结温度下深60 cm处各测温点温度随时间变化曲线，对比同一深度处测温点在不同冻结温度下的变化曲线，发现冻结阶段初期，同一深度处各曲线倾斜程度基本一致，即在冻结初期降温速率受冻结温度影响不大。随着冻结阶段的继续发展，各温度场开始向平稳变化阶段过渡，从图中可以明显看出，冻结温度越低，过渡时间越短，过渡越剧烈；
还可以看出，冻结温度越低，土体处于平稳变化阶段的时间就越短；
三阶段中平稳变化阶段降温速率最小。当进入再降温阶段后，冻结温度较小的点的温度变化出现滞后现象，且随着时间的推移这种现象越加明显。在再降温阶段，不同冷端温度下试样的温度变化速率不同，冻结温度越低降温速率越大。对比数值模拟和试验曲线发现，同一时刻数值模拟所得温度稍稍滞后于试验所得温度，且数值模拟的平稳变化阶段温度略高于试验所得平稳变化阶段温度，这与前期所得结论（试验中所涉及水的冻结温度低于数值模拟所设置的温度）相符。造成滞后的原因可能是在冻结温度设置时，为让数值模拟结果很好的收敛，数值模拟的冻结温度是由环境温度通过20 h的线性变化到达冻结温度的，而在试验时，是直接将温度控制系统的温度调为冻结温度的，即在试验当中，模型箱内温度会在短时间由环境温度降为冻结温度，所以数值模拟温度变化曲线出现滞后现象。除上述2点区别外，数值模拟所得曲线与试验所得曲线其他变化规律基本一致。
图4 不同冻结温度下深60 cm处温度随时间变化曲线Fig.4 Temperature change curve with time at depth of 60 cm at different freezing temperatures
图5 不同冻结温度下土体温度沿深度变化规律Fig.5 The temperature variation of soil with depth under different freezing temperatures
图5为数值模拟得到的含水率为16%、冻结温度为-15℃和-20℃条件下，不同深度处土体温度在0～200 h内的变化规律对比曲线。由图可知，不同冻结温度下同一深度处温度场变化规律一致，在0℃前后同一深度处温度曲线基本平行，同一时刻同一深度处冻结温度为-15℃的温度变化曲线滞后于冻结温度为-20℃的温度变化曲线，且在0℃前的滞后现象较0℃后的滞后现象明显。观察图5我们还可以发现，当冻结时间小于200 h时，深度大于65 cm时，土体的温度基本在零摄氏度左右或者大于0℃，冻结温度对深度大于65 cm的土层影响很小，所以我们可认为，当冻结时间小于200 h时，冻结温度在-20度以内冻结深度不会超过地下65 cm
3.4 不同含水率对土体中温度场的影响
为更好比较不同含水率下土体温度变化情况，取冻结时间为0～300 h内的实验曲线及数值模拟曲线进行研究。图6为试验所得冻结温度为-20℃、含水率分别为14%、16%及18%条件下，深48 cm处各测温点温度随时间变化曲线。观察图6发现，在降温阶段，相同冻结温度不同含水率条件下，同一深度处各点降温速率基本一致；
在平稳变化阶段，随着含水率越的变小，平稳变化阶段持续的时间有逐渐变大的趋势；
在再降温阶段两曲线发展趋势有了明显的不同，含水率小的曲线温度变化出现滞后现象，同一深度处同一时刻含水率越大的点温度越低。
图6 不同含水率下深48 cm处温度随时间变化曲线Fig.6 Temperature change curve with time at depth of 48 cm under different water content
文中利用室内试验和数值模拟对不同环境温度下土体温度场变化进行了研究，结果表明：
（1）冻结过程中土体内温度场变化只表现出降温阶段、平稳变化阶段、再降温阶段，过冷阶段表现不明显。降温阶段温度梯度较大，温度变化速率较快。平稳变化阶段持续的时间与到冷端距离与冻结温度有明显关系，即同一冻结环境温度下，距冷端距离越远，平稳变化阶段持续时间越长；
同一深度处冻结温度越低，平稳阶段持续时间越短。
（2）试验中水的相变发生在-0.7℃左右，即试验中所涉及的水的冻结温度在-0.7℃左右。
（3）在不同的冻结温度下，冻结温度对冻结初期降温速率影响不大，但对平稳变化阶段和再降温阶段影响较大，冻结温度越低，土体处于平稳变化阶段的时间就越短，在再降温阶段，冻结温度较小的点的温度变化出现滞后现象，随着时间的推移这种现象越加明显，且在再降温阶段，不同冷端温度下试样的温度变化速率不同，冻结温度越低降温速率越大。
（4）环境温度变化对土体内温度场的影响较小。在降温阶段，在相同冻结温度不同含水率条件下，同一深度处各点降温速率基本一致；
在平稳变化阶段，随着含水率的变小，平稳变化阶段持续的时间有逐渐变大的趋势；
在再降温阶段两曲线发展趋势有了明显的不同，含水率小的曲线温度变化出现滞后现象，同一深度处同一时刻含水率越大的点温度越低。