何晓莹 林潮宁,2 杜效鹄 袁 力 周兴波 刘晓青

(1.河海大学 水利水电学院,南京 210098;2.河海大学 土木与交通学院,南京 210098;3.水电水利规划设计总院,北京 100120)

基岩内的渗流场和应力场耦合是一个重要课题,该课题被研究人员广泛关注并取得了一些重要研究成果.如黄耀英等[1-3]提出了重力坝应力-渗流的无单元耦合分析,研究了渗流场和应力场耦合对混凝土坝位移的影响等;梁坤[4]探讨了Flac3D 流固耦合计算时初始渗流场的设置方式等.但在目前重力坝的安全分析中,一般不考虑坝基渗流场作用,设计规范对坝基整体渗流场的计算也没有明确规定,如对某水库浆砌石重力坝进行稳定分析和对覆盖层地基重力坝应力和变形影响规律分析时均未考虑渗流场影响[5-7].基于此,本文初步设想在重力坝的计算分析中对坝基渗流作用加以考虑.

弹性模量在8～20 GPa之间、Ⅲ类以上的岩体,通常被认为是力学性质较好的岩石[8].坝基为较好岩石时,渗流场对重力坝坝体位移的影响较小,可以不予考虑.但当基岩不满足上述条件时被界定为软岩,渗流作用将对混凝土重力坝的安全运行产生影响,然而如何影响及其影响程度却鲜有研究.

此外,坝基的力学参数往往随着大坝运行而出现时变的特征,因此设计时采用的力学参数值通常与实际值存在差异.为了反映坝基力学性能的长期变化规律,常采用有限元方法,结合原型观测资料反演推求坝基的不同时期的弹性模量,进而预测大坝运行状态,保证大坝安全.近年来,许多学者将各种智能算法引入混凝土重力坝反演分析中,并取得了一系列的研究成果[9].其中,灰狼算法作为一种新兴的群体智能优化算法,因为能够通过模拟灰狼群体的捕食行为、表达狼群的协作机制而达到优化目的,具有结构简单高效且收敛快的特点,因而被用于混凝土重力坝力学参数的反演分析中[10].同时为了提高反演结果的准确性,也应考虑反演对象自身的工程特性和工作环境,如张向韬等针对高心墙堆石坝发展了考虑各向异性与颗粒破碎的堆石料弹塑性本构模型,并将其运用到反演中[11];陈益峰运用有限元正交数值试验、回归分析与优化结合的方法,探讨了复杂地基条件下坝基力学参数反演的问题[12].针对软岩地基上的重力坝,宋海印探讨了在软岩地基上修建重力坝的可能性和建议[13].本文认为在进行重力坝坝基力学参数反演时,应对软岩地基的特点加以重视.

基于渗流场对混凝土重力坝安全影响研究较少的现状,本文在文献[1-3]的基础上,首先在有限元计算中考虑渗流场作用,通过改变坝基弹性模量,分析不同弹性模量下渗流场对重力坝位移的影响规律,进而研究重力坝安全评价时是否需要考虑渗流场作用的判断条件.针对因不考虑渗流场影响导致对软岩坝基上混凝土坝的反演结果出现与工程实际情况不符的问题,进一步提出将智能算法与有限元法相结合反演渗流场作用下重力坝坝基力学参数的方法,探索分析渗流场对重力坝坝基力学参数反演结果的影响问题.本文依据实测位移资料,考虑渗流场与应力场的耦合,通过产生于自然仿生学并行式灰狼算法[14-15]反演坝基的弹性模量,研究渗流场对坝基力学参数的影响,并且给出进行软基上重力坝稳定分析时需要考虑渗流场作用的条件,研究成果对科学合理地评价大坝的安全性具有参考意义.

1 基本原理

1.1 有限元两场耦合计算基本原理

等效连续介质中的渗流层流运动满足达西定律:

式中:v为渗透流速;k为渗透系数;H=H(x,y)为渗流点的测压管水头.

