寒区路基水-热耦合理论模型及变形规律

Theoretical model of water-heat coupling and deformation law of roadbed in cold area

摘要

为探究寒区路基在水-热耦合作用下的变形机理，以鹤哈高速（G1111）苔青至伊春段（62.44 km）为对象，构建了二维水热迁移耦合模型，结合现场监测与数值模拟方法，揭示了冻融循环过程中路基温度场、水分场与变形场的相互作用规律。基于相变传热与非等温水分迁移理论，采用Abaqus软件进行仿真分析，在路基横断面上自左侧路肩至右侧路肩按1.25 m间距均匀布设9个综合监测点，量化路基表层冻胀与融沉阶段的变形特征。结果表明：冻结期模拟路段路基表层最大冻胀位移可达8.48 mm，横截面变形呈抛物线分布；春季融沉速率（0.12 mm/d）显著高于冻胀期（0.07 mm/d），且11月至次年1月的冻胀速率较1~2月高18.7%；竖向位移极值与地温梯度呈负相关。研究成果为寒区道路抗冻设计与病害防控提供了理论依据。

Abstract

To investigate the deformation mechanisms of subgrade in cold regions under hydro-thermal coupling effects， this study focuses on the TQ-Yichun section （62.44 km） of the G1111 He-ha Expressway. A two-dimensional coupled hydro-thermal migration model was developed to analyze the interactions among temperature fields， moisture fields， and deformation fields during freeze-thaw cycles. Based on the theory of phase change heat transfer and non-isothermal moisture migration， Abaqus software was used for simulation and analysis， and nine integrated monitoring points were deployed along the roadway cross-section at a spacing of 1.25 m to quantify the deformation characteristics of the surface layer of the roadway in the freezing and thawing stages. Key findings include： The maximum frost heave displacement at the subgrade surface during freezing periods reached 8.48 mm， with cross-sectional deformation exhibiting a parabolic distribution. The thaw settlement rate in spring （March-May， 0.12 mm/day） was significantly higher than the frost heave rate during freezing periods （0.07 mm/day）， with frost heave rates between November-January surpassing those of January-February by 18.7%. Vertical deformation extremes showed a negative correlation with geothermal gradients. By introducing dynamic ice-water phase-change correction coefficients， the study improves simulation accuracy and proposes a predictive model for subgrade deformation. These results provide critical insights for optimizing frost-resistant designs in cold-region road engineering， particularly in determining interfacial heat transfer coefficients and insulation layer thickness. The research establishes a theoretical foundation for mitigating subgrade deformation and controlling freeze-thaw-induced pavement deterioration in permafrost environments.

Graphical abstract

关键词

水-热耦合 / 变形规律 / 数值模拟 / 冻土路基 / 有限元

Key words

hydro-thermal coupling / deformation law / numerical simulation / permafrost roadbed / finite element

0 引 言

随着我国基础设施建设的快速发展，寒区公路、铁路等交通线路的规模不断扩大^［1］。寒区路基作为道路工程的重要组成部分，其稳定性直接关系到道路的运营安全和使用寿命^［2］。在水分、温度等因素的作用下，寒区路基容易产生冻胀、融沉现象，导致路面出现坑槽、开裂等问题，严重影响行车安全^［3］。我国东北43°01′N-48°19′N，125°13′E-135°05′E地区为典型的季节性冻土区，冻胀融沉现象是造成该地区路基病害的主要原因^［4］。季节性冻土区冻融循环现象受大气温度周期性影响显著^［5］，路基内部热量分布和水分迁移与土质特性密切相关^［6，7］，需综合考虑多种因素对路基冻胀融沉变形的影响^［8］。本文对揭示寒区路基冻胀、融沉机理，提高路基抗冻性能具有重要意义。

针对季节性冻土中的冻胀融沉现象，国内外学者展开了一系列研究，探讨路基变形与土体冻融循环之间的相互关系。Lu等^［9］通过COMSOL二次开发模块，构建了温度-湿度多场耦合数值模型，结合经典冻融土柱实验验证，揭示了冻融过程中水分迁移规律及其与应力场的耦合机制；Su等^［10］创新性地以非饱和土壤水流动的达西定律为核心，构建以相对饱和度和温度为场函数的黑土水热耦合模型，揭示了路基边坡冻融相变过程中温度场、水分场及含冰量的动态分布规律，阐明了冰水相变驱动下水分迁移与热传导的交互机制；张哲等^［11］采用室内试验的方法，考察了冻融循环次数和应力水平对路基土静剪切强度和塑性变形的影响规律。国内外冻土路基冻胀融沉变形研究已经取得了一定的进展，特别是对冻融过程中的水热参数化研究^［12］，但是仍需进一步提升数值模型的准确性和实用性。

路基中水-热耦合效应客观存在，直接影响路基填土的稳定性及季节冻融层的水热变化，是导致路基出现一系列病害的直接诱因^［13］。在多物理场耦合模型的构建方面，国内外学者提出了冻土的水-热两场至水-热-力三场耦合模型，主要归纳为经验与半经验型、流体动力学模型、刚冰模型和热力学模型^［14］。考虑到水热耦合机制的复杂性，以及多物理场仿真软件可以更好地模拟复杂空间内流体流动和相变潜热等过程^［15］的特点，既有研究多采用COMSOL或Ansys进行水-热-力耦合分析^{［16，17］}。目前，通过建立数学模型和计算模型来预测冻胀融沉的发生、发展及其对路基稳定性影响的研究相对较少^［18］。

