基于非饱和土理论的黄土地区引水渠道的渗流场研究

王 鲁 琦

(华北水利水电大学 岩土工程与水工结构研究院， 河南 郑州 450045)



基于非饱和土理论的黄土地区引水渠道的渗流场研究

王 鲁 琦

(华北水利水电大学 岩土工程与水工结构研究院， 河南 郑州 450045)

摘要：针对黄土地区引水渠道的渗漏问题，提出了两种防渗方案衬砌层加灰土垫层(2 m厚)和衬砌层加灰土垫层(2 m厚)加渠底原土翻夯(3 m厚)，并且利用基于有限元法的非饱和土水模拟软件SEEP/W，对不同厚度黄土区段的引水渠道进行了数值模拟。通过对比实测数据与室内试验数据，对模型中风积黄土的渗透系数函数及土水特征曲线进行了校正。之后，耦合了应力场及位移场，对引水渠道进行全面多工况的分析发现，随着风积黄土层厚度的增大，等水位线密集区逐渐从衬砌层延伸到灰土垫层；随着时间的增加，风积黄土层厚度越小，等水位线越密集。

关键词：SEEP/W软件；黄土湿陷性；非饱和土；多场耦合；渗流场分析

黄土地区引水渠道的渗漏及湿陷破坏是急需解决的问题，其湿陷变形的结果使原土体更加密实，孔隙比减小，干密度增大 。目前，水利工程中为了防止黄土地基浸水后发生破坏性失陷，一般通过混凝土衬砌、土工膜防渗等防渗措施对湿陷性黄土地基进行处理。其实，采用强夯法、水泥垫层或灰土垫层等处理方法既能消除黄土湿陷性，也能起到降低地基土层渗透性的作用[1-5]。黄土地区引水渠道在投入使用之后，风积黄土层饱和会产生渗流通道，不仅引水渠道的渗漏问题会趋于严重，土体的强度也会降低，进而影响处在斜坡上引水渠段的稳定性[6]。考虑到工区的引水渠道地基黄土湿陷性特点，在进行衬砌防渗的基础上，提出灰土垫层及渠底原土翻夯的处理方案。

针对引水工程中所涉及到不同厚度的引水渠道剖面，拟基于非饱和土理论对提出的两种防渗方案进行计算分析，研究多工况下湿陷性黄土的渗流场特性。由于非饱和理论中涉及的土水渗流计算过程过于复杂，并且大量计算后所得出的结果精确度不高[7-9]，SEEP/W软件基于有限元理论可对多孔渗水材料进行稳态问题及时变问题的分析，且通过瞬态分析可得出不同时刻不同点的孔隙水压力分布状况[10]，因而本文采用了数值软件SEEP/W对引水渠道的渗流场进行了多工况的分析。

1基本理论

1.1渗流计算的基本理论

H. Darcy通过在垂直圆管中装砂进行的渗透试验得到：

(1)

式中：v为截面A上的平均流速；J为水力坡降；k为渗透系数；h为测压管水头。

达西定律的使用有一定的限制，当支配层流的黏阻力在惯性力增加的作用下失去主导作用时不再适用；当黏土中为微小流速的渗流时也不再适用。该工程中的引水渠道所涉及的渗流属于层流范围，达西定律适用。

1.2渗流的基本微分方程及有限单元解法

在渗流场中针对某一单元体，根据质量守恒定律所设定的前提可知，水体质量在某一单元体内积累的速率和进出该单元体内水体质量随时间变化的速率相等，可得渗流的水流运动连续性方程为：

(2)

式中：S为贮水率(贮水率表示当水头h变化一个单位时，因水头变化所引起的骨架变形和水自身性质所导致的膨胀压缩，使得单位体积含水层中释放或贮存的水量)，t为时间。

符合达西定律的各向异性连续体中的稳定渗流，在只考虑二维渗流的情况下，用微分方程表示为：

(3)

式中：h为水头，kx、ky为x、y主轴方向的渗透系数。

SEEP/W软件通过Richards方程对模型中非饱和流进行控制，该非饱和渗流的基本微分方程为：

(4)

式中：H为总水头；k(Θ)x、k(Θ)y、k(Θ)z分别为x轴、y轴、z轴方向的水力传导率；Q为流量边界；Θ为土壤体积含水率；t为时间[11]。

有限单元法是通过用有限个单元的集合体来模拟具有连续性质的渗流场。水头函数在单元内部和单元公共边上均具有连续性，对整个渗流场而言，离散后的单元共同构成了渗流场的分片插值函数，该函数仍是连续函数。整个渗流场的水头分布是由各节点的水头值决定的[12]。