考虑介质和水体的压缩性,结合流体在介质中流动的连续性方程,将达西定律代入可以得到需要考虑压缩性的非稳定渗流的微分方程式:

式中:k x、k y分别是x、y向的渗透系数;ρ为水的 密度;g为重力加速度;α为岩基的压缩系数;β为水的压缩性;n为岩基的孔隙率.

根据公式(1)和公式(2)计算坝基渗流场.坝基内的应力场和渗流场相互作用,能使双场耦合并形成一种平衡状态,逐渐形成稳定的渗流场和应力场.本文在两场迭代分析时,利用二维有限元模型,先通过渗流场计算求得各单元节点的水头值,进而求出相应的渗透体积力和等效节点力,再将它叠加到与初始应力相应的荷载项,求得应力值.

两场耦合的计算步骤为:将渗流场的计算结果,即每个节点的渗透压力作为初始值,转化为单元等效节点外荷载,按照弹性分析求出位移,在荷载作用下重复以上步骤,直至所有单元收敛至精度要求,再施加荷载增量,直至所有的荷载增量计算完成.

1.2 灰狼算法基本原理

灰狼算法是一种新兴的群体智能优化算法,基于自然界灰狼种群觅食的仿生学[10],即通过模拟狼群的捕食行为达到优化的目的,结构简单,需要输入的参数少,且收敛快容易实现[16],其中存在能够自适应调整的收敛因子,可实现局部寻优与全局搜索的平衡.

常规灰狼算法反演坝基力学参数的过程是串行计算,同一时间计算机只能进行一次有限元计算,耗费时间较长.为了提高计算效率,利用并行计算的原理,使用并行式的灰狼算法对坝基弹性模量反演[8],流程如图1所示.基于重力坝的实测位移资料建立统计回归模型,将实测值中的水压分量分离出来.把有限元计算位移与观测点实测位移的水压分量的残差加权平方和,作为参数反演优化问题的目标函数(适应度),寻求较接近实际的坝基弹性模量.

2 工程应用

2.1 工程概况

某混凝土重力坝工程的基岩岩性较为软弱,坝基变形是大坝结构整体变形的主要组成部分,故选取4号坝段为研究对象,其坝顶高程153 m,坝底高程41 m,最大坝高112 m.坝顶宽度6 m,坝体上游面84 m高程以上为竖直面,84 m 高程以下坡度为1∶0.3;下游面坡度为1∶0.75,折坡点高程为145 m.该工程的平面布置图如图2所示.

图2 混凝土重力坝平面布置图(单位:m)

如图3所示,该坝段设置4个测点:倒垂线测点IP9和正倒垂线测点PL1-1、PL1-2、PL1-3,安装高程分别为43、88、120和153 m,倒垂线IP9深入建基面以下40 m(A 点).

图3 混凝土 重力坝剖面及测点布置图(单位:m)

2.2 有限元模型

利用hypermesh 建立二维模型,计算范围选取一个典型重力坝坝段和一定范围内的坝基(上游侧与下游侧均延伸200 m,向下延伸200 m).采用四边形网格对模型实体单元划分,模型的节点851个,实体单元778个.采用直角坐标系,以坝轴线为原点,取顺河向为x轴正向,从上游指向下游;y轴铅直向上,从地层指向地表.有限元模型的分区及网格划分情况如图4所示.

图4 混凝土重力坝有限元计算模型图

针对图4所示模型,进行二维渗流场计算.基岩的上下游侧面及底面均为不透水边界;计算大坝渗流场与应力场耦合时,假定上下游侧面约束为x向连杆约束,底部为y向连杆约束.帷幕的渗透系数为6×10-6cm·s-1,坝基渗透系数为1×10-4cm·s-1.

2.3 有限元计算

本文从有限元正分析的角度,采用最直观的坝体位移指标,在是否考虑渗流场作用的两种情况下,设置多组坝基弹性模量,分析坝基弹性模量变化及渗流场作用对坝体位移造成的影响.

2.3.1 方案设计

首先选取坝体弹性模量25 GPa,坝基弹性模量为6 GPa,并将6 GPa视为基本工况,假设其变化幅度为2 GPa,设置5种研究方案,具体见表1.