综上所述，相较于以往研究采用的通用仿真策略，本研究基于Abaqus平台构建寒区路基水-热-力耦合分析框架，考虑路基阴阳坡效应，借助Abaqus子程序（Film、Dflux、Usdfld等）灵活定义相变潜热、未冻水含量-温度关系及冰阻滞效应。通过融合Abaqus非线性求解优势与冻土多场耦合理论，本研究形成的数值模型可为寒区路基冻胀融沉的全生命周期预测及防控设计提供参考价值。

1 水-热耦合理论模型

1.1　温度场模型

外界气温变化导致路基土体经历冻结与融化循环，引发相变导热现象。本文利用显热容法模型处理相变问题^［19］，相应的相变平面非稳态导热模型由下式表示：

\rho \frac{\partial }{\partial t}\left[c_p^{\mathrm{e}\mathrm{f}\mathrm{f}}\left(T\right)T\right]=\frac{\partial }{\partial x}\left(k\left(T\right)\frac{\partial T}{\partial x}\right)（1）

式中：T（x，t）为温度场（℃），随位置x（m）和时间t（s）变化；\rho为路基材料密度（kg/m³），假设不随相变改变；c_p^{\mathrm{e}\mathrm{f}\mathrm{f}}\left(T\right)为等效比热容（J/（kg·℃）），包含显热和潜热；k（T）为导热系数（W/（m·℃）），随相态（冻/融）变化。

1.2　水分迁移模型

冻土路基的变形与水的冻结和融化息息相关，所以分析冻土路基的变形时必须考虑水场分布的影响。应用质量守恒定律与多孔介质流体得到流体连续性方程，结合达西定律，得到土体中水分运动的基本模型如下所示^［20］：

\frac{\partial \epsilon }{\partial t}=D\left(\frac{{\partial }^{\mathrm{2}}\epsilon }{\partial x^{\mathrm{2}}}+\frac{{\partial }^{\mathrm{2}}\epsilon }{\partial y^{\mathrm{2}}}\right)+D_{\mathrm{T}}\left(\frac{{\partial }^{\mathrm{2}}T}{\partial x^{\mathrm{2}}}+\frac{{\partial }^{\mathrm{2}}T}{\partial y^{\mathrm{2}}}\right)+

K-\frac{{\rho }_{\mathrm{I}}}{{\rho }_{\mathrm{w}}}\frac{\partial f_{\mathrm{s}}}{\partial t}（2）

式中：\epsilon为体积含水量；D（\epsilon）描述含水量梯度下的扩散能力（m²/s）；K（\epsilon）描述水力梯度下的流速（m/s）；\rho_I、\rho_w分别为冰和水的密度（kg/m³）；f_{\mathrm{s}}为固相率；D_{\mathrm{T}}为温度梯度下的水分扩散率。

温度梯度差异导致的水分流动和相变使得水分迁移模型复杂化。因此，引入温度梯度下的水分扩散率D_T和固相率f_s，其控制模型如下：

\frac{\partial \epsilon }{\partial t}=D\left(\epsilon \right)\left(\frac{{\partial }^{\mathrm{2}}\epsilon }{\partial x^{\mathrm{2}}}+\frac{{\partial }^{\mathrm{2}}\epsilon }{\partial y^{\mathrm{2}}}\right)+D_{\mathrm{T}}\left(\epsilon \right)\left(\frac{{\partial }^{\mathrm{2}}T}{\partial x^{\mathrm{2}}}+\frac{{\partial }^{\mathrm{2}}T}{\partial y^{\mathrm{2}}}\right)+

\frac{\partial }{\partial y}\left(K\left(\epsilon \right)\right)-\frac{{\rho }_{\mathrm{I}}}{{\rho }_{\mathrm{w}}}\frac{\partial f_{\mathrm{s}}}{\partial t}（3）

冰的阻滞作用主要体现在对液态水迁移路径的物理堵塞和对渗透率的动态抑制。液态水的渗透率k与未冻水含量正相关，且受固相率（f_s，即冰的体积占比）影响。通常将渗透率进行修正，如式（4）所示：

k=k_{\mathrm{0}}\cdot {\left(\frac{{\omega }_{\mathrm{u}}}{{\omega }_{\mathrm{s}\mathrm{a}\mathrm{t}}}\right)}^n\cdot {\left(\mathrm{1}-f_{\mathrm{s}}\right)}^m（4）

S=\frac{{\omega }_{\mathrm{u}}}{{\omega }_{\mathrm{s}\mathrm{a}\mathrm{t}}},\mathrm{0}\le S\le \mathrm{1}（5）

式中：k₀为饱和渗透率；{\omega }_{\mathrm{u}}为未冻水含量；{\omega }_{\mathrm{s}\mathrm{a}\mathrm{t}}为土体饱和含水率；S为饱和度；m为冰阻滞系数。f_s越高，渗透率越低。

水分扩散系数D_w也可通过固相率f_s​直接修正，如下所示：

D_{\mathrm{w}}=D_{\mathrm{w}\mathrm{0}}\cdot \left(\mathrm{1}-f_{\mathrm{s}}\right)（6）

冰的存在（f_s>0）会显著降低液态水的扩散能力。本研究在Abaqus水-热耦合分析中，通过用户子程序（USDFLD）显式计算f_s，并将其作为场变量传递给材料模型，用于修正渗透率k及扩散系数D_w。

1.3　模型简化

在模型简化中，设定特定时间点的导热系数、扩散系数及导水率为定值参数，忽略空间位置差异对参数的影响，同时基于土体各向同性假设，由此构建冻土路基温度场与水分场的耦合计算模型，如下所示^［21］：