1.3SEEP/W中的土水特征曲线的数学模型

土水特征曲线能够反映土壤基质吸力和含水率之间的函数关系。在非饱和土力学中，通过土水特征曲线可以得到非饱和土的渗透系数、强度以及孔隙水的分布变化规律。SEEP/W模块中的土水特征曲线主要使用以下4种基本数学模型[13]：

①AryaandParis模型

(5)

式中：Ψi为单位土体土壤吸力；T为水分的表面张力；β为接触角，取0；ρw为水的密度；ri为单位土体的孔隙半径。

②Kovacs修正模型

(6)

式中：Θw为体积含水率；n为孔隙率；h∞为毛细水上升高度；Ψ为土壤吸力；C4为校正系数；a、m和Ψn为与材料有关的参数。

③FredlundandXing模型

(7)

式中：其余参数见上述模型，a、m、n为拟合参数。

④VanGenuchten模型

(8)

式中参数见上述模型。

需要通过试验获得Θw和Ψ，结合SEEP/W所提供的数学模型进行数据拟合，可以得到符合实际工况的土水特征曲线。

2引水渠道非饱和渗流场模拟

2.1建立模型

本文讨论的引水渠道所处区段经过风积黄土区域，针对该区差异性厚度风积黄土对引水渠道渗漏问题进行研究，特选取了四个典型剖面进行分析。剖面一：风积黄土层厚度为2 m，下覆卵砾石层；剖面二：风积黄土层7 m，下覆卵砾石层；剖面三：风积黄土层13 m，下覆卵砾石层；剖面四：下覆均为风积黄土层。

四个剖面处的引水渠道均进行了初步的防渗处理，即在渠道内采用现浇混凝土板10 cm，两布一膜200 g/0.3 mm/200 g，以及50 cm的砂砾石换填层的处理方法做成衬砌层。此外，渠道内通水后保持水深3.75 m，地下水位距离地表20 m。将引水渠道处的防渗处理概化为模型，见图1。

图1渠道横断面布置图

2.2参数校正

在进行模拟渗流场之前，针对风积黄土层进行了参数校正。由于现场工况的黄土与试验室黄土的应力状态不同，需对测得的黄土参数进行适当的校正，使之接近甚至吻合现场工况的实测数据，从而保证渗流场数值模拟的可靠性。

在暂态下，可以通过人为调整模型中的参数，使得所需关键参数与实际测定数据相似[14]。风积黄土层其饱和导水率初定为5.1×10-7m/s。土水特征曲线能够反映土体的持水能力，即保有在土壤中的依靠基质吸力而克服重力的水量的多少；非饱和渗透系数反映的是在非饱和区域土体所能导水的能力。非饱和渗透系数以及土水特征曲线对于模型校正具有至关重要的作用，是基于非饱和土理论下表征土体非饱和渗流特性的最重要的两个参数[15]。在距离衬砌层底中心线下方0.5 m、1 m、1.5 m处安置TDR自动土壤含水率测定仪，可测定渠底下的30 d内的含水率数据，由此可得到实测数据源。需对有限元模型中TDR位置所对应的节点进行数据采集，因而模型进行剖分网格时，设定网格剖分的尺寸为0.5 m，针对衬砌层处的网格剖分，主要是沿着衬砌层平行剖分，对于在衬砌层下方的网格有一定程度的细化。之后，结合试验室给出的非饱和水力参数，进行参数校正[16]，本文对此不作赘述。在斜率做出细微调整之后，得到修正后风积黄土层的非饱和渗透系数函数曲线和土水特征曲线(图2)。

图2风积黄土层渗透系数、土水特征曲线

针对剖面的初始衬砌处理工况，对所建立的四个剖面的有限元数值模型进行了暂态分析，分析了引水渠道段在30 d内的渗流场变化。在调整材料的非饱和水力参数后，经过模拟值和实测值的对比，模型与实际所测含水率的拟合较好。在此前提下，数值模拟了渠道瞬态总水头、Y向位移以及总应力，作为参照基础。

2.3防渗处理措施

根据设计方案，针对加强引水渠道的防渗处理措施，共包括两种防渗处理方案(见图3)：

(1) 衬砌层+灰土垫层(2 m厚)

(2) 衬砌层+灰土垫层(2 m厚)+渠底原土翻夯(3 m厚)

