工程科学与技术   2021, Vol. 53 Issue (3): 106-114
抽降水引发的弱透水层非线性固结解析解
李传勋, 江留慧     
江苏大学 土木工程与力学学院,江苏 镇江 212013
基金项目: 国家自然科学基金项目(51878320)
摘要: 抽降水会引发含水层下卧的弱透水层的固结沉降,同时,土体非线性压缩、渗透特性对固结变形的影响也较为明显。故抽降水引发的弱透水层的固结模型中有必要考虑土体的非线性固结特性,但考虑非线性固结特性后弱透水层固结模型的解析求解存在着困难。本文在分析现有非线性固结模型的基础上,获得了弱透水层非线性固结的近似解析解。具体过程为:引入描述土体的非线性压缩、渗透特性的经典非线性关系,在弱透水层初始有效应力为恒值的假定下,采用分离变量法推导了压缩指数Cc与渗透指数Ck比值不等于1时,潜水层大面积抽水引发的弱透水层1维非线性固结近似解析解。当Cc/Ck→1时,本文解可退化为相同条件下已有的Cc/Ck=1时1维非线性固结解,初步地验证了其解答的正确性。将本文解答应用于实际工程的沉降分析中,通过理论计算结果与实测结果的对比分析,进一步表明该弱透水层非线性固结解答应用于实际工程的可靠性。最后,以水位单级等速下降为例,利用本文解计算不同影响因素下弱透水层的固结曲线,分析其非线性固结性状。结果表明:Ck值不变时,Cc值越大(Cc/Ck越大),固结速率越缓慢,但相同时间因子下土层沉降量及最终沉降量均越大;Cc值不变时,Ck值越小(Cc/Ck越大),固结速率越慢,相同时间因子下土层沉降量越小,但其并不影响最终沉降量。水位下降速率越快,弱透水层固结速率越快,但最终沉降量相同。水位下降终值hc越大,沉降速率越快,此时Cc/Ck对固结速率的影响越明显。砂土层自然重度γ越大,相同时间因子下,土层沉降量越大,同时超静孔压消散速率越缓慢,但其并不影响土层内的最终超静孔压。
关键词: 水位下降    非线性固结    固结系数变化    解析解    
Analytical Solutions for Nonlinear Consolidation of Aquitard Induced by the Dropping of Groundwater Table
LI Chuanxun, JIANG Liuhui     
Faculty of Civil Eng. and Mechanics, Jiangsu Univ., Zhenjiang 212013, China
Abstract: The dropping of groundwater table by pumping may cause consolidation settlement of aquitard underlying aquifer, and the nonlinear compressibility and permeability of soils have evident influences on the consolidation deformation. Therefore, the nonlinear compressibility and permeability of soils have to be considered in the nonlinear consolidation model of aquitard caused by groundwater pumping, while difficulties exist in solving this consolidation model. Based on the analysis of the existing nonlinear consolidation models, an approximate analytical solution for the nonlinear consolidation of aquitard was derived, in which by adopting the classical nonlinear relationships and assuming a constant initial effective stress at different depths, an approximate analytical solution for one-dimensional nonlinear consolidation of soils caused by large area pumping of phreatic water layer were obtained using the method of separation of variables when the ratio of compressibility index Cc to permeability index Ck was not equal to 1. When Cc/Ck→1, the solution can be reduced to the existing one-dimensional nonlinear consolidation solution for the case of Cc/Ck =1. The solution was applied to the settlement analysis of engineering case, and the comparisons between theoretical and actual measurement results further showed that the solution is reliable for actual problems. Finally, taking the single-stage and constant-rate dropping of groundwater table as an example, the consolidation settlement curves of aquitard under different factors were calculated by using the above solution and the non-linear consolidation behaviors were analyzed. The results showed that when the value of Ck is constant, the larger Cc is (resulting in larger Cc/Ck), the slower the consolidation rate is, but the greater both the settlement at the same time and the final settlement are. When the value of Cc remains constant, the smaller Ck is (resulting in larger Cc/Ck), the slower the consolidation rate is, and the smaller the settlement at the same time during the consolidation process is, however it does not affect the final settlement of aquitard. The faster the dropping rate of groundwater table is, the faster the consolidation rate of aquitard is, while the final settlements under different dropping rate are same. The settlement rate grows with the increase in the final value of groundwater table (hc), whereas the influence of Cc/Ck on the consolidation rate becomes more evident. The larger the specific weight of soils (γ) is, the larger the settlement at the same time factor is, and the slower the dissipation rate of the excess pore pressure is, which indicates that γ has no influence on the final value of the excess pore pressure.
Key words: dropping of groundwater table    nonlinear consolidation    variation of consolidation coefficient    analytical solution    

随着经济发展及城市化进程的加快,人类对地下水资源的需求日趋增长,其后果是由大规模抽取地下水引发的地面沉降问题日益严重[1-3]。大范围地面沉降给人类生活和生产环境、生命财产等造成危害,并严重制约社会和经济的可持续发展。因此,深入研究由抽取地下水引发的地面沉降问题具有重要的理论和实际意义。

计算机技术的发展极大地促进了Biot固结理论在抽降水引发地面沉降研究中的应用。Lewis等[4]最早建立基于Biot耦合固结理论的有限元模型,并分析了威尼斯地面沉降问题。Bear等[5]也基于Biot固结理论建立含水层抽水引起的地面沉降计算模型,分析了降水引起的区域地面沉降问题。Teatini等[6]利用基于Biot固结理论建立的3维有限元模型,分析了Emilia–Romagna地区的地面沉降问题。陈杰等[7]运用2维Biot固结平面有限元程序计算因超采地下水而引起的地面沉降问题。骆祖江等[8]以Biot固结理论为基础,考虑土体的非线性压缩及渗透特性,建立了浅层地下水开采与地面沉降的3维全耦合模型。由此可见,地下水位变化引起的地面沉降采用Biot 3维固结理论建立的完全耦合模型计算已相对完善。但也需注意到Biot 3维固结理论必须借助于相对复杂的数值计算,模型较复杂、参数较多导致数值计算不仅耗时,而且计算结果往往因人而异,很难在实际工程中得到推广应用。

