基于有限元软件的水库黏土心墙坝渗流及边坡稳定研究

何 明

(新疆昌源水务准东供水有限公司,乌鲁木齐 831700)

0 引 言

水库大坝的稳定性直接影响着水库的正常运行,以及周边居民的生命财产安全。在各类坝型中,黏土心墙坝因其良好的防渗性能和经济效益而被广泛应用[1-2]。但随着工程规模的扩大和建设环境的复杂化,黏土心墙坝的渗流特性和边坡稳定性问题逐渐显现出来。渗流是坝体中的常见现象,但过大的渗流量可能导致坝体内部的水压过高,从而影响坝体的稳定性[3-4]。因此,对渗流特性的研究,特别是在特定工况下的渗流场分析,对于保障坝体的稳定性具有重要意义。

传统的渗流和边坡稳定性分析方法大多基于简化的理论模型或试验数据,难以准确描述实际工程中复杂的地质条件和工况变化[5]。近年来,有限元方法因其在处理复杂问题上的优势而得到广泛应用。因此,本文利用有限元软件,对新疆昌吉某水库黏土心墙坝的渗流及其边坡稳定性进行分析,研究成果可为渗流-应力耦合计算的应用提供参考与借鉴。

1 工程概况

新疆昌吉某水库工程大坝采用传统黏土心墙坝,位于木垒县城南方约2km处,工程主要由4个部分组成,隶属于拦河水库。该拦河大坝的主坝体为砂壳黏土心墙坝,整个坝体高42m,长1 100m,顶部宽5m,水库总库容约1 400×104m3。水库主要功能是为灌溉提供水源,同时具有防洪、养鱼和旅游等多重功能,是一座综合性的中型调节水库。

水库恰好位于木垒河的山口出口,从1958年开工,经过几次续建和改建,至1998年形成目前规模。近些年,当地居民积极投身绿化工作,使水库周围的绿化区域达到8.333hm2,为水库区域创造宜人的生态环境。考虑到该水库主要服务于农业灌溉需求,其库区水位会频繁发生波动,使得黏土心墙坝的土体经常处于交替的饱和与非饱和状态,从而影响土体的稳定性和其抗剪强度。因此,对黏土心墙坝进行坝体渗流与稳定性分析十分必要。通过分析,可以准确识别潜在的风险点,并制定有效的应对策略。针对所存在的问题,实施相应的修复措施是确保水库大坝长期安全的关键步骤。

2 有限元软件在水库黏土心墙坝渗流与边坡稳定性研究中的应用

2.1 水库黏土心墙坝渗流与应力耦合计算

渗透稳定问题是渗流对土工建筑物或地基内土颗粒产生渗透力作用,以及土的稳定性受到渗流破坏的问题。由于渗透力而引起的渗流稳定问题在实际工程中经常发生,是土工建筑物或地基发生破坏从而引发工程事故的重要原因之一。在天然环境中,渗流场与应力场常常相互交织和影响,形成一种紧密的耦合关系。渗流分析的主要方法有水力学法、流网法、差分法以及有限元法等。

在对土坝的渗流和应力进行分析时,如果仅从单一的渗流场或应力场角度出发,其结果与真实场景可能会有所偏差[6-7]。因此,许多学者开始在渗流与应力的计算中融入两者之间的相互关系。由于耦合计算具有一定的复杂性,因此需要利用专门的计算程序来实现二者之间的耦合[4,8]。同时,考虑到渗流场与温度场有一定的相似性,其温度-应力的耦合分析模块可被用来满足渗流-应力的分析需求。

因此,本研究采用ANSYS软件作为基础平台,通过融合APDL命令流语言与Fortran编程语言,开发一套专用于坝体渗流与应力耦合的计算程序。通过对试验选取的黏土心墙坝进行分析,以确认耦合计算的重要性和准确性。所有模型中,邓肯-张(Duncan-Chang)模型是常用的非线性弹性模型之一,主要包括E-v模型和E-B模型。该模型无需考虑土体的塑性变形,简化了迭代计算,确保计算的收敛性。同时,可以准确反映土体的应力-应变特性,且参数较少,可通过试验确定,因此在数值计算中得到广泛应用[9-10]。虽然ANSYS软件目前尚无法直接支持邓肯-张模型,但鉴于该模型计算公式的简明特性,研究利用APDL编程语言进行模型的集成实现。

