二维固结过程中土体各向异性影响的研究★

2022-08-24裴美娟李乐晗

山西建筑 2022年17期

裴美娟，杨鹏，李乐晗

(兰州石化职业技术大学土木工程学院，甘肃兰州 730060)

1 概述

土体的固结过程是岩土工程中重要的研究课题，对于实际工程而言，土体固结过程伴随着土体体积的变化，进而导致地基的沉降，会对工程建设造成较大的影响。再者，岩土体具有多相性质，土颗粒之间的孔隙往往被水所填充，而在外力作用下，土体会产生超孔隙水压力，超孔隙水压力的产生会对土体的强度造成强烈的影响，会导致土体发生剪切破坏，而太沙基的渗流固结原理指出了固结排水过程中孔压的变化规律，可以说固结过程不仅仅体现了土体体积的变化，而且更为重要的是能够描述排水过程中土颗粒之间有效应力的变化规律。

土体的固结过程受到了众多学者的关注与研究，不同土体性质有较大差异，其固结排水特性不同，由于土体渗透性具有较多的影响因素，因此，不同因素下，其渗透、固结性质会表现出差异。邓岳保等[1]研究了温度效应下土体的固结蠕动耦合行为，张玉国等[2]给出了外荷载变化条件下组合桩复合地基的固结计算方法，指出不同边界条件会导致孔压消散的速率及分布出现差异。王洁等[3]建立了竖井的非线性固结模型，并引入了井阻变化的影响。当土体处于欠固结状态时，其固结特性与正常固结土体不同[4]，张志龙等[5-6]基于分段线性方法分析了欠固结土体的大变形非线性固结过程，发现该过程中土体的孔压消散具有明显的非线性性质。贺建清等[7]则研究了有机质土体的次固结行为，有机质含量会对固结特性产生影响，有机质的增加会造成次固结系数的上升。循环荷载作用下，汪磊等[8]研究了循环荷载下非饱和土的固结过程，并采用变换方法给出了解析解。土体性质是具有复杂性的，不同因素影响下都具有不同的特性，因此，一些学者也开展了复杂条件下的固结行为研究[9-11]，这些研究为工程应用提供了技术支持。

目前，对于固结过程中土体各向异性的研究较少，黄土受到成因、气候等因素的影响，具有强烈的竖向节理，竖向渗透性往往比水平渗透性大得多，一维渗流固结理论无法描述出二维空间中的各向异性表现，多数的研究基于各向同性假设开展，无法准确刻画出黄土的固结过程，因此，开展各向异性下二维固结理论的研究，对于评价黄土等特殊土的固结特性具有重要的意义。本文通过Galerkin法进行控制方程的离散，获得其有限元形式，并基于Python语言编程开发专有的有限元求解器，完成各向异性下的二维固结计算，能够有效反映出不同各向异性程度下超孔压的变化，并对比分析了不同异性程度下及排水边界条件二维固结模型的孔压消散规律以及分布特性。

2 控制方程的离散

2.1 控制方程及边界条件

在二维空间中，考虑竖向渗透系数kv与水平渗透系数kh不同，同时理论模型依然满足水流连续性原理，孔隙水的排出依然满足达西定律，这样竖向渗流速度与水平渗流速度可以表示为：

(1)

(2)

其中，u(x,y,t)为孔隙水压力函数，是关于空间坐标和时间的函数；γw为水的重度；vh，vv分别为水平向渗流速度和竖向渗流速度。

而土体体积的变化完全由孔隙水的排出量决定，即有式(3)：

(3)

其中，εv为土体体积应变。

引入土体体积压缩系数,同时由太沙基有效应力原理得到二维渗流固结控制方程,如式(4)所示。

(4)

其中，Ch为水平固结系数,Ch=kh/mvγw;Cv为竖向固结系数,Cv=kv/mvγw;mv为土体体积压缩系数。

方程(4)为各向异性条件下的二维渗流固结控制方程,该方程为标准的二维抛物型偏微分方程,属于非稳态物理场问题,需要在边界条件以及初始条件给出的基础上,才能够完成求解。

2.2 控制方程的离散

为了对控制方程进行求解,考虑到空间维度的复杂性,采用有限元方法进行偏微分方程的求解。基于Galerkin法进行控制方程的离散,获得其有限元形式。首先,单元选择三节点三角形单元,此时,单元内任意点处的场函数值可以通过节点值进行插值表达,见式(5)。

u(x,y,t)=[N]·[u(t)]

(5)

其中，[N]为形函数行向量;[u(t)]为t时刻单元3个节点处的函数值,为列向量。

针对方程(4),构建其加权残差,得到对应的等效积分强形式,见式(6)。

(6)

利用分部积分思路或者格林公式进行式(6)的处理,可以进一步得到等效积分弱形式,见式(7)。

(7)