事实上,抽水引起的场地地面沉降范围与发生主要沉降变形的土层厚度相比往往大几倍、几十倍之多,此时采用1维固结理论计算弱透水层沉降是完全适用的。鉴于此,不少学者考虑土体压缩性和渗透性参数与孔隙比的关系,建立了3维地下水渗流和1维固结理论组合的部分耦合模型,对降水引起的地面沉降进行数值模拟[9-11],但仍不可避免数值计算的局限性。

解析解能简捷地预测地面沉降的时空发展规律,便于工程应用。故学者们相继开展了由抽降水引发的地面沉降相关解析解的研究。如:基于太沙基1维固结理论,骆冠勇等[12]推导了承压水大面积减压引起的上覆弱透水层1维线性固结解析解。谢康和等[13]获得了越流系统中承压层水位大面积瞬时下降引发的弱透水层的1维固结解析解。谢海澜等[14]考虑土体固结过程中的非达西渗流特性,给出了承压层水位瞬时下降的弱透水层固结半解析解。张玮鹏等[15]考虑土中渗流存在的起始比降,推导了承压层水位瞬时下降引发的弱透水层固结解析解。夏长青[16]考虑土体成层性,推导了潜水层和承压层水位随时间线性下降情况下成层软土固结解析解。晏凌晗[17]考虑土体流变特性,推导了潜水层水位瞬时下降引发的弱透水层1维流变固结解析解。上述研究虽然考虑了软黏土的成层性、流变性、非达西渗流特性等,但均假定土体在固结过程中其渗透系数和压缩系数保持不变,而实际土体在固结过程中土体渗透系数和压缩系数随有效应力增加而不断变化。

Li[18]基于非线性应力应变关系推导了越流系统中,水位波动变化时,软土层1维非线性固结解析解,但其未能考虑土体渗透性的变化。张云等[19]考虑了土体的非线性变形及渗透特性,建立了以有效应力为变量的1维非线性地面沉降模型,用半解析法进行了求解。黄大中[20]在压缩系数与渗透系数同比例减小(即固结系数保持不变)的假定下得到了潜水层和承压层水位随时间大面积线性下降时的1维非线性固结解析解。但实际上土体在固结过程中压缩系数与渗透系数的变化往往并不成比例。目前,对由抽降水引发的、能反映土体非线性固结特性的弱透水层非线性固结解析解的研究仍较缺乏。

综上所述,十分有必要深入开展抽降水条件下,考虑土体非线性压缩和渗透特性的固结理论研究,尤其是更适合工程技术人员使用的、形式简单的非线性固结解析解。本文引入土体经典的e–lgσ′和e–lgkv非线性关系描述土体的压缩和渗透特性,建立潜水层大面积抽水引发的下卧弱透水层的非线性固结模型,并获得其近似解析解。该解答不仅进一步完善了非线性固结理论体系,而且其简单易用,便于在实际工程中推广应用。

1 问题描述及固结控制方程建立 1.1 问题描述

