基于分数阶本构模型的隧道工后固结沉降

鄢豪杰， 郑俊杰， 潘 雨， 陈鹏飞

(华中科技大学 土木与水利工程学院， 湖北 武汉 430074)

随着我国城市建设的发展，越来越多的城市开始了地铁的建设，为居民出行带来了极大的便利。然而由于地铁线路走向需要综合考虑政治、经济等多方面因素，因此往往不可避免的需要建在软黏土地层中，导致隧道在开挖完成后很长时间内仍然会发生较大的沉降，长期沉降甚至会占到总沉降量的30%～90%[1]，如不能对此部分沉降加以分析和控制，极有可能会引起隧道结构的变形甚至破坏，造成严重的经济损失，因此非常有必要对隧道开挖完成后的长期沉降进行深入研究。

Terzaghi[2]通过对芝加哥地铁的长期监测，提出土体孔隙水压力消散产生的固结导致了隧道的工后沉降。张冬梅[3]将隧道工后固结问题简化为有侧限的单向压缩问题，并利用三元件模型预测了上海地铁1号线的长期沉降，但部分测点存在较大的误差。魏纲[4]针对张冬梅[3]文中的不合理之处，假定盾构隧道开挖引起的初始超孔隙水压力分布为三角形分布、隧道衬砌为不排水边界，对隧道轴线上方土体的工后沉降进行了研究。杨泽飞等[5]运用分层总和法，采用魏新江等[6]提出的盾构施工引起的初始超孔隙水压力分布，计算了盾构扰动区的固结沉降，但同样没有考虑到隧道衬砌的渗透性。

可以看到，现有的研究大部分都假设隧道衬砌为不排水边界，并且没有考虑地下水位的影响，但Mair等[7]指出影响地表长期沉降的因素主要有四点，即隧道开挖引起的初始超孔隙水压力、土体的物理性质、隧道衬砌渗透性和地下水位的影响。因此，本文针对现有研究的不足之处，假设隧道开挖引起的轴线上方初始超孔隙水压力呈三角形分布，压缩层厚度为隧道上方全部覆土深度，隧道衬砌为半透水边界，考虑地下水位线，并且引入能更好的描述土体蠕变的分数阶三元件本构模型来求解隧道轴线上方土体的工后沉降，该方法可以全面考虑Mair等[7]提出的四个因素。

1 分数阶本构模型

元件模型由于其物理意义明确、模型直观的优点被广泛应用于土体固结的研究中。Gemant[8]通过试验得出土体应力与应变之间存在着分数阶导数关系，随后分数阶微积分理论开始被应用于土体固结的研究中。基于此理论，殷德顺等[9]通过对试验数据的拟合发现分数阶本构关系可以更好的描述土体的蠕变行为。目前，分数阶本构模型在土体的研究中已经得到了广泛的应用。

1.1 分数阶微积分理论

在实际工程中应用最广泛的分数阶算子理论是Riemann-Liouville型分数阶微积分算子理论[10]，该理论对f(t)的α阶积分的定义为：

其中Γ(α)为伽玛函数，α为分数阶次。

1.2 分数阶三元件模型

理想弹性固体本构方程为σt=Eεt(式中：σt为模型受到的应力；εt为模型的应变；E为模型的弹性模量)，可将应力看作正比于应变对时间的零阶导数。牛顿黏壶的本构方程为σt=ηdε/dt(式中：η为模型的黏滞系数)，可将应力看作正比于应变对时间的一阶导数。Koeller[11]提出了一种基于分数阶导数的元件，其本构关系介于上述两种元件之间，本构方程如式(1)所示。当α取0时，元件退化为弹性元件；当α取1时，元件退化为牛顿黏壶。将传统三元件模型中的牛顿黏壶替换为本文中的分数阶元件，就得到了分数阶三元件模型，如图1所示。

图1 分数阶三元件模型

(1)

式中：c=E1-αηα，0≤α≤1。

下面推导分数阶三元件模型在拉普拉斯变换域中的本构方程。根据并联模型总应力等于各部分应力之和，串联模型总应变等于各部分应变之和，可得：

(2)

式中：σt为模型受到的应力；εt,1为模型并联部分的应变；εt,2为串联的弹性元件的应变；εt为模型总应变。