表1 坝体材料参数表 (单位:GPa)

2.3.2 计算结果分析

选取43 m 高程测点处的坝体位移进行分析比较.将2018年6月21日的位移值视为基点,x轴为顺河向,顺河向位移指向下游为正,指向上游为负.

图5为43 m 高程处考虑和不考虑渗流场坝体x向(顺河向)位移随时间变化情况.基岩内的水头随时间而变,大坝在水荷载的作用下,顺河向位移随时间而变化,位移受到水位变化的影响.

图5 43 m 高程坝体顺河向位移图

从图5可以看出:位移与水位呈正相关关系,符合重力坝的变形规律;随着坝基弹性模量减小,坝体位移及变幅增大;在同样的坝基条件下,考虑渗流场作用时坝体顺河向的位移比不考虑时的位移明显增大,说明渗流场对坝体位移有影响.从图5来看,在5种方案中,当坝基的弹性模量从4 GPa变化到2 GPa且水位变化时,顺河向位移变化幅度很大.为了更加直观描述这种变化,对几种计算方案进一步细化.选择测点IP9安装高程43 m 的顺河向位移,在已有方案基础上,增加1 GPa和3 GPa两种方案,研究在渗流场作用下,坝基弹性模量变化对位移及其变幅造成的影响.计算结果如图6～7所示.

图6 43 m 高程坝体顺河向位移图(考虑渗流场)

图6和图7分别是考虑渗流场作用和不考虑渗流场作用时,坝体43 m 高程的顺河向位移图,考虑渗流场时的顺河向位移明显大于不考虑渗流场作用时.对水位下降较明显的2020年5月4日到2020年8月4日期间的位移变幅进行分析,见表2和图8.

图7 43 m 高程坝体顺河向位移图(不考虑渗流场)

表2 各坝基弹性模量下顺河向位移变幅 (单位:mm)

图8 坝基弹性模量43 m 高程顺河向位移变幅图

从图8 和表2 可以看出,随着坝基弹性模量减小,顺河向位移的变幅增大.当坝基弹性模量小于3 GPa时,曲线变陡,也就是顺河向位移变幅随着坝基弹性模量的减小而增大的速度加快,并且考虑渗流场的作用时这种趋势更加明显.

当坝基弹性模量处于6～10 GPa时,坝基弹性模量每下降1 GPa,位移变幅平均增长了9.1%;当坝基弹性模量处于3～6 GPa时,坝基弹性模量每下降1 GPa,位移变幅平均增长了17.2%;当坝基弹性模量小于3 GPa时,坝基弹性模量每下降1 GPa,位移变幅平均增长了38.0%.这说明当坝基弹性模量较小时,是否考虑渗流场作用将会造成较大的差异.因此对于本文针对的混凝土重力坝及其坝基情况所做的二维分析可知,当水位变化较大,坝基的弹性模量处于3～6 GPa时,建议考虑渗流场对坝基造成的影响;当坝基的弹性模量小于3 GPa时,应考虑渗流场对坝基造成的影响.

3 坝基弹性模量反演

第2节利用有限元计算分析了渗流场作用对不同坝基弹性模量下的重力坝位移的影响,为了验证并进一步研究影响规律,利用该工程重力坝的实测位移数据,在坝基反演时考虑渗流场作用,并与不考虑渗流场作用的结果进行对比,分析讨论渗流场对于坝基弹性模量反演结果的影响.本文利用Matlab软件,编写并行式灰狼算法程序[17],子种群数量为8个,最大迭代次数为50 次,收敛判断条件为目标函数小于10.

3.1 方案设计

该工程地基大致分为3层,其中第1层为软岩,第2～3层泥岩类岩层埋深较深,性状较好.为简化计算,分成如图4 所示的上、下两层,上层为浅层坝基(第1层),下层为深层坝基(第2～3层),深层坝基的弹性模量按照工程资料选取,反演浅层坝基的弹性模量.各材料具体参数见表3.