图1为潜水层水位下降引发的弱透水层1维固结示意图。如图1所示: $ {\textit{z}}$ 为以弱透水层顶面为坐标原点、向下为正的竖向坐标系,弱透水层厚度为H,底面不透水;假定弱透水层土体初始有效应力沿深度保持不变,其值等于自重应力沿深度的加权平均值,将其记为 ${{\sigma_0 '}} $ ;潜水层发生大面积抽水(降水),h0t=0时刻初始变化水位,最终水位下降值为hc;潜水层水位下降引起弱透水层顶面处总应力和孔压均减小,但孔隙水压力减小更明显,故引发弱透水层中土体有效应力进一步增大,进而引发弱透水土层的固结沉降变形。

图1 潜水层水位下降引发的弱透水层1维固结示意图 Fig. 1 Schematic diagram of one-dimensional consolidation of aquitard caused by the dropping of groundwater table in phreatic aquifer

图2为潜水层中水位变化h(t)与时间t的关系曲线。特殊的水位瞬时下降和单级等速下降模式如虚线所示。

图2 水位下降量与时间的关系 Fig. 2 Relationship between the dropping of groundwater table and time

水位瞬时下降其表达式为:

$h(t){\rm{ = }}{h_{\rm{c}}}$ (1)

水位单级等速下降时表达式为:

$h(t) = \left\{\!\!\!\! {\begin{array}{*{20}{c}} {\dfrac{{{h_{\rm{c}}}}}{{{t_{\rm{c}}}}}t,}\;{t \le {t_{\rm{c}}};} \\ {{h_{\rm{c}}},}\;{t > {t_{\rm{c}}}} \end{array}} \right.$ (2)

式中,tc为水位下降历时。

1.2 固结相关土工参数的推导

引入经典的e–lgσ′和e–lgkv关系描述土体固结过程中的非线性压缩与渗透特性:

$e = {e_0} - {C_{\rm{c}}}\lg \left({{{\sigma '} / {{{\sigma_0 '}}}}} \right)$ (3)
$e = {e_0} + {C_{\rm{k}}}\lg \left({{{{k_{\rm{v}}}} / {{k_{{\rm{v0}}}}}}} \right)$ (4)

式中,CcCk分别为压缩指数和渗透指数,e为孔隙比,e0为初始孔隙比,σ′为土层竖向有效应力,kv为土体的渗透参数, ${\sigma_0 '} $ kv0分别为e–lgσ′和e–lgkv曲线上相对应于e0的土体初始有效应力和初始竖向渗透系数。

对式(3)求导,可得土体体积压缩系数mv为:

${m_{\rm{v}}} = - \frac{1}{{1 + {e_0}}}\frac{{\partial e}}{{\partial \sigma '}} = {m_{\rm{v}}}_0\frac{{{{\sigma_0 '}}}}{{\sigma '}}$ (5)

式中, $ {m_{\rm{v}}}_0 = \dfrac{{{C_{\rm{c}}}}}{{{{\sigma_0 '}}\left({1 + {e_0}} \right)\ln\; 10}}$

式(4)减式(3),可得土体渗透系数kv为:

${k_{\rm{v}}} = {k_{{\rm{v}}0}}{\left({\frac{{{{\sigma_0 '}}}}{{\sigma '}}} \right)^{{{{C_{\rm{c}}}} / {{C_{\rm{k}}}}}}}$ (6)

进而根据土体固结系数 $C_{\rm{v}} $ 定义得到:

${C_{\rm{v}}} = {C_{{\rm{v}}0}}{\left({\frac{{{{\sigma_0 '}}}}{{\sigma '}}} \right)^{{{{C_{\rm{c}}}} / {{C_{\rm{k}}} - 1}}}}$ (7)

式中, ${C_{{\rm{v}}0}} = \dfrac{{{k_{{\rm{v}}0}}}}{{{\gamma _{\rm{w}}}{m_{{\rm{v}}0}}}} = \dfrac{{{{\sigma_0 '}}{k_{{\rm{v}}0}}\left({1 + {e_0}} \right)\ln\; 10}}{{{\gamma _{\rm{w}}}{C_{\rm{c}}}}}$ γw为水的重度。由式(7)可知:当Cc/Ck=1时,固结系数保持不变;Cc/Ck≠1时,固结系数在固结过程中随有效应力的增加而不断改变。

1.3 固结控制方程及求解条件

饱和土体1维非线性固结普遍连续方程为:

$\frac{1}{{{\gamma _{\rm{w}}}}}\frac{\partial }{{\partial {\textit{z}}}}\left({{k_{\rm{v}}}\frac{{\partial u}}{{\partial {\textit{z}}}}} \right) = \frac{1}{{1 + {e_0}}}\frac{{\partial e}}{{\partial t}}$ (8)

式中,u为超静孔隙水压力。

将式(3)对t求偏导可得:

$\frac{{\partial e}}{{\partial t}} = - \frac{{{C_{\rm{c}}}}}{{\sigma '\ln\; 10}}\frac{{\partial \sigma '}}{{\partial t}} = - \frac{{{C_{\rm{c}}}}}{{\sigma '\ln\; 10}}\left(\frac{{\partial \sigma }}{{\partial t}} - \frac{{\partial u}}{{\partial t}}\right)$ (9)

根据有效应力原理,土体有效应力表达为:

$\sigma ' = {\sigma '_0} + q\left(t \right) - u$ (10)

式中, $q\left(t \right) $ 为外界因素改变引起的总应力增量。

将式(7)、(9)和(10)代入式(8)可得1维非线性固结控制方程为:

${\;\;\;\;\;\;\;\;C_{{\rm{v}}0}}\frac{{\sigma '}}{{{{\sigma_0 '}}}}\frac{\partial }{{\partial {\textit{z}}}}\left[ {{{\left({\frac{{{{\sigma_0 '}}}}{{\sigma '}}} \right)}^{{{{C_{\rm{c}}}} / {{C_{\rm{k}}}}}}}\frac{{\partial u}}{{\partial {\textit{z}}}}} \right] = \frac{{\partial u}}{{\partial t}} - \frac{{{\rm{d}}q}}{{{\rm{d}}t}}$ (11)

潜水层水位下降,计算弱透水层顶面处总应力时,水位下降区域内土体重度由饱和重度降为自然重度,导致作用在弱透水层顶面处总应力降低,t时刻总应力增量q为:

$q = (\gamma - {\gamma _{{\rm{sat}}}})h(t)$ (12)

式中,γ为降水后潜水层中砂土的天然重度,γsat为潜水层中土的饱和重度。

弱透水层顶面处水头变化与潜水层水位变化一致,土层的边界条件表达为:

$u\left({0,t} \right) = - {\gamma _{\rm{w}}}h(t),\;t > 0$ (13)
$\frac{{\partial u}}{{\partial {\textit{z}}}} = 0,\;{\textit{z}} = H$ (14)

其初始条件表达为:

${\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;u} = {q_0} = (\gamma - {\gamma _{{\rm{sat}}}}){h_0},\;t = 0$ (15)

黄大中[20]Cc/Ck=1(固结系数在固结过程中保持不变)下,获得该模型解析解。目前,关于Cc/Ck≠1时,潜水层水位下降引发的弱透水层1维非线性固结解析解研究甚少,本文将对此问题展开探究。

2 固结模型的近似解析解 2.1 变量代换后的控制方程及求解条件

引入变量 $w = {\left[ {\left({{{{{\sigma_0 '}} + q - u} / {{{\sigma_0 '}}}}} \right)} \right]^{{{1 - {C_{\rm{c}}}} / {{C_{\rm{k}}}}}}}$ ${p_{\rm{c}}} \!=\! {\gamma _{\rm{w}}}{h_{\rm{c}}}$ $v \!=\! {\left[ {\dfrac{{{{\sigma_0 '}} + q + {\gamma _{\rm{w}}}h(t)}}{{{{\sigma_0 '}}}}} \right]^{1 - {C_{\rm{c}}}/{C_{\rm{k}}}}} \!-\! w$ ${q_{\rm{c}}} \!=\! (\gamma \!- \!{\gamma _{{\rm{sat}}}}){h_{\rm{c}}}$ ${N_{\rm{q}}} = ( {{\sigma_0 '}}{\rm{ + }} $ $ {q_{\rm{c}}} + {p_{\rm{c}}} ) / {{{\sigma_0 '}}}$ ,应用以上变量,式(11)可改写为:

${C_{\rm{v}}}_0w\frac{{{\partial ^2}v}}{{\partial {{\textit{z}}^2}}} = \frac{{\partial v}}{{\partial t}} - R(t)$ (16)

式中, $R(t) = \dfrac{{1 - {C_{\rm{c}}}/{C_{\rm{k}}}}}{{{{\sigma_0 '}}}}{\left({\dfrac{{{{\sigma_0 '}} + q + {\gamma _{\rm{w}}}h(t)}}{{{{\sigma_0 '}}}}} \right)^{ - \frac{{{C_{\rm{c}}}}}{{C{}_{\rm{k}}}}}}\left[ {\dfrac{{{\rm{d}}q}}{{{\rm{d}}t}} + {\gamma _{\rm{w}}}h'(t)} \right]$ h´(t)为h(t)关于t的导数。

控制方程通过变量代换,虽然得到极大的化简,但由于系数w的原因仍然无法求得其精确解。变量w的最小值为1,最大值为 ${\left[ {{{\left({{{\sigma_0 '}} + q + {p_{\rm{c}}}} \right)} / {{{\sigma_0 '}}}}} \right]^{1 - {C_{\rm{c}}}/{C_{\rm{k}}}}}$ ,如果取固结过程中的平均值并令其保持不变,则解析解可求。李传勋等[21]曾采用这种求解思路对软土非线性固结模型进行求解,通过与数值解对比发现,其解答具有较高的可靠性。因此,令 ${w_0} = \dfrac{1}{2}\left[ {1{\rm{ + }}{{\left({\dfrac{{{{\sigma_0 '}} + q + {p_{\rm{c}}}}}{{{{\sigma_0 '}}}}} \right)}^{1 - {C_{\rm{c}}}/{C_{\rm{k}}}}}} \right]$ ,此时非线性固结控制方程式(16)转化为线性偏微分方程,具体为:

${C_{\rm{v}}}_0{w_0}\frac{{{\partial ^2}v}}{{\partial {{\textit{z}}^2}}} = \frac{{\partial v}}{{\partial t}} - R(t)$ (17)

应用以上变量,边界条件和初始条件分别转化为:

$v = 0,\;{\textit{z}} = 0$ (18)
$\frac{{\partial v}}{{\partial {\textit{z}}}} = 0,\;{\textit{z}} = H$ (19)
${\;\;\;\;\;\;\;\;\;v} = {\left[ {{{\left({{\sigma _0'} + {q_0} + {p_0}} \right)} / {{\sigma _0'} }}} \right]^{1 - {C_{\rm{c}}}/{C_{\rm{k}}}}} - 1,\;\;t = 0$ (20)

式中, $ {p_0} = {\gamma _{\rm{w}}}{h_0} $

2.2 超静孔隙水压力解答

根据控制方程及边界条件形式,假定满足控制方程和边界条件的解的形式为:

$v = \sum\limits_{m = 1}^\infty {\sin \left(\frac{{M{\textit{z}}}}{H}\right)} {{\rm{e}}^{ - {\beta _m}t}}\left[ {{B_m} + {C_m}{T_m}(t)} \right]$ (21)

式中:M=(2m–1)π/2,m=1,2,3,…;βmBmCm为待定参数;Tm(t)为时间t的函数,由R(t)确定。

将式(21)代入公式(17)~(20),采用分离变量法求得待定参数为 ${\;\beta _m} = \dfrac{{{M^2}{C_{{\rm{v0}}}}{w_0}}}{{{H^2}}}$ ${B_m} = \dfrac{2}{M} \cdot $ $ \left[ {{{\left({\dfrac{{{\sigma _0'} + {q_0} + {p_0}}}{{{\sigma _0'} }}} \right)}^{1 - {C_{\rm{c}}}/{C_{\rm{k}}}}} - 1} \right]$ ${C_m} = {2 / M}$ Tm(t)表达式为:

${T_m}(t) = \int_0^t {{{\rm{e}}^{{\beta _m}\tau }}R\left(\tau \right){\rm{d}}\tau } $ (22)

式中, $\tau $ 为积分变量。

进而可得超静孔隙水压力的普遍解答为:

$\begin{aligned}[b] u =& {\sigma '_0} + q - {\sigma '_0} \left[ {\left(\dfrac{{{{\sigma_0 '}} + q + {\gamma _{\rm{w}}}h(t)}}{{{{\sigma_0 '}}}}\right)^{1 - {C_{\rm{c}}}/{C_{\rm{k}}}}} - \right.\\ &\left.\displaystyle \sum\limits_{m = 1}^\infty {\sin \left(\dfrac{{M{\textit{z}}}}{H}\right)} {{\rm{e}}^{ - {\beta _m}t}}\left[ {{B_m} + {C_m}{T_m}(t)} \right] \right]^{\frac{{{C_{\rm{k}}}}}{{{C_{\rm{k}}} - {C_{\rm{c}}}}}}\end{aligned}$ (23)

在式(23)基础上,可求得土层按孔压定义的平均固结度Up为:

${U_{\rm{p}}} = \dfrac{{\displaystyle \int_0^H {(\sigma ' - {{\sigma_0 '}}){\rm{d}}{\textit{z}}} }}{{\displaystyle \int_0^H {({{\sigma' _f}} - {{\sigma_0 '}}){\rm{d}}{\textit{z}}} }}=\dfrac{{qH - \displaystyle \int_0^H {u{\rm{d}}{\textit{z}}} }}{{({p_{\rm{c}}} + {q_{\rm{c}}})H}}$ (24)

t时刻土层发生的沉降变形值St为:

${S\!_{\rm{t}}} = \int_0^H {\frac{{{e_0} - {{e}}}}{{1 + {e_0}}}} {\rm{d}}{\textit{z}} = \frac{{{C_{\rm{c}}}}}{{1 + {e_0}}}\int_0^H {\lg \frac{{\sigma '}}{{{{\sigma_0 '}}}}} {\rm{d}}{\textit{z}}$ (25)

土层的最终沉降变形值S为:

${S\!_\infty } = \frac{{{C_{\rm{c}}}}}{{1 + {e_0}}}\int_0^H {\lg \frac{{{{\sigma'_f }}}}{{{{\sigma'_0 }}}}} {\rm{d}}{\textit{z}}{\rm{ = }}\frac{{{C_{\rm{c}}}\lg\; {N_{\rm{q}}}}}{{1 + {e_0}}}H$ (26)

进而得到按变形定义的土层平均固结度Us为:

${U_{\rm{s}}} = \frac{{{S\!_{\rm{t}}}}}{{{S\!_\infty }}} = \frac{{\displaystyle \int_0^H {\lg \dfrac{{\sigma '}}{{{{\sigma'_0 }}}}} {\rm{d}}{\textit{z}}}}{{\lg\; {{N_{\rm{q}}}} \cdot H}}$ (27)
3 特殊降水模式下的解答 3.1 水位瞬时下降模式

图2中虚水平线代表t>0时,潜水层水位发生瞬时下降,水位下降深度为hc。可得该降水模式下总应力增量和求解条件分别为:

$q = \left({\gamma - {\gamma _{{\rm{sat}}}}} \right){h_{\rm{c}}}{\rm{ = }}{q_{\rm{c}}}$ (28)
${\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\left\{\!\!\!\! \begin{array}{l} u = - {\gamma _{\rm{w}}}{h_{\rm{c}}} = - {p_{\rm{c}}},\;{\textit{z}} = 0; \\ \dfrac{{\partial u}}{{\partial {\textit{z}}}} = 0,\;{\textit{z}} = H; \\ u = {q_{\rm{c}}},\;t = 0 \end{array} \right.}$ (29)

将式(29)代入待定参数Bm和式(22),可解得:

    ${B_m} = \dfrac{2}{M}\left({{N_{\rm{q}}^{1 - {C_{\rm{c}}}/{C_{\rm{k}}}}} - 1} \right)$ ${T_m}(t) = 0$

$B_m $ $C_m $ $\;\beta_m $ 代入式(23)得到弱透水层中超静孔隙水压力为:

${\;\;\;\;\;\;\;\;\;\;\;\;\begin{aligned}[b]u =& {\sigma '_0} + {q_{\rm{c}}} - {\sigma '_0}\bigg[ {N_{\rm{q}}^{1 - \frac{{{C_{\rm{c}}}}}{{{C_{\rm{k}}}}}}} - \left({{N_{\rm{q}}^{1 - \frac{{{C_{\rm{c}}}}}{{{C_{\rm{k}}}}}}} - 1} \right)\cdot\bigg. \\ &\bigg.\displaystyle \sum\limits_{m = 1}^\infty {\dfrac{2}{M}} \sin \left(\dfrac{{M{\textit{z}}}}{H}\right){{\rm{e}}^{ - {M^2}{w_0}{T_{\rm{v}}}}} \bigg]^{\frac{{{C_{\rm{k}}}}}{{{C_{\rm{k}}} - {C_{\rm{c}}}}}}\end{aligned}}$ (30)

式中, $ {T_{\rm{v}}} = {{{C_{{\rm{v}}0}}t} / {{H^2}}}$

3.2 水位单级等速下降模式

图2中虚斜线所示,当ttc时,潜水层水位发生等速下降,t>tc时水位下降值恒定为hc。可得该降水方式下总应力增量和求解条件分别为:

$ q=\left\{\!\!\!\! {\begin{array}{cc}\dfrac{{q}_{\rm{c}}}{{t}_{\rm{c}}}t, \; t\le {t}_{\rm{c}}; \\ {q}_{\rm{c}}, \; t>{t}_{\rm{c}}\end{array}} \right.$ (31)
$ {\;\;\;\;\;\;\;\;\;\;\;\;\;\left\{ \!\!\!\!\begin{array}{l} u = \left\{ \!\!\!\!{\begin{array}{*{20}{c}} { - \dfrac{{{p_{\rm{c}}}}}{{{t_{\rm{c}}}}}t,}\;{t \le {t_{\rm{c}}},}\;{\textit{z}} = 0\\ { - {p_{\rm{c}}}},\;{t > {t_{\rm{c}}},}\;{\textit{z}} = 0 \end{array}}; \right.\\ \dfrac{{\partial u}}{{\partial {\textit{z}}}} = 0,\;{\textit{z}} = H;\\ u = 0,\;t = 0 \end{array} \right.}$ (32)

同理,可求得超静孔压表达式为:

$u = {\sigma '_0} + q - {\sigma '_0}{\left[ {\left( {\frac{{{{\sigma '}_0} + q + {\gamma _{\rm{w}}}h(t)}}{{{{\sigma '}_0}}}} \right)^{1 - {C_{\rm{c}}}/{C_{\rm{k}}}}} - \sum\limits_{m = 1}^\infty {\frac{2}{M}} \sin \left(\frac{{M\textit{z}}}{H}\right){{\rm{e}}^{ - {M^2}{w_0}{T_{\rm{v}}}}}{T_m}(t) \right]^{\frac{{{C_{\rm{k}}}}}{{{C_{\rm{k}}} - {C_{\rm{c}}}}}}}$ (33)
${T_m}(t)=\left\{\!\!\!\! {\begin{array}{*{20}{l}} {\left(1 - \dfrac{{{C_{\rm{c}}}}}{{{C_{\rm{k}}}}}\right){{\rm{e}}^{\frac{{ - {M^2}{w_0}{T_{{\rm{vc}}}}}}{{{N_{\rm{q}}} - 1}}}}\left[ \dfrac{{{{\left(1 + \dfrac{{{N_{\rm{q}}} - 1}}{{{T_{{\rm{vc}}}}}}{T_{\rm{v}}}\right)}^{1 - \frac{{{C_{\rm{c}}}}}{{{C_{\rm{k}}}}}}} - 1}}{{1 - {C_{\rm{c}}}/{C_{\rm{k}}}}} + \displaystyle \sum\limits_{k = 1}^\infty {\dfrac{{{{\left({\dfrac{{{M^2}w{}_0{T_{{\rm{vc}}}}}}{{{N_{\rm{q}}} - 1}}} \right)}^k}}}{{k!}}} \dfrac{{{{\left(1 + \dfrac{{{N_{\rm{q}}} - 1}}{{{T_{{\rm{vc}}}}}}{T_{\rm{v}}}\right)}^{k + 1 - \frac{{{C_{\rm{c}}}}}{{{C_{\rm{k}}}}}}} - 1}}{{k + 1 - {C_{\rm{c}}}/{C_{\rm{k}}}}} \right],}\;{{T_{\rm{v}}} \le {T_{{\rm{vc}}}};} \\ {\left(1 - \dfrac{{{C_{\rm{c}}}}}{{{C_{\rm{k}}}}}\right){{\rm{e}}^{\frac{{ - {M^2}{w_0}{T_{{\rm{vc}}}}}}{{{N_{\rm{q}}} - 1}}}}\left[ \dfrac{{{N_{\rm{q}}}^{1 - \frac{{{C_{\rm{c}}}}}{{{C_{\rm{k}}}}}} - 1}}{{1 - {C_{\rm{c}}}/{C_{\rm{k}}}}} + \displaystyle \sum\limits_{k = 1}^\infty {\dfrac{{{{\left({\dfrac{{{M^2}w{}_0{T_{{\rm{vc}}}}}}{{{N_{\rm{q}}} - 1}}} \right)}^k}}}{{k!}}} \dfrac{{{N^{k + 1 - \frac{{{C_{\rm{c}}}}}{{{C_{\rm{k}}}}}}_{\rm{q}}} - 1}}{{k + 1 - {C_{\rm{c}}}/{C_{\rm{k}}}}} \right],}\;{{T_{\rm{v}}} > {T_{{\rm{vc}}}}} \end{array}} \right. $ (34)

式中, $ {T_{\rm{v}}}_{\rm{c}} = {{{C_{{\rm{v0}}}}{t_{\rm{c}}}}/ {{H^2}}} $

3.3 Cc/Ck→1条件下超静孔压解的退化

Cc/Ck→1的条件下,水位单级等速下降时的超静孔压解可退化为文献[20]中相同条件下Cc/Ck=1时解的形式。超静孔压解u的退化过程为:

$\begin{aligned}[b] \mathop {\lim }\limits_{{C_{\rm{c}}}/{C_{\rm{k}}} \to 1} u \!=\!& \mathop {\lim }\limits_{{C_{\rm{c}}}/{C_{\rm{k}}} \to 1} {{\sigma' _0}} + q -\\ &{{\sigma _0'}}{\left[ {{{\left( {\frac{{{{\sigma' _0}} + q + {\gamma _{\rm{w}}}h(t)}}{{{{\sigma '_0}}}}} \right)}^{1 - {C_{\rm{c}}}/{C_{\rm{k}}}}} - v} \right]^{\frac{1}{{1 - {C_{\rm{c}}}/{C_{\rm{k}}}}}}}=\\ &{{\sigma' _0}} \!+\! q \!-\! {{\sigma '_0}}{{\rm{e}}^{{{\mathop {\lim }\limits_{{C_{\rm{c}}}/{C_{\rm{k}}} \to 1} }}\frac{{{{\left( {\frac{{{{\sigma '_0}} + q + {\gamma _{\rm{w}}}h(t)}}{{{{\sigma '_0}}}}} \right)}^{1 - {C_{\rm{c}}}/{C_{\rm{k}}}}} \!-\! 1 \!-\! v}}{{1 - {C_{\rm{c}}}/{C_{\rm{k}}}}}}}=\\ &{{\sigma' _0}} + q - {{\sigma '_0}}{{\rm{e}}^w} \\[-10pt] \end{aligned} $ (35)

式中,w ${\sigma '_0} + q - {\sigma '_0}{{\rm{e}}^w} $ 分别为文献[20]中的变量和超静孔压解。

4 工程实例对比

为验证本文解答的正确性,将其与某工程降水引发的地面沉降观测结果、分层总和法计算结果及数值模拟结果[22]进行对比。

简化文献[22]中第五章表5–2和图5–22所示工程实例土层信息,将杂填土、粉土、黏土、粉砂和细砂土层视为土质均一的潜水层,天然重度取各土层天然重度的厚度加权平均值为19.637 kN/m3,这种简化对于弱透水层的沉降计算无影响。初始水位为5.5 m;隔水层粉质黏土视为弱透水层,弱透水层底面不透水;土层性质详见表1,降水引起的有效应力增量如表2所示。在距离基坑5.2 m处,水位下降6.4 m,忽略潜水层水位下降区域内饱和土体重度变化,根据表2和式(26)可分别计算潜水层和弱透水层沉降量。

表1 土性参数 Tab. 1 Parameters of soil characteristic

表2 有效应力增量计算 Tab. 2 Calculation of effective stress increment

图5 Ck值对沉降曲线的影响 Fig. 5 Influence of Ck on settlement curves

将计算结果与文献[22]对比,如表3所示。可以看出,本文解与分层总和法计算结果相比较更趋近于实测沉降值,可能原因主要体现两个方面:其一,本文认为水位下降区域最终有效应力增加值呈线性分布,其更能体现土体中的真实受力特性,故本文相较于分层总和法中有效应力保持恒值,更能体现土体的真实状况;其二,本文解弱透水层沉降计算考虑了土体压缩性与初始有效应力间的非线性关系。本文计算的沉降充分考虑土层中的实际有效应力,其结果更能反映土体的真实性状。从计算的总沉降量与实测结果对比说明了本文的弱透水层沉降计算解答是合理有效的。

表3 不同计算方法下沉降量对比 Tab. 3 Comparison of settlement results with different calculation methods

5 算例及固结性状分析

某越流系统,砂土层饱和重度γsat=20 kN/m3,自然重度γ=18 kN/m3,水的重度按γw=10 kN/m3,弱透水层土层厚度H=10 m,顶面透水,底面不透水,弱透水层自重应力的平均值 ${\sigma '_0} $ =50 kPa。依据郑刚等[23]反演得到的弱透水层土性参数计算表,取弱透水层压缩指数Cc=0.02,土体初始孔隙比e0=0.7,与初始孔隙比e0相对应的渗透系数kv0=4×10–9 m/s。Tv<0.05时潜水层水位等速下降,而后水位下降值恒定为hc=5 m。

基于上述算例,利用本文解分别计算绘制在不同Cc/Ck值、水位下降速率Tvc、水位下降终值hc及砂土层自然重度γ下,弱透水层固结沉降曲线,分析潜水层水位单级等速下降时弱透水层的1维非线性固结性状。

5.1 Cc/Ck对固结性状的影响

图3为不同Cc/Ck下,弱透水层平均固结度与时间因子的关系曲线。图45分别为改变CcCk值时弱透水层沉降量与时间因子的关系曲线。如图35所示:Cc/Ck越大,相同时间因子下弱透水层的平均固结度越小,即固结速率随Cc/Ck增大而减慢;当Ck值保持不变时,Cc值越大,Cc/Ck越大,虽然此时固结速率会减慢,但相同时间因子下弱透水层的沉降量及最终沉降量仍随Cc/Ck增大而变大;当Cc值保持不变时,Ck值越小,Cc/Ck越大,固结完成前,相同时间因子下土层沉降量越小,但固结完成时的最终沉降量相同。该固结性状可由固结系数表达式(7)和沉降表达式(25)得到解释:Cc/Ck越大,固结系数越小,固结速率减慢,最终沉降量与Cc值成正比,与Ck值无关;Cc值越大,Cc/Ck越大,固结速率减慢,但相同时间因子下沉降量越大;当其他参数一定时,增大Ck,固结速率及沉降速率加快,但不改变沉降量终值。

图3 Cc/Ck对平均固结度Up曲线的影响 Fig. 3 Influence of Cc/Ck on the curves of average consolidation degree Up

图4 Cc值对沉降曲线的影响 Fig. 4 Influence of Cc on settlement curves

5.2 水位下降速率对固结性状的影响

图6为水位下降速率与沉降量的关系曲线。由图6可知:随着Tvc增大(水位下降速率减慢),同一时间因子下沉降量随之减小,但最终沉降量相同;在水位下降历时发生改变而其他参数不变的条件下,水位瞬时下降时弱透水层沉降速率最快。该固结性状可由沉降量表达式(25)解释,水位下降越快,在相同时间因子下,弱透水层有效应力越大,相应的沉降量越大。

图6 水位下降速率对沉降曲线的影响 Fig. 6 Influence of the dropping rate of groundwater table on settlement curves

5.3 不同Cc/Ck下水位下降终值对固结性状的影响

图78分别为不同水位下降终值时,σ′/ ${\sigma '_0} $ 沿深度分布和沉降量随时间变化曲线。由图78可得:相同Cc/Ck和同一深度处的有效应力σ′随hc增大而增大;hc越小,Cc/Ck值对有效应力σ′的影响越小。相同Cc/Ck和时间因子下,hc越大,沉降量越大;相同hc和时间因子下,Ck越大,Cc/Ck越小,沉降量越小,但最终沉降量相同;Cc/Ck对沉降曲线影响随hc减小而减弱,且hc=1不同,Cc/Ck下沉降曲线几乎重合,甚至可忽略。由式(7)和(25)可知,Cc/Ck越小,固结系数越大,固结速率跟沉降速率越快;hc越大,相同时间因子下,有效应力越大,相应的沉降量越大;hc越小,σ′/ ${\sigma '_0} $ 越趋近于1,Cc/Ck对固结速率的影响越小,同时Ck不影响沉降量终值,故而出现上述性状。

图7 Cc/Ck下水位下降终值对 ${{{\sigma '}}} $ / ${{{\sigma_{\bf{0}} '}}} $ 沿深度分布的影响 Fig. 7 Influence of the final dropping values of groundwater table on distribution of ${{{\sigma '}}} $ / ${{{\sigma_{\bf{0}} '}}} $ along depth with different Cc/Ck

图8 Cc/Ck下水位下降终值对沉降曲线的影响 Fig. 8 Influence of the final dropping values of groundwater table on settlement curves with different Cc/Ck

5.4 砂土层自然重度 ${{\gamma}} $ 对固结速率的影响

图910分别为不同砂土自然重度γ下,弱透水层超静孔压随时间和深度变化曲线。

图9 砂土层自然重度 ${{\gamma}} $ 对超静孔压消散的影响 Fig. 9 Influence of specific density of sand layer on dissipation of excess pore pressure

图10 砂土层自然重度对超静孔压沿深度分布的影响 Fig. 10 Influence of specific density of sand layer on distribution of excess pore pressure with depth

图910可发现:相同深度处的超静孔隙水压力随γ增大而增大,说明孔压消散速率随γ增大而减慢;相同条件下,相较于Cc/Ck>1,Cc/Ck<1时,对应的超静孔压值更小,即Cc/Ck越小,孔压消散速率越快。

图1112分别为不同砂土层自然重度γ下,Cv/Cv0沿深度分布曲线和沉降量随时间因子变化曲线。由式(7)和式(25)可知:Cc/Ck<1时,砂土层自然重度γ越大,相同时间因子下有效应力σ′越大,固结系数越大,固结速率越快,相应的沉降量越大;Cc/Ck>1时,砂土层自然重度γ越大,相同时间因子下有效应力σ′越大,固结系数越小,固结速率越慢,但相应的沉降量越大;γ越大,Nq越大,最终沉降量越大。

图11 砂土层自然重度对Cv/Cv0沿深度分布的影响 Fig. 11 Influence of specific density of sand layer on distribution of Cv/Cv0 with depth

图12 砂土层自然重度 ${{\gamma}} $ 对沉降曲线的影响 Fig. 12 Influence of specific density of sand layer on settlement curves

6 结 论

1)基于弱透水层初始有效应力沿深度均匀分布的假定,推导了考虑土体非线性压缩和渗透特性的由潜水层水位随时间大面积均匀下降引发的弱透水层1维非线性固结近似解析解,且当Cc/Ck→1时,超静孔压解可退化为相同模型下Cc/Ck=1时的现有解。

2)Cc/Ck越大,弱透水层固结速率越缓慢。Ck不变时,Cc值越大,固结速率越慢,但沉降速率越快;Cc不变时,Ck越小,固结速率越慢,固结过程中沉降速率越慢,但不影响最终沉降量。

3)水位下降越快,沉降速率越快,但最终沉降量相同。

4)水位下降终值hc越大,沉降速率越快,Cc/Ck对固结速率的影响越显著。

5)Cc/Ck<1时,砂土层自然重度γ越大,相同时间因子下固结系数越大,固结速率越快,相应的沉降量越大。Cc/Ck>1时,砂土层自然重度γ越大,相同时间因子下固结系数越小,固结速率越慢,相应的沉降量越大。

参考文献
[1]
Xu Y S,Shen S L,Cai Z Y,et al. The state of land subsidence and prediction approaches due to groundwater withdrawal in China[J]. Natural Hazards, 2008, 45(1): 123-135. DOI:10.1007/s11069-007-9168-4
[2]
Cui Zhendong,Tang Yiqun. Domestic and international recent situation and research of land subsidence disasters[J]. Northwestern Seismological Journal, 2007, 29(3): 275-278. [崔振东,唐益群. 国内外地面沉降现状与研究[J]. 西北地震工程学报, 2007, 29(3): 275-278.]
[3]
Chen C X,Pei S P,Jiao J J. Land subsidence caused by groundwater exploitation in Suzhou City,China[J]. Hydrogeology Journal, 2003, 11(2): 275-287. DOI:10.1007/s10040-002-0225-5
[4]
Lewis R W,Schrefler B. A fully coupled consolidation model of subsidence of Venice[J]. Water Resources Research, 1978, 14(2): 223-230. DOI:10.1029/WR014i002p00223
[5]
Bear J,Corapcioglu M Y. Mathematical model for regional land subsidence due to pumping:2.Integrated aquifer subsidence equations for vertical and horizontal displacements[J]. Water Resources Research, 1981, 17(4): 947-958. DOI:10.1029/WR017i004p00947
[6]
Teatini P,Ferronato M,Gambolati G,et al. Groundwater pumping and land subsidence in the Emilia–Romagna coastland,Italy:Modeling the past occurrence and the future trend[J]. Water Resources Research, 2006, 42(1): W01406.
[7]
Chen Jie,Zhu Guorong,Gu Aming,et al. Application of Biot consolidation theory to calculation of land subsidence[J]. Hydrogeology and Engineering Geology, 2003, 30(2): 28-31. [陈杰,朱国荣,顾阿明,等. Biot固结理论在地面沉降计算中的应用[J]. 水文地质工程地质, 2003, 30(2): 28-31. DOI:10.3969/j.issn.1000-3665.2003.02.007]
[8]
Luo Zujiang,Zeng Feng,Li Ying. Study on three-dimensional full coupling model of groundwater exploitation and land-subsidence control[J]. Journal of Jilin University(Earth Science Edition), 2009, 39(6): 1080-1086. [骆祖江,曾峰,李颖. 地下水开采与地面沉降控制三维全耦合模型研究[J]. 吉林大学学报(地球科学报), 2009, 39(6): 1080-1086.]
[9]
Shen S L,Xu Y S. Numerical evaluation of land subsidence induced by groundwater pumping in Shanghai[J]. Canadian Geotechnical Journal, 2011, 48(9): 1378-1392. DOI:10.1139/t11-049
[10]
Li Qinfen,Fang Zheng,Wang Hanmei. A mathematieal model and forecast of groundwater workable reserves for Shanghai[J]. Shanghai Geology, 2000, 21(2): 36-43. [李勤奋,方正,王寒梅. 上海市地下水可开采量模型计算及预测[J]. 上海地质, 2000, 21(2): 36-43.]
[11]
Xu Yeshuang,Shen Shuilong,Tang Cuiping,et al. Three-dimensional analysis of land subsidence based on groundwater flow model[J]. Rock and Soil Mechanics, 2005, 26(Supp1): 109-112. [许烨霜,沈水龙,唐翠萍,等. 基于地下水渗流方程的三维地面沉降模型[J]. 岩土力学, 2005, 26(增刊1): 109-112.]
[12]
Luo Guanyong,Pan Hong,Cao Hong,et al. Analysis of settlements caused by decompression of confined water[J]. Rock and Soil Mechanics, 2004, 25(11): 196-200. [骆冠勇,潘泓,曹洪,等. 承压水减压引起的沉降分析[J]. 岩土力学, 2004, 25(11): 196-200.]
[13]
Xie Kanghe,Tao Liwei,Wang Yulin,et al. One-dimensional consolidation solution and analysis for aquitard in leakage system[J]. Journal of Shenyang University of Technology, 2012, 34(5): 581-585. [谢康和,陶立为,王玉林,等. 越流系统中弱透水层的一维固结解及分析[J]. 沈阳工业大学学报, 2012, 34(5): 581-585.]
[14]
Xie Hailan,Wu Qiang,Zhao Zengmin,et al. Consolidation computation of aquitard considering non-Darcy flow[J]. Rock and Soil Mechanics, 2007, 28(5): 1061-1065. [谢海澜,武强,赵增敏,等. 考虑非达西流的弱透水层固结计算[J]. 岩土力学, 2007, 28(5): 1061-1065. DOI:10.3969/j.issn.1000-7598.2007.05.040]
[15]
Zhang Weipeng,Xie Kanghe,Lyu Wenxiao,et al. Analytical solution to one-dimensional consolidation of saturated soil with threshold gradient induced by groundwater pumping[J]. Journal of Central South University(Science and Technology), 2016, 47(3): 875-881. [张玮鹏,谢康和,吕文晓,等. 抽水引起的有起始比降饱和土固结解析解[J]. 中南大学学报(自然科学版), 2016, 47(3): 875-881. DOI:10.11817/j.issn.1672-7207.2016.03.021]
[16]
Xia Changqing.Studies on the nonlinear large strain consolidation theory for lavered structured and rheological soft soils[D].Hangzhou:Zhejiang University,2018.
夏长青.成层结构性软土非线性大应变流变固结理论研究[D].杭州:浙江大学,2018.
[17]
Yan Linghan.Analytical and experimental study on one-dimensional rheological consolidation of soft induced by groundwater level drop[D].Hangzhou:Zhejiang University,2017.
晏凌晗.水位下降引发的软土层一维流变固结解析与试验研究[D].杭州:浙江大学,2017.
[18]
Li J. A nonlinear elastic solution for 1–D subsidence due to aquifer storage and recovery applications[J]. Hydrogeology Journal, 2003, 11(6): 646-658. DOI:10.1007/s10040-003-0283-3
[19]
Zhang Yun,Xue Yuqun. Sensitivity of parameters in one-dimensional nonlinear model for land subsidence[J]. Hydrogeology and Engineering Geology, 2005, 32(3): 1-4. [张云,薛禹群. 一维非线性地面沉降模型参数敏感性分析[J]. 水文地质工程地质, 2005, 32(3): 1-4. DOI:10.3969/j.issn.1000-3665.2005.03.001]
[20]
Huang Dazhong.Studies on theory of coupled consolidation for soil laver induced by groundwater level variation[D].Hangzhou:Zhejiang University,2014.
黄大中.水位变化引发的土层耦合固结变形理论研究[D].杭州:浙江大学,2014.
[21]
Li Chuanxun,Wang Su. An analytical solution for one-dimensional nonlinear consolidation of soft soil[J]. Rock and Soil Mechanics, 2018, 39(10): 3548-3554. [李传勋,王素. 软土一维非线性固结近似解析解[J]. 岩土力学, 2018, 39(10): 3548-3554.]
[22]
Shang Chuanchuan.The research on land subsidence caused by engineering dewatering[D].Beijing:China University of Geosciences,2010.
尚川川.工程降水引发的地面沉降研究[D].北京:中国地质大学,2010.
[23]
Zheng Gang,Zeng Chaofeng,Xue Xiuli. Settlement mechanism of soils induced by local pressure-relief of confined aquifer and parameter analysis[J]. Chinese Journal of Geotechnical Engineering, 2014, 36(5): 802-817. [郑刚,曾超峰,薛秀丽. 承压含水层局部降压引起土体沉降机理及参数分析[J]. 岩土工程学报, 2014, 36(5): 802-817. DOI:10.11779/CJGE201405002]