将式(2)进行拉普拉斯变换可得：

(3)

式中：σs，εs,1，εs,2，εs分别为σt，εt,1，εt,2，εt的拉普拉斯变换式。

由式(3)得到分数阶三元件模型本构方程的拉普拉斯变换式为：

(4)

式(4)仅适用于无侧限压缩，对于有侧限压缩的情况，根据张冬梅[3]得出的无侧限和有侧限情况下的对应关系可得有侧限情况下分数阶三元件模型本构方程的拉普拉斯变换式为：

(5)

2 公式推导

2.1 基本假设

张冬梅[3]认为盾构工后沉降主要是由隧道上方1 m范围内土体的压缩变形引起的，因此只考虑隧道上方1 m范围内的土体扰动，并且将隧道上方1 m范围以外的土体自重视为均布荷载，但事实上，盾构引起的初始超孔隙水压力应一直延续到地下水位处[4]。因此，本文取隧道拱顶至地下水位间的土体进行计算，并将地下水位之上的土体自重视为均布荷载，同时作如下假设：

(1)不考虑隧道下方土体的变形，同时隧道结构变形相对土体沉降非常小，因此也不予考虑；

(2)土体压缩层厚度取为隧道上方地下水位以下全部覆土厚度，外荷载只考虑地下水位以上土体的自重；

(3)隧道上方初始超孔隙水压力分布形式为三角形分布；

(4)地表视为完全透水边界，隧道衬砌视为半透水边界。

基于上述假设的计算示意图如图2所示。

图2 计算示意

2.2 公式推导

设隧道衬砌厚度为b；所取土层厚度为h；土体重度为γ；土体渗透系数为k；衬砌渗透系数为k1；隧道开挖引起的拱顶和地表初始超孔隙水压力分别为q2，q1，当q1取0时即退化为本文假设的三角形分布；均布荷载q为地下水位以上土体压力。

采用太沙基一维固结理论中的假设，则土体渗流方程为：

(6)

式中：ε(z,t)为深度z处的土体在t时刻的应变；σ′(z,t)为深度z处的土体在t时刻有效应力的增量；γw为水的重度。

将式(6)进行拉普拉斯变换得：

(7)

由式(5)知分数阶三元件模型本构方程的拉普拉斯变换式为：

(8)

将式(8)代入式(7)得：

(9)

解式(9)得：

(10)

衬砌视为半透水边界，则：

(11)

设深度为z处初始超孔隙水压力为pz，根据有效应力原理由式(11)得：

(12)

将式(12)进行拉普拉斯变换得：

(13)

将式(10)对z求导可得：

(14)

将式(14)代入式(13)得：

(15)

对于本节假设有：

(16)

将式(16)代入式(15)同时令z=h可得底面的边界条件为：

(17)

地表为全透水边界，则z=0时有：

(18)

将式(18)进行拉普拉斯变换并带入式(10)可得地表的边界条件为：

A1+A2=(q1+q)/s

(19)

联立式(17)(19)可得：

(20)

代入式(10)可得深度z处土体的有效应力为：

(21)

则土体的沉降在变换域内的表达式为：

(22)

将式(22)进行拉普拉斯数值反演即可得到任意时刻隧道工后土体的沉降值，本文采用文献[12]中的方法来反演，并通过欧拉变换来加快收敛速度。具体方法如下：

(23)

式中：a为计算参数，取为6。

对式(23)中的求和部分进行欧拉变换：

(24)

式中：Ek通过如下递推式求得：

k=ndif-1,ndif-2,…,1

(25)

式中：ndif为欧拉求和的项数。

3 对盾构隧道工后固结沉降的预测

上海地铁1号线全长14.57 km，所处土层为淤泥质粉质黏土和淤泥质黏土，区间盾构隧道外径6.2 m，内径5.5 m，覆土厚7 m，地下水位1 m，土体加权平均重度γ=18 kN/m3，渗透系数为k=5.4×10-9m/s，泊松比为ν=0.35，包鹤立[13]根据上海地铁允许渗流速度得出衬砌的渗透系数为k1=3.89×10-11m/s，张冬梅[3]和魏纲[4]分别使用应力路径法和应力释放理论得出上海地铁1号线区间盾构隧道拱顶的初始超孔隙水压力为30.3，32.8 kPa，参考两人的结果将拱顶初始超孔隙水压力取为30 kPa。本节将通过2.2节中推导的公式对上海地铁1号线土体的工后长期固结沉降进行预测，为了方便与张冬梅[3]的计算结果进行对比，选取P14，N12两个测点进行计算。