大坝土体中,渗流与应力的相互耦合主要体现在两个方面:①应力场会修改土体的渗透系数,进而调整渗流场。②渗流场的渗透力也会影响应力场的分布情况,通过二者的相互耦合,实现对大坝渗流的研究。试验通过迭代法(间接耦合法)获得这两种场在耦合下的分布情况,其中渗流场对应力场的影响可以分为静水压力(渗透压力)和动水压力(渗透体积力)[11-12]。同时,渗流力能够被直接作用于坝体与应力的接触面上,在计算时将其作为面力被考虑在节点荷载中。渗透体积力则可以看作与自重类似,在计算时按照体力被计入节点荷载中。基于上述原理,结合APDL和Fortran语言,开发了渗流-应力的耦合迭代计算程序,其流程见图1。

利用图1中的模块化程序设计,工作人员只需要初步掌握有限元知识,便可得到建模效果。为了提高计算效率,研究选取大坝最高断面进行有限元建模,而非三维分析。在考虑流固耦合时,利用邓肯-张E-v模型进行应力计算[13]。试验模拟的渗流边界条件设定如下:大坝基底作为不渗透边界,其左侧界定为上游水位,右侧则为无水位状态,下游边界可能设为溢流边界。关于应力边界条件的设定如下:大坝基底被施加X轴和Y轴方向的位移限制,而大坝上下游两侧则只受X轴方向的位移限制。模型基于笛卡尔直角坐标系构建,其中X轴朝向大坝下游,Y轴则向上。根据大坝最高断面,建立大坝平面有限元模型,在建立模型时,选用四面体网格划分。离散后,通过计算得到的有限元网格见图2。

图1 渗流-应力耦合程序计算流程

图2 坝体断面有限元网格示意图

2.2 融合邓肯-张模型的强度折减法设计

边坡稳定分析是黏土心墙坝和土石坝设计的核心环节,它可以对坝坡是否失稳以及可能存在的破坏形式进行评估,进而为整体结构与设计提供关键的参考。尽管黏土心墙坝边坡稳定的研究不断深入,但大多数研究仍主要关注应力场对边坡稳定性的影响,较少涉及渗流与应力之间的耦合作用,因此相关研究可能存在盲区。

为了更精确地评估边坡稳定性,有些学者使用有限元强度折减法,对应力进行计算。但该模型会由于屈服面的尖锐而引起计算上的复杂性、收敛速度减慢或不收敛等问题。鉴于此,试验结合邓肯-张模型的有限元强度折减法,并基于自主开发的模块化渗流与应力耦合计算程序,对边坡稳定性进行深入分析。将强度储备概念的安全系数作为基础,利用有限元强度折减法,可以得到边坡稳定安全系数。对黏聚力、内摩擦角进行折减,得到新的表达式,公式如下:

式中:Fs为安全系数,即折减系数;φf、cf分别为折减计算后得到的内摩擦角与材料内聚力;φ为内摩擦角;c为黏聚力。

试验引入邓肯-张模型,并通过Mohr-Coulomb定律,推导求得边坡受到破坏时的主要应力差值,然后可以得到新的切线模量。具体计算公式如下:

式中:Ei为初始弹性模量;Rf为破坏比;(σ1-σ3)f为破坏时的主应力差;P为渗透压力。

有些学者的研究指出,仅仅对材料的强度或相关参数进行调整,并不能完全反映强度折减法的深层含义。为了更加充分体现土体的变形特性,不仅需要考虑初始切线模量的折减,同时还需要时刻关注初始切线模量k的调整与折减。这种计算方式提供了一个更为综合的视角,能够更深入地理解强度折减对边坡稳定性性质的影响,得到新的切线模量计算表达,公式如下:

式中:k为初始切线模量。

3 黏土心墙坝渗流及边坡稳定研究

为了对研究提出的有限元分析方法检测时的可靠性与有效性进行验证,研究引入CNPM3D软件,将其与试验所构建的ANSYS软件系统运行下不同坝体料区的最大渗透坡降进行分析,具体结果见表1。

表1 运行期坝体各料区的最大渗透坡降

由表1可知,在不同工况中,坝体黏土心墙和防渗帷幕的渗透坡降显著,而其他部分如反滤层和坝壳风化料的渗透坡降则相对较小。在校核洪水位工况下,各料区的渗透坡降最大,均满足稳定性要求。此外,不考虑流固耦合时,ANSYS软件得到的结果与CNPM3D软件十分近似。而考虑流固耦合时,ANSYS的结果明显低于未考虑耦合的情况。这是因为在应力场的作用下,材料的渗透性增强,导致心墙和帷幕对水头的削弱能力较弱。在渗径不变的情况下,上下游水头的差值会变小,并且可以计算得到相对较小的渗透坡降值。