\rho C_{\mathrm{p}}\frac{\partial T}{\partial t}=\lambda \left(\frac{{\partial }^{\mathrm{2}}T}{\partial x^{\mathrm{2}}}+\frac{{\partial }^{\mathrm{2}}T}{\partial y^{\mathrm{2}}}\right)+q_{\mathrm{v}}+\rho Q\frac{\partial f_{\mathrm{s}}}{\partial t}（7）

式中：T为物体的瞬态温度（℃）；t为过程进行的时间（s）；\rho为土体的密度（kg/m³）；C_p为土体的定压比热（J/（kg·℃））；q_v为内热源强度（W/m³）；Q为土壤冻融过程中的相变潜热（J/kg）。

模型（2）和（7）构建了冻土路基温度-水分耦合的理论框架。针对基于质能守恒定律的水热耦合问题构建数值求解框架，选用CPE4PTHT等支持温度-孔隙压力耦合的专用单元类型进行非结构化网格划分，确保梯度敏感区域的网格分辨率满足收敛需求。在材料定义阶段，通过定义子程序（USDFLD）精确表征流体渗透率、动态黏度及导热系数等关键参数随温度-压力场的非线性演化规律。数值求解过程采用耦合温度-孔隙压力分析模块，通过稳态/瞬态分析步设置，赋予初始温度场与孔隙压力场分布，并结合水力边界（Dirichlet压力边界、Neumann流量边界）与热边界（恒温条件、对流换热条件）构建封闭物理系统。对于复杂本构行为，则通过UMATHT热力学子程序实现自定义传热模型扩展。最终通过Abaqus后处理模块量化评估温度-压力时空分布特征，并基于系统总内能（ALLIE）与流体耗散（ALLFD）等守恒量指标验证数值解的物理合理性，其相对误差控制在1×10⁻⁶量级以内，证实了该建模策略在耦合传热-渗流问题中的计算鲁棒性与工程适用性。

2 数值模型构建

2.1　路基几何模型

本文着重分析路基横断面中二维水-热耦合场的分布特性，并假设其在长度方向上具有无限延展性。同时，公路路基采用分层压实工艺形成具备特定几何形态的路基断面，以满足路面结构的荷载传递需求。因各土层特性相异，需结合工程实际状况模拟路基各土层的材料组成。以鹤哈高速（G1111）苔青至伊春段（62.44 km）为对象进行建模与分析，该区域的最大冻结深度在2.0 m左右^［22］。根据中国《公路路基设计规范》（JTG D30-2015）规定，季冻区路基最小防冻厚度不小于当地最大冻结深度的1.2~1.5倍，由于该区域冻结深度达2 m，路基主体厚度需≥2.4 m。本文取路基主体厚度为2.5 m，路基主体填料为砂砾土；季节性活动层填土换填为亚黏土，厚度为3.0 m；路基顶面宽度为10 m；边坡1∶1.5。路基横断面及计算点位见图1。

2.2　网格划分

计算区域及网格划分如图2所示，有限元网格划分采用四结点平面应变四边形单元，共划分1 724个单元。计算中需考虑填土和砂土在长期自重荷载下的应力分布和变形特性，因此，其材料模型采用摩尔-库伦模型（Mohr-Coulomb model）。设置以下边界条件：地表边界设置为排水自由边界，允许水分自由渗透和排出；在季节活动层的底层边界，设定为不排水边界，阻止水分的垂直流动；对于模型的侧面边界，左侧和右侧分别施加固定横向位移的约束，以模拟实际的侧向限制；在模型底部，应用固定支撑条件，限制底部的水平和垂直位移，确保底部边界在模拟过程中保持稳定状态。

2.3　填土材料模型参数

路基填土结构如图1所示。砂砾土的干密度为1 870 kg/m³；亚黏土的干密度为2 025 kg/m³。通过热探针实验测定实际土样的热参数，数据见表1。

为确保路基稳定，需研究填筑土体水分迁移规律，通过数值模拟预测参数变化，制定水力耦合设计标准。土壤导水率与扩散率的精确测量是建立渗流变形数学模型的基础，为分析水分运动提供理论依据，数据见表2。

2.4　初始条件的确定

基于非稳态热传导理论框架，本研究采用区域大气月均值构建路基结构表层区热初始条件。针对原位热力学参数缺失现状，通过建立热边界将表层完全解冻相态的温度场与施工季节气候要素建立量化关联，形成具有工程适用性的温度边界设定方法^［23］。观测地区的气候数据如表3所示：

2.5　边界条件

路基温度场的边界条件复杂，现有研究多采用第一类边界条件（Dirichlet条件）。但该方法不仅需要对路基表层温度实施持续观测，还存在参数采集困难及忽略实际环境变量等问题，致使模型预测结果与实际工况存在偏差^［24］。针对季节性冻土区，需综合考虑太阳辐射、地面辐射和对流换热等因素，采用Neumann边界条件（第二类边界条件）与Robin边界条件（第三类边界条件）相结合的综合研究模式。基于此，本研究利用综合边界条件对温度场进行分析。

（1）阴阳坡边界设置

季冻区路基的阴阳坡因太阳辐射差异导致温度场与水分场分布不对称，显著影响冻胀融沉行为。路基阴阳坡为斜面，其太阳辐射吸收效率与平面有所差异，为此引入坡面系数ε^［24］。

阴坡：

\epsilon =\left(\mathrm{c}\mathrm{t}\mathrm{g}\alpha -\mathrm{c}\mathrm{t}\mathrm{g}\gamma \cdot \left|\mathrm{c}\mathrm{o}\mathrm{s}\left(\delta -\theta \right)\right|\right)\mathrm{s}\mathrm{i}\mathrm{n}\alpha（8）