表3 坝体坝基材料参数表(基本情况)

以图3～4为算例进行初步分析.为了验证方法的正确性与可行性,先假定一个坝基弹性模量,通过有限元计算得出位移值,再根据位移值去反演坝基的弹性模量,观察两者是否一致.预设浅层坝基的弹性模量为3 GPa,坝基弹性模量取值范围1～6 GPa,将考虑渗流作用时有限元计算的位移值作为基准值,设置两种反演计算方案:

(1)方案1:考虑渗流场与应力场耦合作用(施加水荷载和自重,水荷载是加载在上游坝面的面荷载,自重为加载在所有结点上的体荷载,并考虑坝基的渗流场作用),反演浅层坝基弹性模量;

(2)方案2:不考虑渗流场的作用(施加水荷载、扬压力和自重,水荷载为加载在上游坝面的面荷载,扬压力为作用在建基面上的节点力,自重为加载在所有结点上的体荷载),反演浅层坝基弹性模量.

方案1反演的数值为3.0 GPa,与预设的结果相同,说明算法适合用于反演混凝土重力坝坝基的弹性模量.对于方案2,在同样的基准值下,不考虑渗流场作用时反演的浅层坝基弹性模量为1.8 GPa,明显小于设计值,说明此时坝基渗流场对重力坝的位移有较为明显的影响.因此,根据算例初步分析若坝基条件较差,不考虑渗流场的作用会导致反演的坝基弹性模量偏小,与工程实际存在较大差异.

3.2 工程实例反演

将算例方法运用到工程实例中.利用工程实测值,即倒垂线测点IP9和正倒垂线测点PL1-1、PL1-2、PL1-3的位移测值进行浅层坝基的弹性模量反演.反演的计算位移值分别为各测点的相对位移,即计算IP9相对于点A 的相对位移、PL1-1相对于IP9的相对位移、PL1-2 相对于PL1-1 的相对位移、PL1-3 相对于PL1-2的相对位移.计算时段为2018年6月21日到2022年3月26日.

为了研究渗流场作用对反演结果造成的影响,依然以3.1中的方案1和方案2为反演计算方案.由于位移监测数据较多,本文仅将接近坝基面的IP9的位移监测数据列于表4,以便进行反演分析和与有限元计算结果的比较.

表4 2018-2022年IP9测值 (单位:mm)

对监测数据进行初步筛选,剔除异常值.反演时把有限元计算位移与观测点实测位移的水压分量的残差加权平方和作为参数反演优化问题的目标函数,寻求接近实际的坝基弹性模量,反演结果见表5.

表5 工程实例反演结果(单位:GPa)

将测点IP9的实测值的水压分量与预测值进行对比,结果如图9所示.

图9 实测值(水压分量)与反演值对比图(IP9)

从图9可以看出,实测值与反演值总体吻合,规律形似,且量级相同,说明反演的效果较好,反演结果可信.从表5可以看出,考虑渗流场时反演的坝基弹性模量明显大于不考虑渗流场时.这是由于当坝基弹性模量较小的时候,渗流场的影响较为明显,因此当坝基弹性模量相同时,考虑渗流场作用时的坝体位移较大.这也说明了为什么对相同的监测值,考虑渗流场作用时反演得到的坝基弹性模量较大.

综上所述,当坝基较软时,不能忽略渗流场的作用,否则反演得到的结果会小于实际的坝基弹性模量.

4 结论

1)基于灰狼算法的坝基弹性模量有限元反演方法对解决工程实际问题是可行有效的.

2)不考虑渗流场作用时,反演得到的重力坝力学参数偏小.当坝基为软岩时,不考虑渗流场的反演结果与实际情况不符,这可能对大坝安全评价造成偏差.

3)针对本文研究的某重力坝工程,当其坝基的弹性模量处于3～6 GPa时,建议考虑渗流场的作用;当坝基的弹性模量小于3 GPa时,应考虑渗流场对坝基造成的影响.本文研究所针对的是二维的简化坝基,对包含复杂坝基的实际工程情况应做更深入的探讨和研究.