灰土垫层处理(图3(a))：此方法是在渠道按照设计开挖轮廓线开挖完成后，不立即通水进行预浸水处理方法，而是对浅部2 m深度内的土体进行换填处理，施工2 m厚度的灰土垫层，在灰土垫层施工完成后再通水运行，并保持渠内水深3.75 m，采用有限元的计算方法研究分析渠道通过此种处理方法后，渠道地基的渗流场、位移场和应力场的分布规律。

灰土垫层+原土翻夯处理(图3(b))：此方法是在灰土垫层处理方法上的一种加强处理手段，即在渠道按照设计开挖轮廓线开挖完成后，对超出设计断面外3 m厚度的土体进行翻夯处理，再施做2 m厚度的灰土垫层，在灰土垫层施工完成后再通水运行，并保持渠内水深3.75 m，针对最易产生渗漏问题的渠底进行了原土翻夯的处理，采用有限元的计算方法分析运行阶段渠道的岩土体的渗流场、位移场、应力场的分布规律。

通过SEEP/W及SIGMA/W模块，针对两种方案进行了不同工况的计算。

(1) 针对衬砌层+灰土垫层： 正常工作；衬砌层失效(给定衬砌层一个大的渗透系数)+垫层正常：计算1 d、10 d、30 d、90 d的渗流场和应力场。

(2) 针对防渗层+原土垫层： 正常工作；衬砌层失效(给定衬砌层一个大的渗透系数)+垫层正常+原土翻夯正常：计算1 d、10 d、30 d、90 d的渗流场和应力场。

进行渗流场模拟计算中所涉及的引水渠道区段岩土体物理力学参数见表1。

图3　防渗处理示意图

3引水渠道数值模拟计算结果分析

3.1渗流场分析

预浸水处理后，由于剖面1、剖面2、剖面3均存在卵砾石层，会形成强渗流通道；且渠底均存在风积黄土层，考虑其湿陷性，施做2 m厚度的具有良好隔水性能的灰土垫层后渠内通水，仅有微量的水透过灰土垫层向其它土层渗漏。随着风积黄土层厚度的增加，断面的渗流场改变愈加显著，通过稳定流和非稳定流以及衬砌层破坏后时间变化的总水头云图进行分析比较可知，施做灰土垫层有利于改善渠堤的渗流影响。

采用2 m厚的灰土垫层+原土翻夯处理后，与仅用灰土垫层进行加固的方案相比，该种方案下的原土翻夯对渗流有一定的阻挡作用，除剖面3外，其它断面在90 d衬砌破坏的情况下，使得可见水位线稳定在原土翻夯层范围内。

由于卵砾石层与风积黄土层的厚度不同，衬砌层破坏后的总水头云图不同，随着风积黄土层厚度的增大，等水位线密集区逐渐从衬砌层延伸到灰土垫层；随着时间的增加，风积黄土层厚度越小，等水位线越密集。特选取图4灰土垫层+原土翻夯方案下渠道防渗层失效后90 d总水头云图作为示例。

3.2位移场及应力场分析

位移场：分析三种处理方法渠底最大沉降可知，预浸水处理后，渠堤的整体沉降变形均较大，渠底处得最大位移达到33 mm，最小位移也有13 mm，而施做灰土垫层或者灰土垫层+原土翻夯处理对整体的沉降并没有明显的的抑制作用，渠底处位移基本不变。

图4灰土垫层+原土翻夯方案下渠道防渗层失效后90 d总水头云图

应力场：对比三种处理方法渠堤总应力分布特点可知，直接采取预浸水处理时，当渠道直接作用在卵砾石层或者直接作用于黄土层时，渠道底部和渠坡中部均会出现不同程度的应力集中现象，这种应力集中作用可能引起渠坡的破坏变形，应该引起足够的重视。当采用灰土垫层方法对渠道处理时，渠道底部和渠坡中部仍会出现不同程度的应力集中现象，但由于灰土垫层的强度较原始土层大，起到了一定的支撑作用，使得应力集中的范围缩小，应力集中数值也较不处理时有所下降，总应力随着深度的增加均匀增大，这对渠堤整体受力是有利的。当采用灰土垫层+原土翻夯的方法处理渠道后，应力的分布及数值变化较仅考虑灰土垫层处理时变化不大，说明增加原土翻夯对总应力的改善效果有限。

4结论

(1) 对于黄土地区的引水渠道，通过灰土垫层处理，对渠堤整体受力有利；在灰土垫层基础上再经原土翻夯处理后，对渠道总应力及地基沉降的改善效果不佳，但是能够有效防治渠道水流的渗漏。

