THM coupled model for simulating frost heave based on a new water film pressure criterion
-
摘要: 冻胀融沉是多年冻土区的主要病害，其受水分场、温度场、应力场的复杂耦合作用影响。基于水膜理论提出了冻土未冻水膜压力作为冰透镜体生成的判据，并重新对水分迁移驱动作用进行描述，建立了以温度、土体孔隙比为变量的全耦合模型。通过考虑已冻区冰基质的影响，推导了涵盖原位冻胀与冰分凝两部分的冻胀量计算公式。基于Matlab和COMSOL Multiphysics的联立平台，提出了模型冰透镜体实时分布的数值求解方法，实现了冻土温度、水分、应力、冰透镜体分布的全耦合数值求解。通过与室内土柱冻结试验及现有水热力模型（Thermal-Hydraulic-Mechanical，又称THM模型）冻胀量结果进行对比分析，证明了温度、含水率与冻胀计算上的可靠稳定。最后通过探讨温度梯度、上覆压力、渗透系数与压缩模量对土柱冻结的影响发现，温度梯度能显著增加土体冻胀量，上覆压力会导致更多水分向冻结锋面迁移，但对冻胀量起着抑制作用，渗透系数与压缩模量均对冻胀量产生正影响。为冻胀理论研究与数值实现提供了新思路。Abstract: The frost heave and thaw settlement are the main frost damage in cold areas, which are the complex coupling process of water, temperature and stress fields. In this study, a coupled thermal-hydraulic-mechanical model is developed based on the water film theory, in which the temperature and void ratio of soils are the input variables. The novelty of this model is that the frozen water film pressure is used as the criterion for the generation of ice lens. The driving force of water migration is newly defined, and the frost heave includes the pristine frost heave and the amount of ice segregation. The fully coupled model is numerically solved based on the Matlab and COMSOL Multiphysics, generating the results of soil temperature, moisture, stress and the layered ice lens. The simulated results are then compared with those of the laboratory freezing tests, which shows that they match quite well and verify the validity of the proposed model. The simulation indicates that temperature gradient can promote the frost heave, and the overburden pressure can attract more water to the freezing front but decrease the amount of the frost heave. In addition, both the hydraulic conductivity and the compressive modulus have positive effects on the frost heave. The proposed model provides a new approach to understand the frost heave.
-
0. 引言
中国现阶段正处于水下盾构隧道建设的高峰期，水下隧道的建设总体呈现出大断面、大埋深、高水压和地质条件日益复杂的趋势。不良地质段的沉降控制是海底隧道施工的核心技术之一，它不仅对施工安全至关重要，而且对隧道运营阶段的结构受力及长期安全性也是极其重要的[1]。断层破碎带是海底盾构隧道建设过程中常见的不良地质段，在断层破碎带地层界面附近，地层特性和隧道上覆荷载均发生显著变化，极易造成隧道的不均匀沉降，从而导致环缝张开、错台、螺栓应力增大乃至接缝渗漏水等现象，严重影响盾构隧道的结构稳定和运营安全[2]。
目前常用的隧道上覆荷载计算方法有全土柱法、太沙基理论[3]、普氏理论[4]和规范公式[5]等。黎春林[6]在Terzaghi[7]基础上推导了一种可以考虑地层损失和管片刚度的松动土压力计算公式。王浩等[8]基于普氏平衡拱理论，提出适用于3洞小净距隧道的围岩压力计算方法。颉永斌等[9]用筒仓理论结合地层应力分布特征来求解断层破碎带内隧道的上覆荷载。由于断层破碎带形态各异，填充物复杂，目前断层破碎带地层隧道上覆荷载仍难以全面准确描述。
不同学者[10-12]将既有盾构隧道简化为无限长弹性直梁或Euler-Bernoulli梁放置于Winkler地基模型或Pasternak地基模型上，研究既有隧道的力学响应。然而Euler-Bernoulli梁忽略了实际工程常发生的剪切变形。为更合理地描述隧道的实际变形模式，Wu等[13]针对隧道的剪切变形，基于Liao等[14]接头有一定影响范围思路同样引入剪切刚度影响系数推导隧道管片等效抗剪刚度，将隧道视为Timoshenko梁，从而考虑隧道纵向剪切效应。Zhang等[15]在Wu等[13]基础上将既有隧道简化为能够考虑地层的压缩及剪切效应的kerr三参数地基上Timoshenko梁，并依据高斯曲线推导出新隧道对既有隧道的纵向响应的解析解。目前Timoshenko梁在多段非均布复杂上覆荷载作用下难以直接求解，通常以有限差分法进行数值计算。
为此，本文将隧道上覆荷载概括为4种特征形态断层破碎带松动土压力的线性组合，提出了归一化上覆荷载理论公式，并将复杂上覆荷载在5%的误差内简化为三次泰勒多项式，结合Timoshenko梁理论推导建立了断层破碎带海底盾构隧道管片纵向沉降及内力系列解析解，并通过现场实测结果验证了本文理论模型和解析解的正确性。
1. 断层破碎带隧道纵向上覆荷载分析
1.1 断层破碎带形态分析
断层破碎带是海底隧道建设中经常遇到的不良地质条件，由于断层破碎带内岩体构造发育，渗透性较大，使得海底盾构隧道的上覆荷载比较复杂。如图 1所示，将隧道沿纵向分为4段：ab段为正常段，bc段为相交段，cd段为延伸段，de段为正常段。不同的相交段和延伸段组合可模拟不同形态的断层破碎带。
断层破碎带的规模根据地质构造的不同，从几米到几千米不等。对于任意断层破碎带形态都可以看成以下4种特征形态的线性组合，如图 2所示。其中相交段为形态一、三的线性组合，延伸段为形态二、四的线性组合。
![]()
图 2 断层破碎带地层典型特征形态Figure 2. Typical characteristics and morphologies of strata in fault fracture zone1.2 隧道纵向上覆荷载计算
对断层破碎带范围内海底隧道任意横剖面分析，隧道的横向上覆荷载均可视为多层土体作用下隧道的松动土压力，力学模型如图 3所示。土体上作用水压力p0，忽略水的渗流作用。从上至下3种土层厚度分别为ha1，ha2，ha3。
取第一层正常围岩松动区内深度y处宽为2Bt，厚度为dy的土条单元，根据土条竖向力的平衡，可得
2Btσy+2Btγsdy=2Bt(σy+dσy)+2τsdy。 (1) 式中：γs为正常围岩的重度；σy为深度y处的垂直方向应力；τs为作用于正常围岩地层土条上向上的摩擦应力；Bt为隧道开挖影响宽度的一半，t=1，2，B1为相交段范围内隧道开挖影响宽度的一半，B2为其余段隧道开挖影响宽度的一半，τs和Bt由下式计算：
τs=cs+σxtanφs=cs+ksσytanφs， (2) Bt={rcot(π /4 + φp/22) (t=1)rcot(π /4 + φs/22) (t=2)。 (3) 式中：σx为深度y+dy/2处的水平方向应力；cs为正常围岩黏聚力；φs为正常围岩内摩擦角；φp为断层破碎带内摩擦角；r为隧道的外半径；ks为正常围岩静止侧压力系数，按经验公式k0=1−sinφs进行计算[16]。
将式(2)，(3)代入式(1)，结合边界条件y=0，σy=p0，并代入y=ha1可以得到第一层土体底部荷载：
q(ha1)=Btγs−cskstanφs(1−e−kstanφs⋅ha1/Bt)+p0e−kstanφs⋅ha1/Bt。 (4) 同理，将上层土体底部荷载当作边界条件代入下层土体中进行计算，即可求得3层土体作用下隧道松动土压力：
q(x)=Btγs−cskstanφs(1−e−kstanφs⋅ha3/Bt)+(Btγp−cpkptanφp(1−e−kptanφp⋅ha2/Bt)+(Btγp−cpkptanφp(1−e−kptanφp⋅ha2/Bt))+p0e−kstanφs⋅ha1/Bt)⋅e−kptanφp⋅ha2/Bt)e−kstanφs⋅ha3/Bt。 (5) 式中：cp为断层破碎带的黏聚力；γp为断层破碎带的重度；kp为断层破碎带静止侧压力系数，按经验公式k0=1−sinφp进行计算[16]。
式(5)即为断层破碎带地层海底盾构隧道归一化上覆荷载理论公式。4种特征形态下海底盾构隧道纵向上覆荷载，如表 1所示。
表 1 4种特征形态下海底盾构隧道纵向上覆荷载Table 1. Longitudinal overlying loads of subsea shield tunnel for four characteristic forms形态 纵向上覆荷载 Bt取值 形态一 q1(x)=q(x)|ha1 = h1, ha1 = h1, ha3 = 0 Bt=B1 形态二 q2(x)=q(x)|ha1 = 0, ha2 = h3, ha3 = h4 Bt=B2 形态三 q3(x)=q(x)|ha1=ha3=0,ha2=h5 Bt=B1 形态四 q4(x)=q(x)|ha1=h6,ha2=h7,ha3=h8 Bt=B2 对于正常段ab、de，海底盾构隧道的纵向上覆荷载q0(x)：
q0(x) = B2γs−cskstanφs(1−e−kstanφs⋅ha1/B2) + p0e−kstanφs⋅H/B2。 (6) 结合正常段、相交段、延伸段的纵向上覆荷载，对于不同形态断层破碎带，纵向上覆荷载如表 2所示，其中H为隧道埋深，bp为断层破碎带宽度，α为断层破碎带与隧道轴线的夹角。
表 2 不同倾角断层破碎带海底隧道纵向上覆荷载Table 2. Longitudinal overlying loads of subsea tunnel in fault fracture zone with different inclination angles倾角 上覆荷载 α=90∘ q0(x)−∞ < x < 0 q3(x)0⩽x < bp q0(x)bp⩽x < + ∞ tanα>Hbp q0(x)−∞ < x < 0 q1(x)0⩽x < Htanα q3(x)Htanα⩽x⩽bp q2(x)bp < x < bp + Htanα q0(x)bp + Htanα⩽x<+∞ tanα=Hbp q0(x)−∞ < x < 0 q1(x)0⩽x < bp q3(x)x=bp q2(x)bp < x⩽2bp q0(x)2bp < x < + ∞ tanα<Hbp q0(x)−∞ < x < 0 q1(x)0⩽x < bp q4(x)bp⩽x<Htanα q2(x)Htanα⩽x⩽bp + Htanα q0(x)bp + Htanα < x < +∞ 2. 盾构隧道纵向沉降解析解
2.1 隧道纵向沉降理论模型
为了更加准确计算盾构隧道的纵向沉降与内力, 考虑管片接头对隧道整体剪切刚度的影响，将隧道管片结构视为置于Winkler地基上的Timoshenko梁，其纵向力学理论模型如图 4所示。其中fj(x)为盾构隧道纵向上覆荷载，由表 3可知，隧道纵向上覆荷载fj(x)共分为四类，即j=1~4。k2为断层破碎带地层地基基床系数；k1为正常地层地基基床系数。根据隧道纵向上覆荷载的分段及地层的变化，盾构隧道沿纵向可以分为多段Timoshenko梁。
![]()
图 4 盾构隧道管片结构纵向理论模型Figure 4. Longitudinal theoretical model for segment structure of shield tunnel表 3 4种破碎带形态及其对应上覆荷载Table 3. Four types of fracture zones and their corresponding overlying loads破碎带形态 α=90∘ tanα>Hbp tanα=Hbp tanα<Hbp 上覆荷载 f1(x) f2(x) f3(x) f4(x) 根据Timoshenko梁理论, 各段梁弯矩Mij、剪力Qij、转角θij满足下面的关系：
Mij=−(EI)eqdθij(x)dx， (7) Qij=(kGA)eq(dwij(x)dx−θij(x))， (8) θij=dwij(x)dx + 1(kGA)eq∫[2fj(x)r−2kirwij(x)]dx。 (9) 式中：(EI)eq为隧道纵向等效抗弯刚度；(kGA)eq为隧道纵向等效抗剪刚度；i为梁段号，ki为第i段梁的地基基床系数；Mij为第j类上覆荷载作用在隧道上i段梁产生的弯矩。4类上覆荷载对应破碎带形态如表 3所示。
在纵向上覆荷载fj(x)作用下，可以得到在Winkler地基上Timoshenko梁关于纵向沉降wij(x)的微分方程：
d4wij(x)dx4−2kir(kGA)eqd2wij(x)dx2 + 2kir(EI)eqwij(x)=2rfj(x)(EI)eq−2r(kGA)eqd2fj(x)dx2。 (10) 式(10)对应的齐次方程为
d4wij(x)dx4−Aiad2wij(x)dx2 + Aibwij(x)=0， (11) 式中，Aia=−2kir(kGA)eq，Aib=2kir(EI)eq。
其齐次方程通解为
wij(x)=eαix[ci1cos(βix) + ci2sin(βix)] + e−αix[ci3cos(βix) + ci4sin(βix)]， (12) 式中，αi = √2√Aib−Aia2，βi=4Aib−A2ib√2√Aib−Aia，ci1，ci2，ci3，ci4为待定系数。
2.2 隧道纵向上覆荷载简化
由表 3可知，fj(x)在其作用范围内任意开区间(x0，xn)具有直到n+1阶的导数。xt为开区间(x0，xn)内的一点，对于任意x∈(x0，xn)，由泰勒公式得：
fj(x)=fj(xt)0！+f′j(xt)1！(x−xt)+f″j(xt)2！(x−xt)2+⋯+fnj(xt)n！(x−xt)n+Rn(x)， (13) Rn(x)=f(n+1)j(ξ)(n+1)!(x−xt)n+1。 (14) 式中，ξ为xt与x之间的某个值。当x→xt时，误差|Rn(x)|是比(x−xt)n+1高阶的无穷小，x越趋近于xt，误差就越小。利用泰勒公式这个性质，可以提出高精度泰勒公式分段展开法，即对不同的x进行不同xt处高精度n次多项式展开。
本文将纵向上覆荷载fj(x)展开为三次泰勒多项式，其与fj(x)的误差为
R3(x)=f(4)j(ξ)4!(x−xt)4。 (15) 对于有界函数fj(x)，当x∈(x0，xn)时，∣f(4)j(x)∣ ⩽B，则有下列不等式：
|R3(x)|⩽B4!(x−xt)4。 (16) 若本文最大误差取为5%，可以得到
|R3(x)|⩽0.05fj(x)。 (17) 联立求解式(16)，(17)，可以求出fj(x)在xt处最大误差为5%的展开区间(x0，xn)。即fj(x)可以在5%误差内简化为三次多项式函数：
fj(x)≈fj(xt)0！+f′j(xt)1！(x−xt)+f″j(xt)2！(x−xt)2+f‴j(xt)3！(x−xt)3，x∈(x0,xn)。 (18) 此时纵向沉降微分方程的特解为
w∗ij(x) = fj(x)ki。 (19) 即Timoshenko梁纵向沉降微分方程通解wij(x)为
wij(x)=eαix[ci1cos(βix) + ci2sin(βix)] + e−αix[ci3cos(βix) + ci4sin(βix)] + fj(x)ki。 (20) 假设隧道两端自由，可得边界条件：
limx→ + ∞Qij(x)=limx→ - ∞Qij(x)=0， (21) limx→ + ∞Mij(x)=limx→ - ∞Mij(x)=0。 (22) 假定隧道沉降曲线连续可微可导，在各段梁节点处隧道的弯矩、剪力、转角、纵向变形必然相等，可以列出以下方程：
limx→xpwij(x)=limx→xpwij(x)， (23) limx→xpθij(x) = limx→xpθij(x)， (24) limx→xpMij(x) = limx→xpMij(x)， (25) limx→xpQij(x) = limx→xpQij(x)， (26) 式中，xp为各连接点横坐标。
由式(21)~(26)及式(7)~(9)，(12)可以求出通解中待定系数ci1，ci2，ci3，ci4。
2.3 相关参数的确定
隧道纵向等效抗弯刚度(EI)eq及纵向等效剪切刚度(kGA)eq是决定隧道变形及内力分布的重要参数，其等效纵向抗弯刚度[17]表达式为
{{\text{(}}EI{\text{)}}_{{\text{eq}}}} = \frac{{{\text{co}}{{\text{s}}^{\text{3}}}\varphi }}{{{\text{cos}}\phi {\text{ + }}\left( {\phi {\text{ + }}\frac{{\text{π }}}{{\text{2}}}} \right){\text{sin}{\rm{ \mathsf{ φ} }} }}}{E_{\text{c}}}{I_{\text{c}}}。 (27) 式中：Ec为管片结构弹性模量；Ic为管片环纵向惯性矩；\phi 为中性轴位置的角度，由下式确定：
\phi {\text{ + cot}}\phi = {\text{π }}\left( {{\text{0}}{\text{.5 + }}\frac{{n{k_{\text{b}}}{l_{\text{s}}}}}{{{E_{\text{c}}}{A_{\text{c}}}}}} \right)。 (28) 式中：kb接头螺栓的平均线刚度，{k_{\text{b}}}{\text{ = }}{E_{\text{b}}}{A_{\text{b}}}{\text{/}}{l_{\text{b}}}，Eb为螺栓的弹性模量，Ab为螺栓的横截面积, lb为螺栓长度；n为螺栓个数；ls为管片环宽；Ac为隧道管片横截面积。
盾构隧道等效剪切刚度(kGA)eq公式[13]为
{(kGA)_{{\text{eq}}}} = \zeta \frac{{{l_{\text{s}}}}}{{\frac{{{l_{\text{b}}}}}{{n{\kappa _{\text{b}}}{G_{\text{b}}}{A_{\text{b}}}}} + \frac{{{l_{\text{s}}} - {l_{\text{b}}}}}{{{\kappa _{\text{c}}}{G_{\text{c}}}{A_{\text{c}}}}}}}。 (29) 式中：\zeta 为修正系数，在本文其值取1；{\kappa _{\text{b}}}及{\kappa _{\text{c}}}分别为螺栓及管片环Timoshenko剪切系数，对于圆形截面螺栓，{\kappa _{\text{b}}}取0.9，对于环形隧道管片环结构，{\kappa _{\text{c}}}取0.5；Gb及Gc为分别螺栓剪切刚度及隧道管片剪切刚度，两者与各自弹性模量的关系为{G_{\text{b}}} = \frac{{{E_{\text{b}}}}}{{2(1 + {\nu _{\text{b}}})}}及{G_{\text{c}}} = \frac{{{E_{\text{c}}}}}{{{\text{2(1 + }}{\nu _{\text{c}}}{\text{)}}}}，其中，{\nu _{\text{b}}}及{\nu _{\text{c}}}分别为螺栓及管片的泊松比。
Winkler地基基床系数k表达式[18]为
{k_{\text{s}}} = \frac{{{\text{1}}{\text{.3}}{E_{\text{s}}}}}{{{\text{2}}r{\text{(1}} - \nu _{\text{s}}^{\text{2}}{\text{)}}}}\sqrt[{{\text{12}}}]{{\frac{{{E_{\text{s}}}{{{\text{(2}}r{\text{)}}}^{\text{4}}}}}{{{{{\text{(}}EI{\text{)}}}_{{\text{eq}}}}}}}}\text{，} (30) {k_{\text{p}}} = \frac{{1.3{E_{\text{p}}}}}{{2r(1 - \nu _{\text{p}}^{\text{2}})}}\sqrt[{12}]{{\frac{{{E_{\text{p}}}{{(2r)}^4}}}{{{{(EI)}_{{\text{eq}}}}}}}}。 (31) 式中：Es为正常围岩地层弹性模量；\nu _{\text{s}}^{}为正常围岩地层泊松比；Ep为断层破碎带地层弹性模量；\nu _{\text{p}}^{}为断层破碎带地层泊松比。
3. 工程实例
3.1 工程概况
珠三角水资源配置工程狮子洋段海底输水隧洞长2252 m，埋深40 m，采用泥水平衡盾构法施工，最大水深为27 m。管片采用C50钢筋混凝土管片错缝拼装，外径8.3 m，内径7.5 m，管片环宽度1.6 m。隧洞穿越f114断层破碎带，该断面水深10 m。断层宽度28.5 m，倾角约为40°。表 4给出断层破碎带与周边岩层的力学参数，其中 \gamma 为重度， \varphi 为内摩擦角，c为黏聚力，E为弹性模量， \nu 为泊松比，k为静止侧压力系数。断层破碎带的形态及其与盾构掘进面相对位置关系如图 6所示。
表 4 断层破碎带及围岩力学参数Table 4. Mechanical parameters of fault fracture zone and surrounding rock类型 \gamma {\text{/}}{\text{(kN}} \cdot {{\text{m}}^{ - {\text{3}}}}{\text{)}} \varphi /(°) c{\text{/}}kPa E{\text{/}}GPa \nu k 正常围岩 26 30 100 1.6 0.33 0.49 断层破碎带 20 25 30 0.5 0.40 0.67 由表 3可知，隧道上覆荷载为第四类f4(x)。根据破碎带的上覆荷载及地层的变化将隧道分为5段Timoshenko梁，其范围分别为：(-∞, 0]，(0, bp]，(bp，H/tan \alpha ]，(H/tan \alpha ，bp+H/tan \alpha ]，(bp+H/tan \alpha ，+∞)。纵向沉降微分方程式(20)相关的系数保留两位小数后分别为：c11=-3.49×10-9；c12=7.04×10-9；c13=0；c14=0；c21=4.01×10-10；c22=-2.54×10-10；c23=-14.76；c24= 17.53；c31=-1.26×10-23；c32=-3.22×10-23；c33=2.0×109；c34=8.75×108；c41=-1.14×10-30；c42=-8.73×10-32；c43=1.76×1011；c44=3.60×1011；c51=0；c52=0；c53=7.22×1018；c54=-5.95×1017。
3.2 监测点布设
根据工程地质条件，选取5个代表管片，分别为处于断层破碎带范围内的第1529、1534和1538环，断层破碎带之外的第1551、1561环。采用全分布式光纤传感器对管片进行监测，如图 7所示。
3.3 现场监测试验结果对比分析
采用偏压短柱法[19]通过实测钢筋应变反算盾构管片内力。如图 8所示，将隧洞衬砌考虑成偏心受压短柱，建立力矩平衡方程，计算相关参数见表 5。
M = {n_{\text{c}}}{\text{(}}{P_{\text{2}}} - {P_{\text{1}}}{\text{)}}\left( {\frac{b}{{\text{2}}} - a'} \right){\text{ + }}\frac{{\text{1}}}{{{\text{12}}}}b{\text{(}}{\sigma _{{\text{c1}}}}{\text{ + }}{\sigma _{{\text{c2}}}}{\text{)}}{A_{\text{c}}}。 (32) 表 5 螺栓及钢筋相关参数Table 5. Relevant parameters of bolt and reinforcementa'/mm Ee/GPa n1 n2 Ec/GPa 50 210 20 10 3.15×104 nb Eb/GPa lb/mm Ab/m2 {\alpha _{\text{b}}}{\text{/}}(°) 28 210 730 7.069×10-4 35 注：a'为钢筋保护层厚度；Ee为钢筋弹性模量；n1为管片内侧环向主筋数量；n2为管片外侧环向主筋数量；Ec为管片弹性模量；nb为纵向螺栓个数；Eb为螺栓弹性模量；lb为螺栓长度，Ab为螺栓横截面面积；{\alpha _{\text{b}}}为螺栓倾角。 式中：{P_{\text{1}}}为管片外侧单根钢筋集中力；{P_{\text{2}}}为管片内侧单根钢筋集中力；{\sigma _{{\text{c1}}}}为管片内侧钢筋位置处混凝土应力；{\sigma _{{\text{c2}}}}为管片外侧钢筋位置处混凝土应力；{n_{\text{c}}}为管片内(外)侧环向主筋数量；b为管片厚度；a'为钢筋保护层厚度；Ac为管片横截面面积。
对第1534、1538环管片结构的内力进行分析，管片环宽度为1.6 m，同一环管片上理论计算弯矩值会在一定的范围波动，其与监测弯矩对比结果如图 9，10所示。
![]()
图 9 1529环管片实测弯矩与计算弯矩Figure 9. Measured and calculated bending moments of segments of No. 1529 ring![]()
图 10 1538环管片实测弯矩与计算弯矩Figure 10. Measured and calculated bending moments of segments of No. 1538 ring由图 9，10可知：现场实测得到的隧洞弯矩与理论计算较为接近。1529环根据现场实测钢筋应变反算管片弯矩值在212.2~228.5 kN，在监测后期管片的弯矩值基本稳定在220 kN·m左右。1538环根据现场实测钢筋应变反算管片弯矩值在175.4~193.6 kN·m，1538环管片弯矩值稳定在180 kN·m左右。
4. 参数敏感性分析
4.1 断层破碎带倾角
假定隧道埋深40 m，断层破碎带宽度为30 m，水深10 m，分别对断层破碎带倾角为20°，40°，60°和80°条件下隧道纵向沉降进行分析，如图 11所示。在断层破碎带区域隧道沉降先急剧增大，后缓慢增加、最后急剧减小，即曲线整体呈现三段变化规律：激增段、缓增段和骤减段。在断层破碎带与正常围岩界面附近隧道不均匀沉降急剧增加。
![]()
图 11 不同断层破碎带倾角隧道纵向沉降图Figure 11. Longitudinal settlements of tunnel with different inclination angles of fault fracture zone4.2 断层破碎带宽度
假定隧道埋深40 m，断层破碎带倾角为60°，分别对断层破碎带宽度为20，30，40，50，60 m条件下隧道纵向沉降进行分析，如图 12所示。随着断层破碎带宽度的增大，断层破碎带范围内隧道的最大纵向沉降逐渐增大并趋于稳定，最大沉降出现位置逐渐远离断层破碎带中部。对于8.3 m外径海底盾构隧道，倾斜断层破碎带对延伸侧1倍破碎带宽度范围内隧道的沉降影响较大，当断层破碎带宽度超过50 m时，隧道的最大纵向沉降基本保持不变。
![]()
图 12 不同断层破碎带宽度隧道纵向沉降图Figure 12. Longitudinal settlements of tunnel with different widths of fault fracture zone不同断层破碎带宽度对隧道纵向弯矩影响如图 13所示。随着断层破碎带宽度增加，在断层破碎带和正常围岩界面附近，隧道纵向弯矩基本不变，在断层破碎带内部隧道的纵向弯矩逐渐减小。对于8.3 m外径海底盾构隧道，当破碎带宽度超过40 m，断层破碎带中部几乎无弯矩。
![]()
图 13 不同断层破碎带宽度隧道弯矩图Figure 13. Bending moments of tunnel with different widths of fault fracture zone4.3 断层破碎带弹性模量
假定隧道埋深40 m，断层破碎带宽度为30 m，倾角为60°，分别对断层破碎带与正常围岩弹性模量比值为1∶2，1∶4，1∶6，1∶8进行分析，如图 14所示。随着断层破碎带与正常围岩弹性模量比值的减小，隧道沉降趋势基本不变，但不均匀沉降迅速增大。对于外径为8.3 m海底盾构隧道，当宽度为30 m的断层破碎带与正常围岩弹性模量比值小于1∶4时，隧道的沉降量超过预警值10 mm[20]。
![]()
图 14 不同断层破碎带弹性模量隧道纵向沉降图Figure 14. Longitudinal settlements of tunnel with different elastic moduli of fault fracture zone5. 结论
本文针对断层破碎带海底盾构隧道纵向沉降问题，基于太沙基松动土压力理论，结合泰勒级数展开公式，对断层破碎带海底盾构隧道纵向沉降及内力进行分析，得到以下3点结论。
(1) 通过分析海底盾构隧道穿越断层破碎带地层全过程，提出了可考虑断层破碎带宽度和倾角的统一荷载模型和理论公式。
(2) 针对复杂上覆荷载条件下Timoshenko梁难以直接求解问题，利用泰勒级数展开公式提出了高精度泰勒公式分段展开法，结合Timoshenko梁理论，推导建立了可用于复杂荷载作用下海底盾构隧道管片纵向沉降及内力系列解析解。工程实例表明本文理论推导的系列解析解能够较为准确地计算管片结构内力。
(3) 计算结果表明，断层破碎带内隧道的沉降曲线表现为激增段、缓增段和骤减段。断层破碎带对延伸侧1倍破碎带宽度范围内隧道的沉降影响较大；当断层破碎带宽度超过6倍洞径时，隧道的沉降峰值基本保持不变。
-
![]()
图 7 不同时刻冰透镜体分布图（白色代表冰透镜体）
Figure 7. Distribution of ice lens at different time (white color indicates ice lens)
![]()
图 9 不同时刻未冻水膜压力随位置变化图(试验4)
Figure 9. Charge of water film pressure with position at different time (#4)
![]()
图 12 不同模拟工况下含水率结果图
Figure 12. Results of water content under different simulation conditions
![]()
图 13 不同模拟工况下冻胀量与冻深结果图
Figure 13. Results of frost heave and depth under different simulation conditions
表 2 冻胀模型输入参数
Table 2 Input parameters in THM model
参数 取值 参数 取值 参数 取值 试验1：k0 /(m·s-1) 3.0×10-10 试验3：T0/℃ -0.3 /(W·(m·K)-1) 0.58 试验2：k0 /(m·s-1) 2.5×10-10 试验4：T0/℃ -0.3 Cs/(kJ·(m3·K) -1) 2160 试验3：k0 /(m·s-1) 2.5×10-10 试验1：/(kg·m-3) 2720 Ci/(kJ·(m3·K) -1) 1874 试验4：k0 /(m·s-1) 2.5×10-10 试验2：/(kg·m-3) 2700 Cw/(kJ·(m3·K) -1) 4180 试验1：Es/MPa 2.5 试验3：/(kg·m-3) 2700 L/(kJ·kg-1) 334.56 试验2：Es/MPa 1.2 试验4：/(kg·m-3) 2768 1 试验3：Es/MPa 0.8 /(kg·m-3) 917 b 3 试验4：Es/MPa 12 /(kg·m-3) 1000 α -5 试验1：T0/℃ 0 /(W·(m·K)-1) 1.2 β -8 试验2：T0/℃ -0.5 /(W·(m·K) -1) 2.22 /kPa 12 
下载: 导出CSV
表 3 不同试验的土性与边界条件参数
Table 3 Soil properties and boundary conditions of different tests
试验编号 顶部温度/℃ 底部温度/℃ 初始温度/℃ 底部状态 初始含水率/% 土的种类 干密度/(g·cm-3) 上覆压力/ kPa 土体高度/cm 冻结时间/h 来源 #1 -3℃ 3℃ 3℃ 补水 17 青海-西藏黏土 1.70 100 10 72 文献[16] #2 -1.7℃ 1℃ 1℃ 封闭 25.71 内蒙古粉质黏土 1.63 0 15 120 文献[26] #3 -2.3℃ 0.7℃ 0.7℃ 补水 18.9 黏土 1.70 0 10 72 文献[28] #4 1→-1 20 ℃/h
-1→-8 0.1 ℃/h1℃ 1℃ 补水 21.7 内蒙古粉质黏土 1.70 500 5 72 文献[24]
下载: 导出CSV
表 4 不同模拟工况参数表
Table 4 Parameters of different simulation conditions
工况 顶温/℃ 上覆压力/ kPa 渗透系数/ (10-10 m·s-1) 压缩模量/ MPa #1 -3 100 3 2.5 #5 -1.5 100 3 2.5 #6 -4.5 100 3 2.5 #7 -3 50 3 2.5 #8 -3 200 3 2.5 #9 -3 100 30 2.5 #10 -3 100 1 2.5 #11 -3 100 3 0.5 #12 -3 100 3 5.0
下载: 导出CSV
-
[1] TABER S. Frost heaving[J]. The Journal of Geology, 1929, 37(5): 428-461. doi: 10.1086/623637
[2] MILLER R D. Freezing and heaving of saturated and unsaturated soils[J]. Highway Research Record, 1972, 393: 1-11
[3] GILPIN R R. A model for the prediction of ice lensing and frost heave in soils[J]. Water Resources Research, 1980, 16(5): 918-930. doi: 10.1029/WR016i005p00918
[4] FOWLER A C. Secondary frost heave in freezing soils[J]. SIAM Journal on Applied Mathematics, 1989, 49(4): 991-1008. doi: 10.1137/0149060
[5] REMPEL A W, WETTLAUFER J S, WORSTER M G. Premelting dynamics in a continuum model of frost heave[J]. Journal of Fluid Mechanics, 2004, 498: 227-244. doi: 10.1017/S0022112003006761
[6] KONRAD J M, MORGENSTERN N R. A mechanistic theory of ice lens formation in fine-grained soils[J]. Canadian Geotechnical Journal, 1980, 17(4): 473-486. doi: 10.1139/t80-056
[7] HARLAN R L. Analysis of coupled heat-fluid transport in partially frozen soil[J]. Water Resources Research, 1973, 9(5): 1314-1323. doi: 10.1029/WR009i005p01314
[8] 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. doi: 10.1139/t78-058
[9] 原国红. 季节冻土水分迁移的机理及数值模拟[D]. 长春: 吉林大学, 2006. YUAN Guohong. The Mechanism and Numerical Simulation of Water Transfer in Seasonal Freezing Soil[D]. Changchun: Jilin University, 2006. (in Chinese)
[10] 李智明. 冻土水热力场耦合机理研究与工程应用[D]. 哈尔滨: 哈尔滨工业大学, 2017. LI Zhiming. Study on Mechanisum of Moisture-Heat-Stress Coupling for Frozen Soil and Engineering Application[D]. Harbin: Harbin Institute of Technology, 2017. (in Chinese)
[11] 白青波, 李旭, 田亚护, 等. 冻土水热耦合方程及数值模拟研究[J]. 岩土工程学报, 2015, 37(增刊2): 131-136. http://manu31.magtech.com.cn/Jwk_ytgcxb/CN/abstract/abstract16270.shtml BAI Qingbo, LI Xu, TIAN Yahu, et al. Equations and numerical simulation for coupled water and heat transfer in frozen soil[J]. Chinese Journal of Geotechnical Engineering, 2015, 37(S2): 131-136. (in Chinese) http://manu31.magtech.com.cn/Jwk_ytgcxb/CN/abstract/abstract16270.shtml
[12] THOMAS H R, CLEALL P, LI Y C, et al. Modelling of cryogenic processes in permafrost and seasonally frozen soils[J]. Géotechnique, 2009, 59(3): 173-184. doi: 10.1680/geot.2009.59.3.173
[13] ZHOU J Z, LI D. Numerical analysis of coupled water, heat and stress in saturated freezing soil[J]. Cold Regions Science & Technology, 2012, 72: 43-49.
[14] HUANG X, RUDOLPH D L. Coupled model for water, vapour, heat, stress and strain fields in variably saturated freezing soils[J]. Advances in Water Resources, 2021, 154: 103945. doi: 10.1016/j.advwatres.2021.103945
[15] LAI Y M, PEI W S, ZHANG M Y, et al. Study on theory model of hydro-thermal–mechanical interaction process in saturated freezing silty soil[J]. International Journal of Heat & Mass Transfer, 2014, 78: 805-819.
[16] MING F, LI D Q. Experimental and theoretical investigations on frost heave in porous media[J]. Mathematical Problems in Engineering, 2015: 1-9.
[17] 曾桂军, 张明义, 李振萍, 等. 正冻土中冰透镜体形成力学判据的分析讨论[J]. 冰川冻土, 2015, 37(1): 192-201. ZENG Guijun, ZHANG Mingyi, LI Zhenping, et al. Review of mechanical criterion for formation office lens in freezing soil[J]. Journal of Glaciology and Geocryology, 2015, 37(1): 192-201. (in Chinese)
[18] NIXON, DERICK J. Discrete ice lens theory for frost heave in soils[J]. Canadian Geotechnical Journal, 1991, 28(6): 843-859. doi: 10.1139/t91-102
[19] TICE A R, ANDERSON D M, BANIN A. The prediction of unfrozen water contents in frozen soils from liquid limit determinations[C]// Symposium on Frost Action on Roads. Oslo, 1973.
[20] 季雨坤. 冰透镜体生长机制及水热力耦合冻胀特性研究[D]. 徐州: 中国矿业大学, 2019. JI Yu-kun. Ice Lens Growth Mechanism and Hydro-Thermal-Mechanical Coupling Research on Frost Heave[D]. Xuzhou: China University of Mining and Technology, 2019. (in Chinese)
[21] ZHOU Y, ZHOU G Q. Intermittent freezing mode to reduce frost heave in freezing soils-experiments and mechanism analysis[J]. Canadian Geotechnical Journal, 2012, 49(6): 686-693. doi: 10.1139/t2012-028
[22] SHENG D C, ZHANG S, YU Z, et al. Assessing frost susceptibility of soils using PCHeave[J]. Cold regions science and technology, 2013, 95(11): 27-38.
[23] TENG J D, LIU J L, ZHANG S, et al. Frost heave in coarse-grained soils: experimental evidence and numerical modelling[J]. Géotechnique, 2022: 1-12
[24] 曹宏章. 饱和颗粒土冻结过程中的多场耦合研究[D]. 北京: 中国科学院研究生院, 2006. CAO Hongzhang. Research on Fields Coupling in Saturated Granular Soil Freezing Process[D]. Beijing: University of Chinese Academy of Sciences, 2006. (in Chinese)
[25] TENG J D, YAN H, LIANG S, et al. Generalising the kozeny-carman equation to frozen soils[J]. Journal of Hydrology, 2020, 594: 125885.
[26] 徐学祖, 王家澄, 张立新. 冻土物理学[M]. 北京: 科学出版社, 2001. XU Xuezu, WANG Jiacheng, ZHANG Lixin. Permafrost Physics[M]. Beijing: Science Press, 2001. (in Chinese)
[27] O'NEILL K, MILLER R D. Exploration of a rigid ice model of frost heave[J]. Water Resources Research, 1985, 21(3): 281-296.
[28] 周家作. 土在冻融过程中水、热、力的相互作用研究[D]. 北京: 中国科学院研究生院, 2012. ZHOU Jiazuo. Study on the Interaction of Water, Heat and Force in Soil during Freezing and Thawing[D]. Beijing: Graduate University of Chinese Academy of Sciences, 2012. (in Chinese)
-
期刊类型引用(3)
1. 邹弈，朱碧堂，吴颖彪，周宇航. 断层破碎带黏性夹泥岩体与基岩接触界面大型直剪试验研究. 华东交通大学学报. 2025(01): 45-51 . 
百度学术
2. 杨会林. 隧道穿越碎裂大涌水岩体变形破坏机理及处治对策. 铁道建筑技术. 2024(07): 121-124 . 百度学术
3. 王志岗，陶连金，史明，刘建功，石城，张海祥，郭文远，谢霖. 穿越活动断裂带的城市地铁隧道抗断性能及损伤状态量化指标研究. 应用基础与工程科学学报. 2024(06): 1772-1786 . 百度学术
其他类型引用(6)
-
其他相关附件
计量
- 文章访问数: 417
- HTML全文浏览量: 90
- PDF下载量: 124
- 被引次数: 9