阳坡：

\epsilon =\left(\mathrm{c}\mathrm{t}\mathrm{g}\alpha +\mathrm{c}\mathrm{t}\mathrm{g}\gamma \cdot \left|\mathrm{c}\mathrm{o}\mathrm{s}\left(\delta -\theta \right)\right|\right)\mathrm{s}\mathrm{i}\mathrm{n}\alpha（9）

式中：\alpha为边坡角度；\gamma为太阳高度角；\delta为路基走向角；\theta为太阳方位角。

在设置路基阴阳坡边界条件时，通过子程序（FILM）定义不同坡向的对流换热系数和环境温度。为设置阴阳坡的太阳辐射差异，通过子程序（DFLUX）引入坡面系数ε，量化坡度角对太阳辐射吸收效率的影响。

（2）确定上边界条件

本文采用第二类和第三类边界条件，考虑太阳辐射、对流换热和地面有效辐射的影响。太阳辐射主要影响大气温度变化，这里采用下式描述大气温度变化见下式：

T_{\mathrm{a}}=\overline{T_{\mathrm{a}}}+

T_{\mathrm{m}}\left[\mathrm{0.96}\mathrm{s}\mathrm{i}\mathrm{n}\omega \left(t-t_{\mathrm{0}}\right)+\mathrm{0.14}\mathrm{s}\mathrm{i}\mathrm{n}\mathrm{2}\omega \left(t-t_{\mathrm{0}}\right)\right]（10）

式中：T_a为大气温度；\overline{T_{\mathrm{a}}}为平均大气温度；T_m为温度振幅；\omega为角频率；t为时间；t₀为参考时间点。

对流换热则采用式（11）描述：

q_{\mathrm{h}}=h_{\mathrm{c}}\cdot \left(T_{\mathrm{a}}-T_{\mathrm{p}}\right)（11）

式中：q_h为对流换热量；h_c为对流换热系数；T_p为物体表面的温度。

（3）确定下边界条件

基于岩土热传导理论体系，研究区热稳定层深度经原位勘察确认为6 m，该层位土壤呈现显著的热惰性特征，其热通量波动幅度不超过±0.5 ℃/m。此外，为了计算方便，在数值模型中采用隔热边界设定下边界。

（4）确定风速边界条件

风速是控制地表热量和水分交换的重要因素，尤其在气候极端或环境敏感区域。在冻土地区，风速的变化对地表热力和水分状况的影响尤为显著。风速通过增强对流传热系数，改变路基表面热交换速率，进而影响其内部温度场。设置初始温度梯度为3 ℃，低风速条件下沥青混凝土路面的对流换热系数与风速间的回归公式如下所示^［25］：

h\left(v\right)=\mathrm{27.4}+\mathrm{29.8}v,R^{\mathrm{2}}=\mathrm{0.987}（12）

式中：h为对流传热系数（W/（m²·K））；v为风速（m/s）。

由国家气象局数据中心知，伊春地区公路路段的累计月平均风速为3.18 m/s，大风时风速可达24 m/s以上^［26］。通过Abaqus对风速边界条件进行非对称施加，需在Interaction模块中借用子程序（FILM）来定义对流换热系数h与风速v的复杂关系。根据伊春气象站（50774）所得地表风速数据，模拟路段风速如表4所示。

3 数据分析

3.1　温度场

分析温度场数据可知，冻土路基会随着季节的转换经历冻结与融化的周期性循环过程，这一过程受到温度波动的直接影响。本研究将这一过程划分为4个主要阶段：冻结阶段、冻结稳定阶段、融化阶段以及融化稳定阶段，并呈周期性重复。模拟计算结果显示，11月份气温降至冰点以下，路基开始进入冻结阶段；8月份后气温持续走低，至次年1月份触底，此时路基完全处于冻结稳定阶段；到了4月份，气温逐渐回暖至冰点以上，路基开始进入融化阶段；随着温度的进一步升高，直至7月份，季节活动层中的冻结土壤已经消融，路基处于融化稳定阶段。鉴于路基物理场数据量较大，分析全年各月份温度场分布时，将选取各阶段的典型状态进行探讨。

在研究区域，9月份路基表面的月平均气温达14.52 ℃。由图3可知，9月份路基基床层阳坡温度高于阴坡，填筑层温度场呈现出明显的层状分布。10月份，大气温度在0 ℃上下浮动，加之太阳辐射的减弱，路基中水分从外缘开始逐渐向内层凝固，内部仍保有未冻结的土壤。进入11月份，伊春市月平均气温降至-10 ℃左右，路基填筑层内温度变化明显滞后于基床层，表现出与外界降温趋势相异的升温状态，如图3所示（彩图参见电子版，以下同）。然而，由于路基温度场的时变覆盖效应，导致其透气透水性（覆盖条件）的改变，影响了土体内水分的迁移、集聚以及冻胀等过程，进而使路基内部保留了一定范围的未冻结土体区域。随着大气温度的持续降低，未冻结土体的区域不断缩小，并逐步向路基底部转移。

经观测，每年1月份，研究区域的平均大气温度降至年度最低水平，同时日均有效日照时间不足5 h，无论是阴坡还是阳坡区域的温度都降至一年中最低点。如图4所示，路基边缘区域的温度场表现出明显的层状分布特征。自2月份起，气温开始缓慢上升，同时太阳辐射强度也逐渐增强。4月份，气温上升至0 ℃以上，路基的各个边界开始吸收热量，标志着路基进入融化阶段。如图4所示，路基的融化过程是从边界处向下逐渐展开的。根据4月份的温度场显示，阳坡的融化速率快于阴坡，表现出不同坡面在融化过程中的温度变化差异。