试验对最不利的工况进行分析,即选取校核洪水位工况,对边坡稳定性的计算展开分析。在考虑与不考虑坝体建造耦合的情况下,坝体Y向位移与Y向应力的分布见图3。

图3 坝体Y向位移等值线分布(单位:m)

由图3可知,当引入坝体渗流与所受到的应力产生耦合效应时,水库坝体在Y方向上的位移会显著增加。这一现象可能源于考虑耦合作用后,渗透水压力在Y垂直方向上的分量开始增加,导致大坝土体呈明显的沉降趋势。因此,将耦合效应考虑在内时,大坝体在Y方向上的位移会使这一现象更为明显。大坝体Y方向上应力等值线分布见图4。

由图4可知,在对耦合与不耦合的情况进行分析时,大坝体在Y方向上的应力均遵循一种规律。这意味着即使坝体受到来自上方荷载的影响,在高程方向上的应力从上到下会开始逐渐增加,其中最大的应力会出现在基岩底部。

图4 大坝体Y方向上应力等值线分布(单位:MPa)

当不考虑耦合的情况时,心墙坝Y方向上的应力会显著低于坝体两侧的应力,呈现出强烈的拱效应。通过分析可知,这一现象可能是因为心墙坝本身的形变模量较邻接坝壳材料小,进而导致心墙坝的沉降量远超坝壳材料,使得心墙的一部分自重转移到坝壳材料上,导致心墙内部所承受的应力小于坝壳材料。

当坝体受到应力的分布考虑耦合效应时,水库心墙坝区域呈现出明显的拱形效应。不同于未考虑耦合作用的情况,当将耦合效应考虑在内时,在相同的高程上,水库坝体上游侧在Y方向上受到的应力会大于下游侧所受到的力。这种差别可能与考虑渗流水压力时的耦合效应有关:上游水位较高,受到渗流水的压力与应力会更大,使得上游侧土体承受的压力超过下游侧的土体。

针对特征点的垂直位移与安全系数之间的关系进行分析,具体结果见图5。

图5 坝体特征点的垂直位移与安全系数的关系

由图5可知,在对校核洪水工况进行分析时,一旦将坝体耦合效应考虑在内,得到的折减系数会从1.2开始增加至1.3,此时特征关键点的位移会开始显著的偏离常态。为了获得更为精准的数据,试验将该区间的折减系数间隔缩小,并将折减系数的间隔逐渐细化,设置为0.025。在细化过程中发现,当坝体的安全系数介于1.225～1.250时,水库大坝会突然出现位移上的跳跃。因此,当坝体处于耦合效应时,水库大坝边坡稳定性的试验安全系数设置为1.225,坝体有更好的稳定性。

另一方面,若忽略坝体所受到的耦合效应,当折减系数的数值在1.30～1.45之间变动时,坝体特征关键点的位移也会发生显著变化。采用同样的0.025系数间隔,进一步精确折减系数,结果发现,当安全系数波动范围在1.300～1.325时,特征点表现出显著位移跳跃。因此,在不考虑耦合效应的情况,安全系数为1.325时,大坝的稳定性最为优越。

但总体来看,考虑耦合效应后,得出的大坝边坡稳定性安全系数相较于忽略渗透应力耦合的情况略有下降。这可能是因为渗透水压力的竖向分力会导致坝体土体有下沉的倾向,进一步增强了坝体边坡土体的下滑趋势。同时也意味着当坝体处于渗流情况,且考虑耦合效应时,渗流会对水库大坝的稳定性产生不利的影响。

4 结 论

本文利用有限元软件,对水库黏土心墙坝的渗流与边坡稳定性进行了分析。结果显示,在不同工况下,黏土心墙坝与防渗帷幕的渗透坡降数值十分显著,而反滤层与坝壳风化料的渗透坡降则较小。考虑耦合情况时,坝体的边坡稳定性安全系数为1.225;不考虑耦合情况时,坝体的边坡稳定性安全系数为1.325。通过对渗流场的准确分析,可以预测边坡可能的滑动面和滑动力,为坝体的设计和施工提供技术支撑。