由于形函数是空间坐标的一次函数,因此求导后为常数,仅与单元节点的坐标值有关,同时,将式(5)代入到式(7)后,可以得到空间离散后的形式,见式(8)。

(8)

由于关于空间的积分是可以直接进行计算,这样将空间项积分记作[B]矩阵,而含有时间项的积分记作[S]矩阵,则得到空间离散后的形式，见式(9)。

(9)

由于该问题属于非稳态问题,利用后向差分进行时间项的离散,得到式(10)。

(10)

式(10)即为离散后的有限元形式，为了实现数值求解，基于Python语言开发数值模块完成计算，采用面向对象的编程范式，构建控制方程实例的类，而单元刚度矩阵的计算以及整体刚度矩阵的集成作为类的方法进行设计，边界条件、初始条件作为属性进行输入，边界条件的引入在每一次迭代计算时进行更新。

3 数值结果分析

3.1 算例设置及参数

为了对比各向异性对于土体固结的影响，设置不同组别下的参数组，具体信息见表1。控制竖向固结系数保持不变，而水平向固结系数为竖向固结系数乘以比例因子，进行各项异性的表达，其中参数组别A中，竖向固结系数等于水平向固结系数，为各向同性，而其他组别，随着比值r的降低，水平向固结系数是随之降低的，典型的黄土含有大量的竖向裂缝，竖向渗透性较强，而水平渗透性低。

表1 固结系数取值

模型长高为5 m×3 m，如图1所示。模型具有4个边界，其中边界1和边界2为不透水边界，为第二类齐次边界条件，顶部地面与大气连通，孔压始终为0，边界4为第一类齐次边界条件。为了考虑边界条件对于各向异性的影响，边界3分别考虑为不透水和透水两种情况，即模型1中边界3为不透水边界，而模型2中边界3与外界连通，孔压为0。采用自编程开发的有限元程序进行计算，网格划分如图2所示，模型节点个数为1 834，模型单元数量为3 506，孔压场为标量场，节点处自由度为1。计算总时间为1e4 s，时间步数为1 000步，模型初始孔隙水压力为20 kPa。

3.2 结果分析与对比

图3为模型1在参数A下末时刻的孔隙水压力分布云图，在该边界条件下，当水平向固结系数降低到参数E时，其孔压分布云图见图4，对比发现，两者区别较小，且分布规律较为一致。在模型1中，由于左侧边界为不透水边界，而渗透系数在水平向的降低，主要影响到水平向渗透性的降低，但是由于边界的限制，两侧边界不透水，因此，在该种条件下土体渗透性的各向异性对于固结过程的影响是较小的。图5为参数A结果与参数E结果的偏差孔压分布图，可以看到，各向异性条件下偏差孔压的幅值较小，在0.02 kPa之内。图6为模型在5组参数下的平均孔隙水压力变化曲线，可以看到各向异性在此种边界条件下，对于模型孔压消散的影响较小，与各向同性下的计算结果偏差不大，平均孔压曲线基本重合。

当边界3为透水边界，其边界条件为第一类齐次边界，图7为模型2在参数A条件下第500时间步的孔隙水压力云图，图8为在参数E下的孔压分布，可以看出，随着水平向渗透性的降低，孔压的消散速度也变慢，图8孔压幅值要高于图7的，反映出边界条件会对各向异性的描述产生一定程度的影响。

为了对比该边界条件下，水平向固结系数降低所带来的影响，分别给出了模型2的参数A与参数E的第100步、第500步及结束时的偏差孔压分布云图，分别见图9～图11。在第100步时，固结系数的降低使得模型整体孔隙水压力偏大，并且在远离排水边界的位置幅值较大，最大达到了6 kPa。随着固结过程的继续发生，孔隙水的持续排出，孔压偏差较大的位置逐步向中间位置移动，而当达到末时刻时，偏差孔压主要集中在右侧边界部分。造成这一种变化规律，主要是受到了右侧边界性质的影响，并且土体水平向固结系数的降低，表示土体水平方向排水能力的降低，在边界条件一定的情况下，土体的各向异性性质会造成固结过程产生差异，尤其是在边界透水情况下，这一差异更为明显。图12给出了在不同参数组下模型固结过程中平均孔隙水压力随时间变化曲线，可以看出，随着水平向固结系数的降低，模型的排水速度逐渐降低，表明在此时边界条件下，各向异性的性质对于模型的整体的孔压分布是有明显影响的。

综上，对于二维空间中的土体固结过程而言，土体的各向异性是会对固结过程中的孔压变化产生一定程度影响的，但是还需要考虑边界条件的影响，如果该方向边界不排水，那么该方向的渗透性对于孔压消散的影响较小。当两个方向均含有排水边界时，各向异性表现得较为明显，与各向同性的结果相比，孔压偏差主要出现在远离排水边界位置处，并且伴随着固结的进行逐步向排水边界发展。