图5为一年内路基顶面中点（即图1的A点）沿纵向垂直方向的温度变化趋势图。如图所示，一年中路基顶面以下2.5 m范围内每月纵向温度都有明显变化，从路基顶面以下2.5 m开始，每月纵向温度变化值逐渐减小，进入基底冻土层后，温度值变化微小，趋于稳定。

3.2　水分场

在季节性冻结区域的路面基层土壤经历冻融循环时，因土壤含水量的变化造成其热特性发生改变，同时土壤各层温度发生变化导致土水势的差异，从而影响路面基层的水分分布情况^［27］。饱和度S是表征土体孔隙中液态水填充程度的关键参数，直接影响渗透率与水分迁移能力，可量化冰阻滞作用对水分迁移的影响。本研究基于冻土路基内部温度与水分长期监测数据集，通过仿真量化了水分场的时空演变特征。

如图6所示，12月路基温度下降形成1.3 m深的冻结面，未冻水减少并聚集在冻结面附近。次年2月，冻结面下移至2.5 m，冻结过程中，冰的形成会占据孔隙空间，导致有效饱和度降低，进而抑制液态水的迁移。4月，气温回升，浅层冻结区融化，未冻水含量增加。6~8月，气温升高，土体温度在0 ℃以上，水分迁移少，未冻水含量稳定。7月，路基含水率从顶部至活动层底部逐渐增大，水分主要受重力影响扩散。

伊春地区位于黑龙江省小兴安岭腹地，属于高纬度岛状冻土区，岛状多年冻土在模拟对象区域有密度不同的分布，因对岛状多年冻土认识不清，调查深度不够，处理措施不当造成公路路基沉陷变形，影响行车质量的现象时有发生^［28］。考虑到岛状多年冻土的影响，本研究在Abaqus模型构建中将水分饱和带设置在地表埋深3 m（即路基深度5.5 m）以下。由图7可知，3~5月随着气温上升，季节活动层中冻结区开始融化，导致路基的含水率迅速上升。此时，路基的承载能力下降，容易出现弹簧现象（路基软化），若路基排水系统不完善，水分不能及时排出，饱和度会保持在较高水平，导致路基不稳定；6~8月大气气温升高，水分蒸发加强，路基饱和度逐渐降低，夏季降雨使路基饱和度短暂上升，但由于高温促进水分蒸发，饱和度通常会迅速回落；9~11月随着太阳辐射量降低以及秋季雨水增多，路基饱和度略有上升，但总体趋于稳定；12月至次年2月随着温度的持续下降，路基中的自由水开始冻结，体积膨胀，冻胀作用导致路基填土的孔隙减小，路基上部的水分因冻胀作用向上移动，导致上层路基水分饱和度增加。同时，冻结锋面以下的未冻水继续向冻结锋面移动，导致冻结锋面附近的水分饱和度较高。

由图6、图7可知，一年中路基浅层土体的未冻水含量随季节变化大，深层变化小，说明近路面含水率受气温影响明显，而深层含水率相对稳定。

3.3　位移场

本研究围绕水热耦合条件下路基表面的不均匀沉降现象，考虑寒区路基材料在低温下的非线性热力学和力学特性，通过Abaqus模拟各月份位移场，对模型位移场的变化规律及长期稳定性进行对比分析。

本研究计算的冻土路基变形包括：10月至次年3月的冻胀变形和4~9月的融沉变形。从图8可以看出，季节冻土区路基的变形过程主要发生在路基表面、路基阴阳坡两侧以及季节活动层等区域。

图9为路基从8月建成起到次年7月路基表面产生的竖向位移：在8、9月，由于路基内部未形成冻结土体，仅受路面车辆载荷的影响，冻土路基的变形较为轻微；10月，随着大气温度降低，路基基床层顶面填土中水分聚集凝结，开始形成零星的冻结土核，伴随冻结土核的扩散，路基顶面出现轻微冻胀变形，此时路基处于冻结阶段初期；从11~2月冻土路基表面产生的冻胀变形不断增大，相较于路基阳坡，路基右侧阴坡的冻胀变形较大。

3月，冻结阶段进入末期，路基中冻结土体的范围达到最大程度。4月，随着外界气温回升，路基中的冻结土体开始自上而下融化，路基表面出现略微沉降，此时路基处于融沉阶段初期。冻结状况下路基材料的热传导性较差，热量传递过程缓慢，导致路基内部温度变化滞后，融化过程所需的相变潜热不足，减缓了深层路基冻结土体的融化速率。3~4月正处于冻结阶段到融沉阶段的过渡期，受滞后性影响，4月路基表面的变形量相比于3月略微降低。5~7月，冻土路基表面的融沉变形持续增大，直至达到最大融沉深度。在路基自身重量及车辆载荷的双重影响下，路基中段的融沉变形最为显著。

图9直观地表达了不同月份路表的变形程度，但无法体现各个横向路段上的非均匀变形趋势。如图10所示：从8~10月份，阳坡由于接收更多的太阳辐射，其温度通常比阴坡高，融沉速率相较于阴坡更快，沉降值更大；从11月份到次年4月份，阴坡的热量交换比阳坡慢，这种热量交换的不平衡导致阴坡的冻土含水量较高，阴坡由于温度较低，水分迁移速度较慢，水分容易在阴坡积累，形成较高的水分饱和度。在冻结过程中，土壤中的水结冰，体积膨胀，导致阴坡更容易发生冻胀现象，从而引起更大的沉降；由表3可知，5~7月份，这3个月份的太阳辐射量全年最大，阴阳坡的温差减小，太阳辐射对于两个坡面的影响差异不再显著，导致两个坡面的温度趋于均衡，但阴坡的融化速率可能仍然低于阳坡。在此期间两个坡面的冻土都已经完成了大部分融化过程，沉降速率趋于一致。

