«上一篇
文章快速检索     高级检索
下一篇»
  哈尔滨工程大学学报  2017, Vol. 38 Issue (8): 1231-1237  DOI: 10.11990/jheu.201607027
0

引用本文  

李刚, 张金利, 杨庆, 等. 基于非线性渗流的大变形固结有限元分析[J]. 哈尔滨工程大学学报, 2017, 38(8), 1231-1237. DOI: 10.11990/jheu.201607027.
LI Gang, ZHANG Jinli, YANG Qing, et al. Finite element analysis of large-strain consolidation with nonlinear flow[J]. Journal of Harbin Engineering University, 2017, 38(8), 1231-1237. DOI: 10.11990/jheu.201607027.

基金项目

国家自然科学基金项目(41572252)

通信作者

杨庆, E-mail:qyang@dlut.edu.cn

作者简介

李刚(1983-), 男, 讲师, 博士;
杨庆(1964-), 男, 教授, 博士生导师

文章历史

收稿日期:2016-07-11
网络出版日期:2017-04-28
基于非线性渗流的大变形固结有限元分析
李刚1,2, 张金利2,3, 杨庆2,3, 蒋明镜4    
1. 西京学院 土木工程学院, 西安 710123;
2. 大连理工大学 建设工程学部土木工程学院 岩土工程研究所, 大连 116024;
3. 大连理工大学 海岸和近海工程国家重点实验室, 大连 116024;
4. 同济大学 地下建筑与工程系, 上海 200092
摘要:为了研究土体在非线性渗流下的固结效应,基于大变形固结理论,考虑非线性渗流(耦合非Darcy渗流与变渗透系数)作用,建立了二维固结控制方程与有限元方程。采用自编程序对某一黏性土体的固结过程进行了计算。通过大量变动参数的计算与分析,探讨了非Darcy渗流与渗透系数随固结变化的影响。计算结果表明:与Darcy渗流相比,非Darcy渗流下的沉降小,超孔压消散缓慢,具有延迟效应。其中,非Darcy渗流试验常数m影响较大,线性渗流起始水力梯度iL影响较小。渗透系数变化对固结的影响与非Darcy渗流相似,压缩指数与渗透系数指数比值αc的影响较大。在非Darcy渗流与渗透系数随固结变化的耦合作用下,沉降显著减小,超孔压消散趋于缓慢,延迟效应更为明显。
关键词大变形    非Darcy渗流    变渗透系数    有限元分析    超孔压    固结    沉降    压缩指数    
Finite element analysis of large-strain consolidation with nonlinear flow
LI Gang1,2, ZHANG Jinli2,3, YANG Qing2,3, JIANG Mingjing4    
1. School of Civil Engineering, Xijing University, Xi'an 710123, China;
2. Institute of Geotechnical Engineering, School of Civil Engineering, Faculty of Infrastructure Engineering, Dalian University of Technology, Dalian 116024, China;
3. State Key Laboratory of Coastal and Offshore Engineering, Dalian University of Technology, Dalian 116024, China;
4. Department of Geotechnical Engineering, Tongji University, Shanghai 200092, China
Abstract: To study consolidation behaviors under non-Darcy flow, finite element equations and governing equations of plane consolidation were deduced with nonlinear flow (coupled non-Darcy flow with variable permeability coefficient) based on large-strain consolidation theory. The consolidation process of clays was studied by using a program that was developed based on finite element equations. The influence of non-Darcy flow and variable permeability coefficient on consolidation was analyzed by changing the parameters. Calculation results indicate that compared with the circumstances that consider Darcy flow, the settlement is smaller and the dissipation of excess pore water pressures is slower when non-Darcy flow is considered. The influence of parameter m (test constant of non-Darcy flow) on consolidation is greater than iL (initial hydraulic gradient of linear flow). The influence of permeability coefficient is similar to that of non-Darcy flow, and the influence of αc (ratio between compression index and permeability coefficient index) is remarkable. The settlement decreased obviously when considering non-Darcy flow and variable permeability coefficient. The dissipation of excess pore water pressure becomes slow, and the delay effects are more evident.
Key words: large strain    non Darcy flow    variable permeability coefficient    finite element analysis    excess pore water pressure    consolidation    settlement    compression index    