(2) 衬砌层破坏后，随着风积黄土层厚度的增大，等水位线密集区逐渐从衬砌层延伸到灰土垫层；随着时间的增加，风积黄土层厚度越小，等水位线越密集。

参考文献：

[1]孙萍萍，张茂省，朱立峰，等.黄土湿陷典型案例及相关问题[J].地质通报，2013,32(6):847-851.

[2]杨校辉，黄雪峰，朱彦鹏，等.大厚度自重湿陷性黄土地基处理深度和湿陷性评价试验研究[J].岩石力学与工程学报，2014，33(5)：1063-1074.

[3]郭斌.湿陷性黄土地区灰土垫层与防渗膜复合地基研究[D].兰州：兰州大学，2014.

[4]马军，王成，李少华.复合土工膜和黄土组合防渗技术在头屯河水库围堰工程中的应用[J].水力发电，2013，39(9)：52-55．

[5]刘兴邦，王文凯，吴裕勇.湖北阳武干渠险段渠道渗漏分析与处理[J].人民长江，2014，45(增2):163-165.

[6]刘博，李江腾，王泽伟，等.非饱和土渗流特性对库岸边坡稳定性的影响[J].中南大学学报:自然科学版，2014，45(2)：515-520.

[7]吴礼舟，黄润秋.非饱和土渗流及其参数影响的数值分析[J].水文地质工程地质，2011，38(1)：94-98.

[8]陈正汉.非饱和土与特殊土力学的基本理论研究[J].岩土工程学报，2014，36(2)：201-272.

[9]李大美，杨小亭.水力学[M]．2版.武汉：武汉大学出版社，2015.

[10]Geo-Slope Int.. Seepage modeling with SEEP/W 2012: an engineering methodology[M]. Geo-Slope, Calgary, Canada，2012.

[11]陈丽刚.基于ABAQUS渗流与应力耦合作用的边坡稳定性分析[D]．郑州：郑州大学，2010.

[12]刘汉东，姜彤，刘海宁，等.岩土工程数值计算方法[M]．郑州：黄河水利出版社，2011．

[13]周奇.SEEP/W在基于非饱和土理论的公路排水设计计算中的应用[J].中外公路，2009，29(6)：42-46.

[14]Ayse Peksezer Sayit Hasan Yazicigil. Assessment of artificial aquifer recharge potential in the Kucuk Menderes River Basin,Turkey[J]. Hydrogeology Journal，2012,20(4):755-766.

[15]张雪东.土水特征曲线及其在非饱和土力学中应用的基本问题研究[D].北京：北京交通大学，2010.

[16]许淑珍，白晓红，马富丽.利用MATLAB拟合压实黄土土水特征曲线的研究[J].太原理工大学学报，2015，46(1)：81-84.

Seepage Filed of a Diversion Channel in Loess Areas Based on the Unsaturated Soil Theory

WANG Luqi

(InstituteofGeotechnicalEngineeringandHydraulicStructureEngineering,NorthChinaUniversityofWaterResourcesandElectricPower,Zhengzhou,He'nan450045,China)

Abstract：To solve the leakage problem of the diversion channels in the loess area,two designs were discussed. Plan 1: lining layer plus spodosol cushion layer(2 meters thick); plan 2: lining layer plus spodosol cushion layer(2 meters thick) plus dynamic compaction soil of the channel bottom (3 meters thick). And then a numerical simulation of different channels in loess area of different depth was developed by using the finite element software SEEP/W. By comparing the measured data and the laboratory test data, the soil-water characteristic curve and the permeability coefficient function of the eolian loess in the model were calibrated. Based on the calibrated model, the stress field and the displacement field were coupled. According to the overall analysis of the channel, the water isoline concentration area extends gradually from lining layer to the spodosol cushion layer along with the increase of thickness of eolian loess; with the increase of time duration, the smaller the thickness of the eolian losses, the denser the water isoline concentration.

Keywords:SEEP/W；loess collapsibility; unsaturated soil; multifield coupling; seepage field analysis

文章编号：1672—1144(2016)01—0192—05

中图分类号：TV139.14

文献标识码：A

作者简介：王鲁琦(1993—)，男，山东荷泽人，硕士研究生，研究方向为工程地质及水文地质。E-mail:377030692@qq.com

收稿日期：2015-08-04修稿日期：2015-09-05

DOI:10.3969/j.issn.1672-1144.2016.01.036