路基表面的平均竖向变形通常为2~4 mm。在接近阴坡的区域，路面竖向变形较为显著，特别是路右侧边缘路段，竖向变形超过7 mm，最大变形量达8.48 mm，详细数据如表5所示。

路面车行道及中央分隔带的竖向位移差值大致为3 mm，呈抛物线分布，靠近路肩两侧路段的竖向位移较大，阴坡侧竖向相对变形最为显著。

根据图11可知，路面宽度区间\left(\mathrm{8.75},\left.\mathrm{10}\right]\right.靠近阴坡一侧路段出现冻胀极值的概率较高，占66.7%，由于阳坡的太阳辐射量大于阴坡，靠近阳坡一侧路段\left(\mathrm{0},\left.\mathrm{1.25}\right]\right.受冻胀影响相对较小，冻胀极值出现概率仅占25%，路肩的冻胀效应相较于路面中部区域更为明显，这表明阴阳坡的冻胀变形受太阳辐射影响较大。融沉阶段，受有效日照时间差异影响，阳坡接受了更多太阳辐射，其地表及地下温度通常高于阴坡，这导致阳坡一侧的沉降速率快于阴坡一侧，从而导致\left(\mathrm{2.5},\left.\mathrm{3.75}\right]\right.和\left(\mathrm{3.75},\left.\mathrm{5}\right]\right.区间内出现融沉变形极小值的概率较高，占75%。融沉变形极值大多集中于靠近中央分隔带两侧2.5 m范围内，左侧靠近阳坡路段出现极值的概率占78.3%，而右侧靠近阴坡路段仅占16.6%，这表明阴阳坡效应加剧了冻融循环的不均匀性，造成路基两侧明显的变形差异。

3.4　路基表面相对变形规律分析

基于线性回归与二次函数模型对路基表面多时相的向位移监测数据进行参数拟合，通过横截面竖向变形特征的解析方程，对冻融循环作用下路基表层的时空演变规律进行定量表征^［29］，如表6所示。

表6中：y为竖向相对变形量，mm；d为横截面上结点距离阳坡路肩边缘的距离，m。2月、3月、6月、8月、9月份时横截面采用抛物线模型拟合精度较高，决定系数R²均超过0.85阈值，说明各结点变形分布与抛物线形态具有显著相关性。

3.5　模型验证

通过埋设温度传感器获取季冻区路基0.1 m埋深的温度时序监测数据，与Abaqus模拟的温度场分布进行对比，结果显示温度场的时空演化特征与实测数据吻合良好，平均误差控制在±0.5 ℃以内，如图12所示。

同时，结合鹤哈高速（G1111）苔青至伊春段（62.44 km）工程现场人工测量的路基沉降量，对比模拟输出的路基变形曲线及空间分布特征，发现路基中心沉降的模拟值与实测值相对误差小于10%，且长期冻胀累积趋势与位移场预测模型一致，如图13所示。

上述验证结果表明，水-热耦合模型能够准确表征季节性冻土区路基温度场动态响应与冻胀变形的耦合机制，为实际工程冻害防控提供依据。

4 结 论

（1）通过对路基温度场的数值模拟分析，发现路基顶面深度6 m以下土体的温度持续处于负值区间，并且呈现出相对稳定的状态，这反映出冻土层的稳定性相对较强。然而，路基基床层及填筑层的温度明显较周边地区高，这一现象揭示了路基工程对冻土稳定性的不利影响。鉴于此，在施工过程中，必须实施相应的保护措施来维护冻土的稳定性。

（2）根据路基水分场数值模拟结果分析得出，在冻胀融沉周期内，路基浅层土体的未冻水含量变化显著，而深层土体变化微乎其微。这表明，靠近路面的路基水分含量受外界温度波动影响较大，与之呈现高度相关性，而深层土体的水分含量则相对稳定，不易受外界气温影响。

（3）对冻胀融沉下的路基位移场进行数据分析，结果表明：冻结阶段路基表面部分冻胀变形产生的竖向位移最大可达8.48 mm。融化阶段路面部分融沉变形产生的沉降值最大可达1.76 mm，融沉变形主要发生于路基表面中部，春季3~5月融沉最快。

（4）本研究通过Abaqus多物理场耦合能力，建立了考虑温度梯度影响的冻土水热力耦合模型，并基于Abaqus子程序模型开发了相应的有限元分析程序。这些模型和分析程序经过实例分析验证，显示出较好的可行性和合理性。该模型的局限性在于简化了冻土饱和度模型，且参数标定需大量实验数据支撑，实际地理与气象条件的复杂性未被完全纳入，导致局部预测精度受限。总体而言，该模型为冻土工程提供了预测工具，但需结合工程经验与数据优化以提升适用性。

参考文献

[1] 陈佩佩, 尚智, 王熳祺, 等. 寒区高速铁路冻土路基的水热耦合响应特性分析[J]. 水文地质工程地质, 2024, 52(3): 152-162.

[2] Chen Pei-pei, Shang Zhi, Wang Firm-qi, et al. Characteristic analysis of hydrothermal coupling response of frozen soil roadbed for high-speed railway in cold region[J]. Hydrogeology Engineering Geology, 2024, 52(3): 152-162.