在描述多孔介质渗流问题时,多数采用Darcy定律。该定律主要用于描述层流问题,当雷诺数较大时,不再适用。黏性土试验结果表明,当水力梯度较小时,土体中的渗流存在偏离Darcy定律的现象,称之为非Darcy渗流[1]。TEH等将非Darcy渗流方程引入到固结方程中,分析了非Darcy渗流对固结性状的影响[2-3]。李传勋等[4-10]基于Hansbo非Darcy渗流模型,考虑变荷载与应力历史的影响,建立了以超孔压为变量的软土固结模型,分析了Hansbo模型参数对地基固结速率的影响。崔莉红等采用渗流-固结-溶质运移耦合模型研究了华北平原弱透水层非Darcy渗流固结过程,并对咸水下移的影响进行了分析[11],指出非Darcy渗流在增加黏性土固结时间的同时,也减缓了咸水下移的速度。李菲菲等基于竖井地基固结理论,假设孔隙水服从指数渗流模式,建立了新的竖井地基固结模型,分析了渗流模型对地基固结的影响[12]。LIU等结合Hansbo非Darcy渗流与Terzaghi固结理论,考虑施工变荷载的影响,用有限体积法重新推导了Terzaghi固结方程,分析了非Darcy渗流与施工速度对黏土固结过程的影响[13]。然而,土体在固结过程中,体积不断减小,孔隙比随之降低,由此导致渗透系数降低,土体排水能力减弱,进而引起渗流速度下降,对土体固结过程产生影响。因此,在固结分析时,应考虑渗透系数随固结变化所产生的影响。本文考虑Hansbo非Darcy渗流与渗透系数变化的耦合影响,建立了二维固结控制方程。通过数值计算分析了非线性渗流对土体沉降与超孔压的影响。

1 渗透系数变化与非线性渗流下的固结控制方程 1.1 考虑渗透系数变化的非线性渗流问题描述

Hansbo[14-15]基于试验给出渗流速度v与水力梯度i之间的关系式为

$ v = \left\{ {\begin{array}{*{20}{c}} {\begin{array}{*{20}{c}} {c{i^m}}\\ {k\left( {i - {i_0}} \right)} \end{array}}&{\begin{array}{*{20}{c}} {i < {i_L}}\\ {i \ge {i_L}} \end{array}} \end{array}} \right. $ (1)

式中:ck分别为曲线段与直线段的渗透系数,cm/s;m为试验常数,取值范围为1.2~1.5;i0为计算起始水力梯度;iL为线性渗流起始水力梯度,取值范围为4~10。

正常固结土的孔隙比e与平均有效应力p′的关系可表示为

$ e = {e_0} - {C_c}\lg \left( {\frac{{p'}}{{{{p'}_0}}}} \right) $ (2)

式中:e0为初始孔隙比;Cc为压缩指数;p0为初始平均有效应力,kPa。

Taylor[16]基于试验结果,给出渗透系数k与孔隙比e的关系为

$ \lg k = \lg {k_0} - \frac{{{e_0} - e}}{{{C_k}}} $ (3)

式中:k0为初始渗透系数,cm/s;e0为初始孔隙比;Ck为渗透系数变化指数。将式(2) 代入式(3) 得

$ k = {k_0}{\left[ {\frac{{p'}}{{{{p'}_0}}}} \right]^{ - {\alpha _c}}} $ (4)

式中:αc=Cc/Ck,取值范围为0.5~1.5。

利用式(1) 中v及其导数在分界点iL处的连续条件,可求得参数c=k/(miLm-1),i0=(m-1)iL/m

将式(4) 代入式(1) 得到考虑渗透系数变化的非Darcy渗流方程为

$ \left\{ \begin{array}{l} v = {k_0}{\left( {\frac{{p'}}{{{{p'}_0}}}} \right)^{ - {\alpha _c}}}{i^m}/\left( {mi_L^{m - 1}} \right),\;\;\;\;\;\;\;\left| i \right| < {i_L}\\ v = {k_0}{\left( {\frac{{p'}}{{{{p'}_0}}}} \right)^{ - {\alpha _c}}}\left( {i - \frac{{m - 1}}{m}{i_L}} \right),\;\;\;\;\;\left| i \right| \ge {i_L} \end{array} \right. $ (5)
1.2 平衡方程

在静荷载下(惯性力为零),平面问题的平衡方程、有效应力原理、本构方程及大变形几何方程可以表示如下