3.1 土体力学参数确定

本节使用BP(Back Propagation)神经网络来反演上海地铁1号线土体的力学参数，待确定的参数有E1m，E2m，ηm，α。通过一定的试算将每个参数取5个水平，然后采用正交设计法[14]设计25组计算参数组合方案作为输出层样本，将25组方案计算出的沉降量作为神经网络的输入层样本，进行神经网络训练。训练完成后将实测沉降数据输入即可反演出土体参数如表1所示。

表1 反演出的土体参数

反演出的地表沉降与实测数据以及张冬梅方法计算结果的对比如图3，4所示。从图中可以看出，本文计算出的结果与实测数据趋势基本一致，并且相比于张冬梅文中的方法更接近实测值。张冬梅计算结果中P14测点前期的沉降值明显比实测数据小，说明该方法低估了土体前期的固结速度，而本文所用的分数阶三元件模型则大大改进了这一缺陷。总体来看，采用分数阶三单元本构模型来模拟隧道上方土体工后长期固结沉降是可行的，并且比普通的三元件本构模型精度更高。

图3 P14测点沉降对比

图4 N12测点沉降对比

3.2 沉降预测

通过反演出的土体力学参数及2.2节中推导出的公式对上海地铁1号线P14，N12两个测点的固结沉降进行预测，结果如图5，6所示。可以看到，两个测点的沉降随时间逐渐增大，并分别在600，400 d左右基本达到稳定，P14测点的最终沉降量约为75 mm，N14测点的最终沉降量约为76 mm。同时可以看到在沉降发展初期N12测点比P14测点沉降速度更大，且最终达到稳定所需时间短，这主要是因为N12测点的黏滞系数较小和分数阶次较大。说明土体性质对沉降的发展有很大的影响，土体黏滞系数越小，分数阶次越大，土体沉降越快。

图5 P14测点沉降曲线

图6 N12测点沉降曲线

4 衬砌渗透系数对工后沉降的影响

现有的大部分研究均假设隧道衬砌为不透水边界，然而隧道衬砌的渗透性对地表长期沉降有着非常重要的影响，实际工程中衬砌渗水现象在盾构隧道中时有发生，造成了严重的安全隐患，因此研究隧道衬砌渗透性对长期沉降的影响有着非常重要的意义。分别取k1/k为0.0001，0.001，0.01，0.1计算P14和N12两个测点的沉降，结果如图7，8所示。

图7 P14测点沉降曲线

图8 N12测点沉降曲线

可以看出，随着k1/k的增大，即隧道衬砌渗透性越大，土体初期沉降速率越大，沉降越早达到稳定。当k1/k小于0.001时，衬砌渗透系数的变化对沉降的影响已不再明显，说明当k1/k小于0.001时隧道衬砌可视为完全不透水边界，大于0.001时应视为半透水边界。而上海地铁1号线实际的k1/k为0.0072，应视为半透水边界，这也验证了本文假设的合理性。

5 结 论

本文在假定隧道开挖引起的初始超孔隙水压力为三角形分布、隧道衬砌为半透水边界的前提下，通过分数阶三元件本构模型，使用太沙基一维固结理论和拉普拉斯变换推导出了盾构隧道开挖引起的工后长期沉降解，可用于预测实际工程中盾构隧道开挖引起的土体长期沉降，并得出了以下结论：

(1)分数阶三元件本构模型可用于模拟隧道开挖引起的土体长期沉降，并且比普通的三元件模型具有更高的精度；

(2)工后沉降受土体黏滞系数和分数阶次的影响很明显，黏滞系数越大，分数阶次越小，土体沉降发展越慢，达到稳定所需时间越长；

(3)隧道衬砌的渗透性对长期沉降的影响不可忽略，衬砌的渗透系数越大，固结初期沉降速率越大，沉降值越大。

(4)当k1/k小于0.001时衬砌可视为完全不透水边界，大于0.001时应视为半透水边界。对于上海地铁1号线，应将衬砌视为半透水边界。