[3] 彭勇. 路堑边坡防护设计方案[J]. 中国科技信息, 2024(14): 93-95.

[4] Peng Yong. Design scheme of road graben slope protection[J]. China Science and Technology Information, 2024(14): 93-95.

[5] 于景飞, 张文月, 李明浩. 季节性冻土区铁路路基冻融变形规律[J]. 长江科学院院报, 2020, 37(5):133-138.

[6] Yu Jing-fei, Zhang Wen-yue, Li Ming-hao. Freezing and thawing deformation law of railway roadbed in seasonal permafrost area[J]. Journal of Changjiang Academy of Sciences, 2020, 37(5): 133-138.

[7] 李晓芳, 官少龙, 滕伟玲. 基于高山草甸区季冻土的冻融变形特性影响因素研究[J]. 公路,2018, 63(4): 47-51.

[8] Li Xiao-fang, Guan Shao-long, Teng Wei-ling. Research on factors affecting freeze-thaw deformation characteristics based on seasonal permafrost in alpine meadow area[J]. Highway, 2018, 63(4): 47-51.

[9] 李万志, 马海玲, 庞昕玮, 等. 气候变暖背景下青海三江源区季节冻土冻融特征研究[J]. 冰川冻土, 2023, 45(4): 1233-1241.

[10] Li Wan-zhi, Ma Hai-ling, Pang Xin-wei, et al. Study on permafrost and thawing characteristics of seasonally frozen soil in the San-jiang-yuan area of Qinghai under the background of climate warming[J]. Glacial Permafrost, 2023, 45(4): 1233-1241.

[11] 赫金良, 李晓宁, 李青龙, 等. 基于优化SVM的粗粒土路基冻胀特性影响因素敏感性分析[J]. 科学技术与工程, 2024, 24(22): 9577-9586.

[12] He Jin-liang, Li Xiao-ning, Li Qing-long, et al. Sensitivity analysis of factors influencing the freezing and expansion characteristics of coarse-grained soil subgrade based on optimised SVM[J]. Science, Technology and Engineering, 2024, 24(22): 9577-9586.

[13] 张玉芝, 杨威, 李锐, 等. 季冻区重载铁路路基冻胀融沉监测及规律分析[J]. 中国铁道科学, 2022, 43(5): 1-10.

[14] Zhang Yu-zhi, Yang Wei, Li Rui, et al. Monitoring and law analysis of freezing and thawing subsidence of heavy railway roadbase in monsoon frozen area[J]. China Railway Science, 2022, 43(5): 1-10.

[15] 甘海龙, 唐先习, 李小博. 季冻区高速铁路路基冻融变形特征[J]. 科学技术与工程, 2022, 22(9): 3760-3767.

[16] Gan Hai-long, Tang Xian-xi, Li Xiao-bo. Freezing and thawing deformation characteristics of high-speed railway roadbase in monsoon frozen area[J]. Science Technology and Engineering, 2022, 22(9): 3760-3767.

[17] Lu B, Zhao W, Li S, et al. Study on seasonal permafrost roadbed deformation based on water-heat coupling characteristics[J]. Buildings, 2024, 14(9): No.2710.

[18] Su A, Hai M, Wang M, et al. Analytical study on water and heat coupling process of black soil roadbed slope in seasonal frozen soil region[J]. Sustainability, 2024, 16(19): No.8427.

[19] 张哲, 付伟, 张军辉, 等. 循环荷载下冻融路基黏土长期塑性行为[J]. 吉林大学学报: 工学版, 2023, 53(6): 1790-1798.

[20] Zhang Zhe, Fu Wei, Zhang Jun-hui, et al. Long-term plastic behavior of frozen-thawed roadbed clay under cyclic loading[J]. Journal of Jilin University(Engineering and Technology Edition), 2023, 53(6): 1790-1798.

[21] 侯雅, 李伟平, 左金清. 土壤冻融过程中的水热参数化方案研究进展[J]. 高原气象, 2024, 44(1): 1-15.

[22] Hou Ya, Li Wei-ping, Zuo Jin-qing. Progress of hydrothermal parameterisation scheme during soil freezing and thawing[J]. Plateau Meteorology, 2024, 44(1): 1-15.

[23] 赵亮, 景立平, 单振东. 冻土冻胀模型研究现状与进展[J]. 自然灾害学报, 2020, 29(4): 43-52.

[24] Zhao Liang, Jing Li-ping, Shan Zhen-dong. Current status and progress of research on permafrost frost expansion modelling[J]. Journal of Natural Hazards, 2020, 29(4): 43-52.

[25] 宋存牛. 冻融过程中土体水热力耦合作用理论和模型研究进展[J]. 冰川冻土, 2010, 32(5): 982-988.

[26] Song Cun-niu. Progress in theory and modelling of soil hydrothermal coupling action during freeze-thaw process[J]. Glacial Permafrost, 2010, 32(5): 982-988.

[27] Zhan Y, Lu Z, Yao H. Numerical analysis of thermo-hydro-mechanical coupling of diversion tunnels in a seasonally frozen region[J]. Journal of Cold Regions Engineering, 2020, 34(3): No.04020018.

[28] 王自力, 李金峰, 滕继东, 等. 基于未冻水膜压力判据的冻土水热力耦合冻胀模型研究[J]. 岩土工程学报, 2023, 45(5): 997-1007.

