基于物质点法的顶管施工模拟分析

陈 涛

(中铁二院工程集团有限责任公司， 成都 610031)

近年来，顶管施工方法因其扰动小、工期短等显著优势，逐渐应用于隧道工程中。鉴于其实用性，很多学者针对开挖效应对管道及土体的影响规律进行了探究，取得了一系列进展与突破。魏纲[1]等人通过现场试验，探究了施工过程中管道的受力性能；房营光[2]等人分别从理论和实测两方面研究了顶管施工扰动区的土体变形情况。同时，随着计算机计算效率的大幅提升和数值模拟软件的不断改进，很多学者尝试采用计算软件模拟顶管施工过程，也取得了很多有益的成果。冯海宁[3]等人通过有限元计算模拟顶管施工，探讨了不同参数对结构受力和土体变形的影响；屠毓敏[4]研究了长距离顶管穿越海堤时的堤面沉降，为工程建设提供了参考。目前的研究多是采用有限元软件进行的，但在分析顶管施工大变形问题时，有限元计算会发生网格畸变，甚至引起数值困难，造成计算难以进行。

针对这一问题，本文采用了一种新的数值计算方法——物质点法，该方法由Sulsky[5]提出，将连续体用一系列质点离散，质点携带了所有物质信息(称为物质点[6])，其运动就代表了连续体的变形过程。物质点法能有效避免由于大变形导致的计算网格畸变问题，因此在模拟大变形破坏时具有独特优势[7-8]。

本文基于广义插值物质点法对顶管施工进行了三维数值模拟，合理调整网格密度并同经验公式进行对比，验证可靠性之后，调整土层强度参数，探究不同土质下顶管施工对土体沉降量的影响规律，据此研究顶管施工的适用性，为工程建设提供参考。

1 物质点法

物质点法由Sulsky从流体力学中引入，采用拉格朗日质点和欧拉网格双重描述，将连续体离散成一系列质点，每个质点携带一块区域的所有信息，称为物质点，其集合就代表整个材料区域。计算网格仅用于动量方程的求解和空间网格的计算，不携带任何物质信息。在每个时间步时，质点和计算网格完全固连，完成动量方程的求解和空间导数的计算。求解后，将计算网格结点的运动量再映射到各质点上，得到这些质点在下一时刻的运动量。由于应力更新可在每个时间步开始时进行，也可在时间步结束时进行，因此根据更新时间的不同，产生了两种相应的形式，分别为USF(初始端应力更新形式)和USL(结束端应力更新形式)，对应的具体计算流程如图1所示。由于计算网格并不携带任何物质信息，因此，完成每个时间步计算后，会用全新的规则网格取代旧网格，从而有效避免了网格畸变问题。

图1 物质点法计算流程图

在更新拉格朗日形式中，质点的特征函数χb(x)定义了质点所占据的空间区域，它是质点当前位置和变形状态的函数，即:

(1)

该特征函数确定了物理量在空间变化的光滑程度。类似地，计算体系的虚功方程可以表示为:

(2)

式中：Ω——物体现时构型所占据的区域；

Ωp——质点p的特征函数在现时构型中的支撑域，或称为质点p在现时构型中所占据的区域。

质量守恒方程为：

(3)

动量方程为：

▽σ+ρb=ρa

(4)

能量方程为：

(5)

式中：ρ——介质密度；

V——速度；

a——加速度；

σ——对称的应力张量；

b——物体所随的体力。

方程中的变量都是时间以及物质点空间位置的函数。

质点的特征函数χb(x)作为试探函数来近似描述系统的物理量，采用背景网格形函数NI(x)为测试函数，并将虚位移δui近似为：

(6)

将式(5)代入式(2)中，并且考虑到δuiI在本质边界Γu上为0，可以得到：

(7)

求解时，首先将粒子状态插值到网格以形成Mi，通过该值获取网格速度vi，网格节点速度vi通过在网格节点上用形函数对周边粒子的动量进行插值，然后通过网格节点质量Mi和动量加权得到：

(8)

(9)

(10)

(11)

根据公式，更新质点信息，进行下一步计算。

2 顶管施工三维模拟

2.1 模型建立

由于土体具有弹塑性性质，因此本构方程选用Drucker-Prager模型，进行不排水条件下模拟开挖。顶管变形相对于土体来说很小，可视为理想弹性体。考虑到模型边界的影响，设置模型尺寸，垂直顶进方向为x轴(-30～30 m)，平行顶进方向为y轴(0～100 m)，土体厚度方向为z轴(-15～0 m)。顶管设置为圆角长方体，外轮廓长10 m，宽4 m，高3 m，结构厚度 0.08 m，中心埋深为7.5 m。依据实际情况，前后左右侧边界施加滑动约束，上侧边界不施加约束，设定为自由表面，底部施加固定约束。选用长度为0.1 m的六面体单元进行网格划分，模型如图2所示。