$ \left\{ \begin{array}{l} {\mathit{\boldsymbol{L}}^{\rm{T}}}\mathit{\boldsymbol{S}} + {f_0} = 0\\ \mathit{\boldsymbol{S}} = \mathit{\boldsymbol{S' + T}}{\mathit{\boldsymbol{p}}_{ow}}\\ \mathit{\boldsymbol{S' = DE}}\\ \mathit{\boldsymbol{E}} = {\mathit{\boldsymbol{E}}_L} + {\mathit{\boldsymbol{E}}_N} \end{array} \right. $ (6)

式中:L为微分算子矩阵,S为Kirchhoff应力矩阵,f0为体力荷载矩阵,S′为有效应力矩阵,T为选项矢量,pow为孔隙水压力,D为本构矩阵,E为Green应变矩阵,EL为几何线性应变矩阵,EN为几何非线性应变矩阵,各矩阵的具体表达形式见文献[17]。

1.3 连续方程

在饱和土体中,单元体内水量的变化率在数值上等于土体积的变化率[18],即

$ \frac{{\partial {v_x}}}{{\partial x}} + \frac{{\partial {v_y}}}{{\partial y}} = \frac{\partial }{{\partial t}}\left( {\frac{{\partial u}}{{\partial x}} + \frac{{\partial v}}{{\partial y}}} \right) $ (7)

式中:vxvy分别为xy方向的渗流速度,uv分别为xy方向的位移。式(5) 可表示为

$ \left\{ \begin{array}{l} {v_\psi } = {S_\psi }{k_{\psi 0}}{\left( {\frac{{p'}}{{{{p'}_0}}}} \right)^{ - {\alpha _c}}}{\left| {{i_\psi }} \right|^m}/\left( {mi_L^{m - 1}} \right),\;\;\;\;\;\;\left| {{i_\psi }} \right| < {i_L}\\ {v_\psi } = {S_\psi }{k_{\psi 0}}{\left( {\frac{{p'}}{{{{p'}_0}}}} \right)^{ - {\alpha _c}}}\left( {\left| {{i_\psi }} \right| - \frac{{m - 1}}{m}{i_L}} \right),\;\;\;\;\left| {{i_\psi }} \right| \ge {i_L} \end{array} \right. $ (8)

式中:下标ψ可取xy;当水力梯度方向与微元体正方向相同时,Sψ为1;方向相反时,Sψ为-1。将式(8) 对ψ求导得