[29] Wang Zi-li, Li Jin-feng, Teng Ji-dong, et al. Study on hydrothermal coupled freezing model of frozen soil based on unfrozen water film pressure criterion[J]. Journal of Geotechnical Engineering, 2023, 45(5): 997-1007.

[30] 李长雨, 马桂霞, 郝光, 等. 季节冻土地区路基冷阻层温度场效应[J]. 吉林大学学报: 地球科学版, 2018, 48(4): 1174-1181.

[31] Li Chang-yu, Ma Gui-xia, Hao Guang, et al. Temperature field effect of cold resistance layer of roadbed in seasonal frozen soil area[J]. Journal of Jilin University(Earth Science Edition), 2018, 48(4): 1174-1181.

[32] 苗祺, 牛富俊, 林战举, 等. 季节冻土区高铁路基冻胀研究进展及展望[J]. 冰川冻土, 2019, 41(3):669-679.

[33] Miao Qi, Niu Fu-jun, Lin Zhan-ju, et al. Progress and prospect of research on freezing and expansion of high-speed railway roadbed in seasonal frozen soil zone[J]. Glacial Permafrost, 2019, 41(3): 669-679.

[34] 毛雪松, 李宁, 王秉纲, 等. 多年冻土路基水-热-力耦合理论模型及数值模拟[J]. 长安大学学报: 自然科学版, 2006(4): 16-19, 62.

[35] Mao Xue-song, Li Ning, Wang Bing-gang, et al. Theoretical model and numerical simulation of water-heat-force coupling of multi-year frozen soil roadbed[J]. Journal of Chang'an University (Natural Science Edition), 2006(4): 16-19, 62.

[36] 赵刚. 寒区高等级公路路堑边坡春季浅层滑塌机理研究[D]. 哈尔滨: 哈尔滨工业大学交通科学与工程学院, 2010.

[37] Zhao Gang. Study on the mechanism of shallow landslide of high-grade highway cutting slope in cold region in spring[D]. Harbin: School of Transportation Science and Engineering, Harbin Institute of Technology, 2010.

[38] 毛雪松, 王秉纲, 胡长顺, 等. 冻土路基水热迁移问题的理论模型及数值模拟[J]. 中外公路, 2006(1): 23-26.

[39] Mao Xue-song, Wang Bing-gang, Hu Chang-shun, et al. Theoretical modelling and numerical simulation of water-heat migration problem in permafrost roadbeds[J]. Chinese and Foreign Highway, 2006(1): 23-26.

[40] 龚强, 晁华, 朱玲, 等. 东北地区地温和冻结深度时空特征的细化分析[J]. 冰川冻土, 2021, 43(6): 1782-1793.

[41] Gong Qiang, Chao Hua, Zhu Ling, et al. Refinement analysis of spatial and temporal characteristics of ground temperature and freezing depth in Northeast China[J]. Glacial Permafrost, 2021, 43(6): 1782-1793.

[42] 杨科. 多年冻土路基温度场变化及数值模拟[J]. 中国科技信息, 2022(14): 81-83.

[43] Yang Ke. Variation and numerical simulation of permafrost subgrade temperature field[J]. China Science and Technology Information, 2022(14): 81-83.

[44] 柳波. 多年冻土地区公路路基不稳定变形研究[D]. 西安: 长安大学公路学院, 2015.

[45] Liu Bo. Study on unstable deformation of highway subgrade in permafrost area[D]. Xi'an: School of Highway, Chang-an University, 2015.

[46] 汪海年, 侯仲杰, 杨玉杰. 多年冻土路基对流换热试验方法[J]. 公路交通科技, 2009, 26(3): 43-48.

[47] Wang Hai-nian, Hou Zhong-jie, Yang Yu-jie. Test method of convective heat transfer in roadbeds with multiyear frozen soil[J]. Highway and Transportation Science and Technology, 2009, 26(3): 43-48.

[48] 李方昊, 范昊明, 石昊, 等. 东北地区极端大风时空分布特征及其对风蚀量的影响[J]. 应用生态学报, 2024, 35(1): 87-94.

[49] Li Fang-hao, Fan Hao-ming, Shi Hao, et al. Temporal and spatial distribution characteristics of extreme high winds in Northeast China and their effects on wind erosion[J]. Journal of Applied Ecology, 2024, 35(1): 87-94.

[50] 毛云程, 李国玉, 张青龙, 等. 季节冻土区黄土路基水分与温度变化规律研究[J]. 冰川冻土, 2014, 36(4): 1011-1016.

[51] Mao Yun-cheng, Li Guo-yu, Zhang Qing-long, et al. Study on the variation of moisture and temperature of loess subgrade in seasonal permafrost area[J]. Journal of Glaciology and Geocryology, 2014, 36(4): 1011-1016.

[52] 徐泉涌. 某公路岛状多年冻土路基处理[J]. 低温建筑技术, 2005(3): 91-92.

[53] Xu Quan-yong. Treatment of island-shaped perennial frozen soil roadbed for a highway[J]. Low Temperature Construction Technology, 2005(3): 91-92.

[54] 纪伦, 张红菊, 郭宏斌, 等. 北方寒冷地区温度场作用下公路路表形变规律分析[J]. 哈尔滨工业大学学报, 2024, 56(9): 113-123.

[55] Ji Lun, Zhang Hong-ju, Guo Hong-bin, et al. Analysis of surface deformation of highway road surface under the action of temperature field in cold areas of northern China[J]. Journal of Harbin Institute of Technology, 2024, 56(9): 113-123.

基金资助

长沙理工大学极端环境绿色长寿道路工程全国重点实验室开放基金项目(kfj230105)

黑龙江省自然科学基金项目(LH2023E011)