图2 模型示意图

为简化分析，采用等效土层模拟，开挖面泥水压力为0.25 MPa，不考虑注浆压力，并将隧道周围土层与浆液的混合体等效为一种弱化土层。同时，设置顶管仅从第2层土中穿过，未与其他土层接触，各土层材料参数，如表1所示。

表1 土层材料参数表

由于运动方程在背景网格上求解，因此背景网格的密度将会对物质点法计算结果的准确性产生影响，本文通过多组试算得出，当计算网格结点数目不少于10×104时，计算精度能够得到保证，即进一步增加网格结点数目计算对计算结果无明显影响。兼顾精度和计算效率，结点数目取为10×104个。

2.2 计算结果及分析

Peck教授在对大量隧道开挖地表沉降实测数据进行分析的基础上，系统地提出了地层损失的概念和估算隧道开挖地表下沉的实用方法，称为Peck经验公式法。为验证物质点法计算的可靠性，本文将一组结果(顶管位于地表以下30 m横截面上的沉降量)与Peck经验公式计算值进行对比，如图3所示。从图中可以看出，两者的数据基本一致，吻合度较好。由此说明，本模型可较好地反映顶管施工对土体扰动的基本规律。

图3 y=30 m地表沉降对比

在模拟顶管掘进过程中，以y=30 m处横截面为参考面，绘制顶管前进到不同位置时y=30 m处的地表横向沉降曲线，如图4所示。从图中可以看出，顶管距离y=30 m处-10 m时，前方土体开始微量隆起，是由泥水压力作用引起的；当顶管距此处-5 m时，地表有轻微塌陷趋势，产生微量沉降；当顶管达到y=30 m处时，地表产生明显沉降，最大沉降值约5 mm；当顶管超出y=30 m处5 m时，地表沉降进一步加大，最大沉降值约 8 mm；此后，随着顶管进一步远离，地表沉降变化不大，逐渐趋于稳定。从模拟结果还可以看出，地表横向受影响范围基本保持一致，约为-15～15 m，其中，顶管中轴位置沉降量最大，向两边逐渐减弱，且大致呈对称性关系。顶管掘进后的土体沉降更为明显，为掘进前沉降量的数倍。

图4 顶管在不同位置L时y=30 m处的沉降量

在实际施工过程中，遇到的土体一般都不是单一均质土层。考虑到这一实际情况，为更好地探究顶管施工的适用性，本文采用不同土质建立了4组新模型进行对比，顶管穿越土层的土质情况分别为全部黏土、全部砂土、上层黏土下层砂土和上层砂土下层黏土。模拟计算结果如图5所示。由模拟结果可知，在不考虑其他因素影响的条件下，相对于黏土而言，砂土沉降量更小，在顶管顶进过程中变形更容易控制。但穿越砂土所需顶推力较大。同时发现，在复合土层中穿越时，砂土层和黏土层的位置关系对整体沉降量影响较小。

图5 穿越不同地层时横向地表沉降曲线

综合模拟情况和实际施工情况分析可知，顶管在土质较好的砂土层中穿越时，控制沉降量较为容易，但所需顶推力较大；在土质较好的黏土层中穿越时，沉降量较大。在复合土层中穿越时，地表沉降量介于两种单一均质土层之间，且复合土层的位置关系对沉降量影响较小，同时顶推力也更易控制。若在土质较差的砂土层和粘土层中采用顶管施工，需提前根据当地环境保护要求和地质条件进行一定的地基加固处理措施，以确保施工安全。

另外，影响顶管适用性的因素还有注浆压力、顶推速度等，在实际工程中也需加以考虑。

3 结论

本文基于物质点法，对顶管施工过程进行模拟，得出了如下结论。

(1)顶管施工过程中，地表横向受影响范围约为-3.7D～3.7D(D为顶管长度)。

(2)针对每个横截面而言，顶管掘进的过程中，土层先因泥水压力而隆起，随后因土层扰动及损失而出现沉降，在顶管远离该截面后土层逐渐趋于稳定。

(3)在土质较好的复合土层中，若顶管穿越砂土层，则沉降控制更为容易，但所需顶推力较大，因此需综合考虑各方面条件，以衡量顶管施工的适用性。