$ \left\{ \begin{array}{l} \frac{{\partial {v_\psi }}}{{\partial \psi }} = \left( {A + B} \right)\frac{{{k_{{\psi _0}}}}}{{{\gamma _w}}}\frac{{{\partial ^2}{p_{{\rm{ow}}}}}}{{\partial {\psi ^2}}},\;\;\;\;\;\;\;\;\;\;\;\;\left| {{i_\psi }} \right| < {i_L}\\ \frac{{\partial {v_\psi }}}{{\partial \psi }} = \left( {C - D + E} \right)\frac{{{k_{{\psi _0}}}}}{{{\gamma _w}}}\frac{{{\partial ^2}{p_{{\rm{ow}}}}}}{{\partial {\psi ^2}}},\;\;\;\;\;\;\left| {{i_\psi }} \right| \ge {i_L} \end{array} \right. $ (9)

式中:

$ A = \frac{{{\alpha _c}}}{{mi_L^{m - 1}\gamma _w^{m - 1}p'\frac{{{\partial ^2}{p_{{\rm{ow}}}}}}{{\partial {\psi ^2}}}}}{\left( {\frac{{p'}}{{{{p'}_0}}}} \right)^{ - {\alpha _c}}}{\left| {\frac{{\partial {p_{{\rm{ow}}}}}}{{\partial \psi }}} \right|^{m + 1}} $
$ B = \frac{1}{{i_L^{m - 1}\gamma _w^{m - 1}}}{\left( {\frac{{p'}}{{{{p'}_0}}}} \right)^{ - {\alpha _c}}}{\left| {\frac{{\partial {p_{{\rm{ow}}}}}}{{\partial \psi }}} \right|^{m - 1}} $
$ C = \frac{{{\alpha _c}}}{{p'\frac{{{\partial ^2}{p_{{\rm{ow}}}}}}{{\partial {\psi ^2}}}}}{\left( {\frac{{p'}}{{{{p'}_0}}}} \right)^{ - {\alpha _c}}}{\left| {\frac{{\partial {p_{{\rm{ow}}}}}}{{\partial \psi }}} \right|^2} $
$ D = \frac{{\left( {m - 1} \right){i_L}{\gamma _w}{\alpha _c}}}{{mp'\frac{{{\partial ^2}{p_{{\rm{ow}}}}}}{{\partial {\psi ^2}}}}}{\left( {\frac{{p'}}{{{{p'}_0}}}} \right)^{ - {\alpha _c}}}\left| {\frac{{\partial {p_{{\rm{ow}}}}}}{{\partial \psi }}} \right| $
$ E = {\left( {\frac{{p'}}{{{{p'}_0}}}} \right)^{ - {\alpha _c}}} $

式(9) 也可表示为

$ \frac{{\partial {v_\psi }}}{{\partial \psi }} = {H_\psi }\frac{{{k_{{\psi _0}}}}}{{{\gamma _w}}}\frac{{{\partial ^2}{p_{{\rm{ow}}}}}}{{\partial {\psi ^2}}} $ (10)

其中,

$ {H_\psi } = \left\{ \begin{array}{l} A + B,\;\;\;\;\;\;\;\;\;\;\left| {{i_\psi }} \right| < {i_L}\\ C - D + E,\;\;\;\;\left| {{i_\psi }} \right| \ge {i_L} \end{array} \right. $ (11)

将式(10) 代入式(7) 中得到

$ {H_x}\frac{{{k_{{\psi _0}}}}}{{{\gamma _w}}}\frac{{{\partial ^2}{p_{{\rm{ow}}}}}}{{\partial {\psi ^2}}} + {H_y}\frac{{{k_{{y_0}}}}}{{{\gamma _w}}}\frac{{{\partial ^2}{p_{{\rm{ow}}}}}}{{\partial {y^2}}} = \frac{\partial }{{\partial t}}\left( {\frac{{\partial u}}{{\partial x}} + \frac{{\partial v}}{{\partial y}}} \right) $ (12)

写成矩阵形式为

$ {\mathit{\boldsymbol{T}}^{\rm{T}}}\mathit{\boldsymbol{LH}}{\mathit{\boldsymbol{k}}_0}{\mathit{\boldsymbol{L}}^{\rm{T}}}\mathit{\boldsymbol{T}}{p_{ow}} - \frac{\partial }{{\partial t}}{\mathit{\boldsymbol{T}}^{\rm{T}}}\mathit{\boldsymbol{E}} = 0 $ (13)

式中:H为渗流控制矩阵,H=$\left[ {\begin{array}{*{20}{l}} {{H_x}}&0\\ 0&{{H_y}} \end{array}} \right]$k0为初始渗透系数矩阵,k0=$\frac{1}{{{\gamma _w}}}\left[ {\begin{array}{*{20}{l}} {{k_{x0}}}&0\\ 0&{{k_{y0}}} \end{array}} \right]$

2 渗透系数变化与非线性渗流下的固结有限元方程 2.1 位移与孔压模式选取

采用八结点四边形等参单元进行离散,如图 1所示。其中,位移结点为8个,孔压结点为4个。单元的位移函数与孔压函数为

$ \left\{ \begin{array}{l} \tilde u = \sum\limits_{i = 1}^8 {{N_i}{u_i}} \\ \tilde v = \sum\limits_{i = 1}^8 {{N_i}{v_i}} \\ \tilde p = \sum\limits_{j = 1}^8 {{N_j}{p_j}} \end{array} \right. $ (14)
图 1 平面八结点等参单元 Fig.1 Plane eight-node isoparametric element

式中:$\tilde u,\tilde v,\tilde p$为单元内任一点的位移与孔压近似值;uivipj为结点值;NiNj分别为位移与孔压形函数,可表示为

$ \left\{ \begin{array}{l} {N_{i\left( {1,3,5,7} \right)}} = \frac{1}{4}\left( {1 + {\xi _0}} \right)\left( {1 + {\eta _0}} \right)\left( {{\xi _0} + {\eta _0} - 1} \right)\\ {N_{i\left( {2,6} \right)}} = \frac{1}{2}\left( {1 - {\xi ^2}} \right)\left( {1 + {\eta _0}} \right)\\ {N_{i\left( {4,8} \right)}} = \frac{1}{2}\left( {1 - {\eta ^2}} \right)\left( {1 + {\xi _0}} \right) \end{array} \right. $ (15)
$ {N_{j\left( {1,2,3,4} \right)}} = \frac{1}{4}\left( {1 + {\xi _0}} \right)\left( {1 + {\eta _0}} \right) $ (16)

式中:ξ0=ξξiη0=ηηi

2.2 平衡方程离散

将有效应力原理、本构方程代入到平衡方程中,结合几何方程得到离散后的平衡方程为

$ {\mathit{\boldsymbol{K}}_{DS}}\Delta {a_e} + {\mathit{\boldsymbol{K}}_C}\Delta p_{ow}^e = - \mathit{\boldsymbol{\bar R}} - {\mathit{\boldsymbol{R}}_{SP}} $ (17)

式中:KDS为刚度矩阵,KC为耦合矩阵,R为等效结点荷载矢量,RSP为等效结点力矢量,Δae、Δpowe分别为单元的结点位移与超孔压增量。具体计算式如下

$ {\mathit{\boldsymbol{K}}_{DS}} = {\mathit{\boldsymbol{K}}_D} + {\mathit{\boldsymbol{K}}_S},{\mathit{\boldsymbol{K}}_D} = \int_{{V_0}} {{\mathit{\boldsymbol{B}}^{\rm{T}}}\mathit{\boldsymbol{DB}}{\rm{d}}{V_0}} , $
$ {\mathit{\boldsymbol{K}}_S} = \int_{{V_0}} {{\mathit{\boldsymbol{G}}^{\rm{T}}}\mathit{\boldsymbol{MG}}{\rm{d}}{V_0}} ,{\mathit{\boldsymbol{K}}_C} = \int_{{V_0}} {{\mathit{\boldsymbol{B}}^{\rm{T}}}\mathit{\boldsymbol{T\bar N}}{\rm{d}}{V_0}} , $
$ \mathit{\boldsymbol{\bar R}} = \int_{{V_o}} {{\mathit{\boldsymbol{N}}^{\rm{T}}}{\mathit{\boldsymbol{f}}_0}{\rm{d}}{V_0}} = \int_{{A_o}} {{\mathit{\boldsymbol{N}}^{\rm{T}}}{\mathit{\boldsymbol{q}}_0}{\rm{d}}{A_0}} , $
$ {\mathit{\boldsymbol{P}}_{SP}} = \int_{{V_0}} {{\mathit{\boldsymbol{B}}^{\rm{T}}}\left( {\mathit{\boldsymbol{S' + T}}{p_{ow}}} \right){\rm{d}}{V_0}} . $
2.3 连续方程离散

采用Galerkin加权残值法与差分法分别对连续方程进行空间与时域离散,得到离散后的连续方程为

$ \mathit{\boldsymbol{K}}_C^{\rm{T}}\Delta {a_e} - \theta \Delta t{\mathit{\boldsymbol{K}}_{SP}}\Delta p_{ow}^e = {\mathit{\boldsymbol{R}}_{PQ}} $ (18)

式中:θ为差分系数,取值范围为0.5~1.0,若取θ=0.5,则属中心差分格式;KSP为渗流矩阵;RPQ为边界流量等效结点力矢量。具体计算式如下

$ \begin{array}{l} {\mathit{\boldsymbol{K}}_{SP}} = \int_{{V_0}} {\mathit{\boldsymbol{B}}_s^{\rm{T}}\mathit{\boldsymbol{Hk}}{\mathit{\boldsymbol{B}}_s}{\rm{d}}{V_0}} \\ {\mathit{\boldsymbol{R}}_{PQ}} = \Delta t\left( {{\mathit{\boldsymbol{K}}_{SP}}\Delta p_{ow}^e + {\mathit{\boldsymbol{R}}_Q}} \right) \end{array} $
3 非线性渗流的计算方法

由有限元方程的推导可知,在考虑非线性渗流时,渗流矩阵KSP中增加了一项渗流控制矩阵H,计算时需要根据应力与孔压来更新控制矩阵H。为了计算KSP,在计算的初始时步,设H为单位阵,在随后的计算中,通过前一时步得到的有效应力与结点孔压来计算H,由此得到KSP

为计算水力梯度iψ,在局部坐标系下,如图 1所示,设j时刻的结点孔压为pij(i表示结点,j表示时刻),则j+1时刻的水力梯度ij+1可按下式计算:

$ \left\{ \begin{array}{l} i_x^{j + 1} = \frac{1}{{{\gamma _w}}}\frac{{\partial p_{ow}^j}}{{\partial x}} = \frac{{\left( {p_{1\# }^j + p_{4\# }^j - p_{2\# }^j - p_{3\# }^j} \right)}}{{4{\gamma _w}a}}\\ i_y^{j + 1} = \frac{1}{{{\gamma _w}}}\frac{{\partial p_{ow}^j}}{{\partial y}} = \frac{{\left( {p_{1\# }^j + p_{2\# }^j - p_{3\# }^j - p_{4\# }^j} \right)}}{{4{\gamma _w}b}} \end{array} \right. $ (19)

当水力梯度iψ确定后,$\frac{{\partial {p_{{\rm{ow}}}}}}{{\partial \psi }}$即可求得,而$\frac{{{\partial ^2}{p_{{\rm{ow}}}}}}{{\partial {\psi ^2}}}$可按下式计算:

$ \left\{ \begin{array}{l} \frac{{{\partial ^2}p_{{\rm{ow}}}^{j + 1}}}{{\partial {x^2}}} = \left( {\frac{{\partial p_{{\rm{ow}}}^{j + 1}}}{{\partial x}} - \frac{{\partial p_{{\rm{ow}}}^j}}{{\partial x}}} \right)/\Delta x\\ \frac{{{\partial ^2}p_{{\rm{ow}}}^{j + 1}}}{{\partial {y^2}}} = \left( {\frac{{\partial p_{{\rm{ow}}}^{j + 1}}}{{\partial y}} - \frac{{\partial p_{{\rm{ow}}}^j}}{{\partial y}}} \right)/\Delta y \end{array} \right. $ (20)
4 计算结果与讨论 4.1 数值解与近似解析解比较

针对前几节所建立的有限元方程,编制了二维有限元程序。为验证本文所建立的有限元方程及计算程序的可靠性,参照文献[19],建立了厚度为20 m(变形模量Es=0.6 MPa,泊松比v=0.4,初始渗透系数k0=10-7 m/s),顶部排水条件下的一维软黏土计算模型,分别对50 kPa与100 kPa荷载下的固结沉降进行了计算,计算结果如图 2所示。由图可见,在两种荷载作用下,数值解与近似解析解的沉降曲线基本吻合,表明本文所建立的有限元方程及程序较为可信。

图 2 数值解与近似解析解的比较 Fig.2 Comparison between numerical solutions and approximate analytic solutions
4.2 非Darcy渗流对固结的影响分析

为分析非Darcy渗流对固结的影响,参照文献[19-20]中算例,针对长度、厚度分别为20 m与10 m的黏土层,建立了二维有限元模型。其中,重度γ=10 kN/m3,静止土压力系数K0=0.53,压缩模量Es=0.6 MPa,泊松比v=0.4,初始渗透系数k0=1.0×10-9 m/s。均布荷载作用宽度为5 m,顶部荷载施加部分以外排水,底部及两侧不排水。载荷分10级施加,每级间隔3 d,分别为15、25、35、42.7、50、55、60、65、70、75 kPa,计算总时间为4 000 d。模型采用八结点四边形等参单元进行剖分,共划分为204个单元,671个结点,网格剖分如图 3所示。

图 3 网格剖分 Fig.3 Mesh divided

为分析非Darcy渗流对固结的影响,令αc=0,对以下两种工况进行分析:1)m=1.2,变动iL(4,7,10);2)iL=4,变动m(1.2,1.3,1.4),计算结果如图 4~7所示。由图 4的沉降过程曲线可见,当渗流符合Darcy定律时,土层在荷载作用下所产生的沉降较大,而非Darcy渗流产生的沉降相对较小,且随iL的增加而略有减小,表明非Darcy渗流对土层固结具有一定的延迟作用,在完成相同固结度时,需要耗时更多。图 5给出500 d时,不同iL下的超孔压分布等值线。由图可见,在500 d时,非Darcy渗流与Darcy渗流下的超孔压分布规律基本一致,在模型左侧加荷面下的超孔压相对较大。对于非Darcy渗流情况,该处超孔压略大,且随iL的增大而略有增加,由此可见,非Darcy渗流延缓了超孔压的消散,导致土体固结沉降减小。

图 4 参数iL对沉降的影响 Fig.4 Influence of iL on settlement
图 5 参数iL对超孔压的影响(500 d) Fig.5 Influence of iL on excess pore water pressure (500 d)
图 6 参数m对沉降的影响 Fig.6 Influence of m on settlement
图 7 参数m对超孔压的影响(500 d) Fig.7 Influence of m on excess pore water pressure (500 d)

当固定参数iL,变动参数m时,非Darcy渗流对固结过程的影响如图 6图 7所示。由图 6的沉降过程曲线可见,参数m对固结的影响规律与参数iL相似,随参数m的增加,沉降减小。与iL相比,参数m的影响相对较大,其对固结的延迟效应也更为明显。图 7给出500 d时的超孔压分布等值线。由图可见,随m的增加,非Darcy渗流现象越明显,超孔压消散越缓慢,延迟效应也越显著。

4.3 渗透系数变化对固结的影响分析

土体在荷载作用下,因排水而压密,渗透系数随固结发展而减小。因此,固定参数iLm(iL=4.0,m=1.0),变动参数αc(0.5,1.0,1.5),分析渗透系数变化对固结过程的影响,计算结果如图 8图 9所示。由图 8沉降过程曲线可见,加载初期,渗透系数变化的影响较小;随固结时间的增加,渗透系数变化对沉降的影响逐步显现,且随αc的增大沉降相应减小。由式(4) 可知,αc取值较大时,尽管渗透系数略有变化,但沉降减小较为明显,与不考虑渗透系数变化情况相比,达到相同沉降时的历时延长,延迟效应明显。图 9给出500 d时的超孔压分布曲线。由图可见,超孔压分布规律基本相似,随αc的增加,超孔压消散变缓,主要是由于渗透系数变小所致。

图 8 参数αc对沉降的影响 Fig.8 Influence of αc on settlement
图 9 参数αc对超孔压的影响(500 d) Fig.9 Influence of αc on excess pore water pressure (500 d)
4.4 非Darcy渗流与渗透系数变化对固结的影响

为分析非Darcy渗流与渗透系数随固结变化的耦合影响,针对以下5种工况进行计算:1)iL=4,m=1.2,αc=0.5;2)iL=10,m=1.2,αc=0.5;3)iL=4,m=1.2,αc=1.5;4)iL=4,m=1.5,αc=0.5;5)iL=4,m=1.5,αc=1.5,将计算结果与Darcy渗流结果进行对比,如图 10图 11所示。由图 10可见,当考虑非Darcy渗流与渗透系数随固结变化的耦合影响时,在100 d以后,沉降曲线出现较大差异,其中,以工况5下的沉降值为最小。尽管参数mαc取值范围变化不大,但两者对沉降的影响较大,而参数iL的影响则相对较小。可以得出,在黏性土固结分析时,若孔压梯度不大,应考虑非Darcy渗流与渗透系数变化的耦合作用,否则,无法真实地反映黏性土的固结过程。图 11给出500 d时,不同工况下的超孔压分布曲线。如图所示,当同时考虑非Darcy渗流与渗透系数变化的作用时,超孔压分布曲线出现显著差异,在模型左侧超孔压消散明显缓慢,以工况5) 最为明显,由超孔压的分布特点较好地印证了上述的沉降分析。

图 10 非Darcy渗流与变渗透系数对沉降的影响 Fig.10 Influence of non-Darcy flow and variable permeability coefficient on settlement
图 11 非Darcy渗流与变渗透系数对超孔压的影响(500 d) Fig.11 Influence of non-Darcy flow and variable permeability coefficient on excess pore water pressure (500 d)
5 结论

1) 在100 d之后,与Darcy渗流相比,非Darcy渗流下的沉降较小,超孔压消散缓慢,且随固结时间的增加,差异越为明显,非Darcy渗流所引起的延迟效应也越来越显著,其中参数m的影响较大,而参数iL的影响较小,可忽略。原因在于土体固结过程中,土颗粒越来越密实,孔隙越来越小,导致渗透性降低;而随着超孔压的消散,水力梯度随之减小,导致渗流速度降低,由此导致上述规律。

2) 参数αc对固结过程具有重要影响,且随αc的增加,沉降逐渐减小,孔压消散趋于缓慢,延迟效应较为明显。

3) 参数mαc对黏土固结过程的影响较大,当m=1.5,αc=1.5时最为显著。建议在进行黏性土固结分析时,应考虑非Darcy渗流与渗透系数随固结变化的耦合作用。

参考文献
[1] POSKITT T J. The consolidation of saturated clay with variable permeability and compressibility[J]. Geotechnique, 1969, 19(2): 234-252. DOI:10.1680/geot.1969.19.2.234 (0)
[2] TEH C I, NIE X Y. Coupled consolidation theory with non-Darcian flow[J]. Computers and geotechnics, 2002, 29(3): 169-209. DOI:10.1016/S0266-352X(01)00022-2 (0)
[3] HANSBO S. Aspects of vertical drain design:Darcian or non-Darcian flow[J]. Geotechnique, 1997, 47(5): 983-992. DOI:10.1680/geot.1997.47.5.983 (0)
[4] 李传勋, 谢康和. 考虑非达西渗流和变荷载影响的软土大变形固结分析[J]. 岩土工程学报, 2015, 37(6): 1002-1009.
LI Chuanxun, XIE Kanghe. Large-strain consolidation of soft clay with non-Darcian flow by considering time-dependent load[J]. Chinese journal of geotechnical engineering, 2015, 37(6): 1002-1009. DOI:10.11779/CJGE201506005 (0)
[5] 李传勋, 徐超, 谢康和. 考虑非达西渗流和应力历史的土体非线性固结研究[J]. 岩土力学, 2017, 38(1): 91-100.
LI Chuanxun, XU Chao, XIE Kanghe. Nonlinear consolidation of clayed soil considering non-Darcy flow and stress history[J]. Rock and soil mechanics, 2017, 38(1): 91-100. (0)
[6] LI C X, XIE K H, WANG K. Analysis of 1D consolidation with non-Darcian flow described by exponent and threshold gradient[J]. Journal of Zhejiang University:science A, 2010, 11(9): 656-667. DOI:10.1631/jzus.A0900787 (0)
[7] LI C X, XIE K H. One-dimensional nonlinear consolidation of soft clay with the non-Darcian flow[J]. Journal of Zhejiang University:science A, 2013, 14(6): 435-446. DOI:10.1631/jzus.A1200343 (0)
[8] LI C X, XIE K H, HU A F. One-dimensional consolidation of double-layered soil with non-Darcian flow described by exponent and threshold gradient[J]. Journal of Central South University, 2012, 19(2): 562-571. DOI:10.1007/s11771-012-1040-3 (0)
[9] 董兴泉, 李传勋, 陈蒙蒙, 等. 考虑非达西渗流的双层软土地基大变形非线性固结分析[J]. 岩土力学, 2016, 37(8): 2321-2331.
DONG Xingquan, LI Chuanxun, CHEN Mengmeng, et al. Analysis of large-strain nonlinear consolidation of double-layer soft clay foundation with considering effect of non-Darcy's flow[J]. Rock and soil mechanics, 2016, 37(8): 2321-2331. (0)
[10] XIE K H, WANG K, WANG Y L, et al. Analytical solution for one-dimensional consolidation of clayey soils with a threshold gradient[J]. Computers and geotechnics, 2010, 37(4): 487-493. DOI:10.1016/j.compgeo.2010.02.001 (0)
[11] 崔莉红, 成建梅, 路万里, 等. 弱透水层低速非达西流咸水下移过程的模拟研究[J]. 水利学报, 2014, 45(7): 875-881.
CUI Lihong, CHENG Jianmei, LU Wanli, et al. Numerical study on saltwater downward migration in aquitard as low velocity non-Darcy flow[J]. Shuili Xuebao, 2014, 45(7): 875-881. (0)
[12] 李菲菲, 谢康和, 邓岳保. 考虑指数流的真空预压竖井地基固结解析解[J]. 中南大学学报:自然科学版, 2015, 46(3): 1075-1081.
LI Feifei, XIE Kanghe, DENG Yuebao. Analytical solution for consolidation by vertical drains with exponential flow under vacuum preloading[J]. Journal of Central South University:science and technology, 2015, 46(3): 1075-1081. (0)
[13] LIU Z Y, SUN L Y, YUE J C. One-dimensional consolidation of saturated clays under time-dependent loadings considering non-Darcy flow[C]//Soil Behavior and Geo-Micromechanics, Shanghai, 2010:1-7. (0)
[14] HANSBO S. Consolidation equation valid for both Darcian and non-Darcian flow[J]. Geotechnique, 2001, 51(1): 51-54. DOI:10.1680/geot.2001.51.1.51 (0)
[15] HANSBO S. Deviation from Darcy's law observed in one-dimensional consolidation[J]. Geotechnique, 2003, 53(6): 601-605. DOI:10.1680/geot.2003.53.6.601 (0)
[16] TAYLOR D W. Fundamentals of soil mechanics[M]. New York: John Wiley and Sons Inc, 1948. (0)
[17] 谢永利, 潘秋元. 物质构型下的大应变固结方程及其矩阵表述[J]. 重庆交通学院学报, 1994, 13(2): 77-82.
XIE Yongli, PAN Qiuyuan. Large-strain consolidation equation and its matrix form based on material configuration[J]. Journal of Chongqing Jiaotong Institute, 1994, 13(2): 77-82. (0)
[18] 邓岳保, 谢康和, 李传勋. 考虑非达西渗流的比奥固结有限元分析[J]. 岩土工程学报, 2012, 34(11): 2058-2065.
DENG Yuebao, XIE Kanghe, LI Chuanxun. Finite element analysis of Biot's consolidation with non-Darcian flow[J]. Chinese journal of geotechnical engineering, 2012, 34(11): 2058-2065. (0)
[19] 蒋明镜, 沈珠江. 饱和软土的弹塑性大变形有限元平面固结分析[J]. 河海大学学报, 1998, 26(1): 73-77.
JIANG Mingjing, SHEN Zhujiang. Finite element analysis of elasto-plastic large-strain consolidation for saturated Cam-Clay soft soils[J]. Journal of Hohai University, 1998, 26(1): 73-77. (0)
[20] 何开胜, 沈珠江, 彭新宣. 两种Lagrangian大变形比奥固结有限元法及其与小变形法的比较[J]. 岩土工程学报, 2000, 22(1): 30-34.
HE Kaisheng, SHEN Zhujiang, PENG Xinxuan. The comparison of large strain method using Total and Updated Lagrangian finite element formulation and small strain method[J]. Chinese journal of geotechnical engineering, 2000, 22(1): 30-34. (0)