参数空间变异性下坝基防渗墙地震反应

doi:10.16183/j.cnki.jsjtu.2021.242

参数空间变异性下坝基防渗墙地震反应

张富有, 周强强^,, 杜鹏程

河海大学土木与交通学院,南京 210024

Seismic Response Analysis of Cut-Off Wall of Dam Foundation Under Spatial Variability of Parameters

ZHANG Fuyou, ZHOU Qiangqiang^,, DU Pengcheng

College of Civil Engineering and Transportation, Hohai University, Nanjing 210024, China

通讯作者: 周强强,男,硕士生,电话(Tel.):15700012815;E-mail:qqzhou@hhu.edu.cn.

责任编辑: 王一凡

收稿日期: 2021-07-6

基金资助:

长江水科学研究联合基金重点支持项目(U2040221)

Received: 2021-07-6

作者简介 About authors

张富有(1972-),男,河南省驻马店市人,副教授,硕士生导师,从事工程抗震和隔震减震方面的研究.

摘要

堆石坝覆盖层土体类型复杂,结构松散且有时存在透镜体,使得覆盖层材料参数存在空间变异性,地震作用下对坝基防渗墙应力变形存在不可忽视的影响.在蒙特卡罗框架下,采用Cholesky分解的高斯空间随机场离散方法,利用Python对Abaqus进行二次开发,实现了覆盖层材料空间变异性的“非侵入式”随机有限元分析.研究结果表明:考虑覆盖层材料空间变异性的情况下防渗墙不同部位的应力变形均有不同程度的超标,且动应力的变异系数大于动位移的变异系数,显示覆盖层材料参数的变异性对防渗墙在地震作用下的应力变化更敏感.因此在设计阶段应考虑这种不确定性,以确保坝基防渗墙的安全运行.

关键词： 覆盖层; 防渗墙; 空间变异性; 随机有限元; 地震反应

Abstract

The overburden of rockfill dams has complex soil types, loose structures, and sometimes lenses, which make the material parameters of the overburden have spatial variability, and the impact of earthquakes on the stress and deformation of the dam foundation cut-off wall cannot be ignored. Under the framework of Monte Carlo, by using the Gaussian spatial random field discrete method of Cholesky decomposition, the secondary development of Abaqus with Python is used to realize the “non-invasive” stochastic finite element analysis of the spatial variability of the overburden material. The research results show that considering the spatial variability of the cover material, the stress and deformation of the different parts of the cut-off wall have different degrees of excess, and the coefficient of variation of dynamic stress is greater than that of the dynamic displacement, indicating the variability of the material parameters of the cover is more sensitive to the stress changes of the impervious wall under the action of earthquake. Therefore, this uncertainty should be considered in the design stage to ensure the safe operation of the dam foundation impervious wall.

Keywords： overburden; cut-off wall; spatial variability; stochastic finite element; seismic response

PDF (4705KB) 元数据多维度评价相关文章导出 EndNote| Ris| Bibtex 收藏本文

本文引用格式

张富有, 周强强, 杜鹏程. 参数空间变异性下坝基防渗墙地震反应[J]. 上海交通大学学报, 2022, 56(5): 684-692 doi:10.16183/j.cnki.jsjtu.2021.242

ZHANG Fuyou, ZHOU Qiangqiang, DU Pengcheng. Seismic Response Analysis of Cut-Off Wall of Dam Foundation Under Spatial Variability of Parameters[J]. Journal of Shanghai Jiaotong University, 2022, 56(5): 684-692 doi:10.16183/j.cnki.jsjtu.2021.242

混凝土防渗墙作为土石坝地基覆盖层防渗的关键结构和大坝防渗的重要组成部分,不仅要承受上覆坝体的压力,还要承受库水压力、地基土侧压力、河谷拱效应以及覆盖层与防渗墙之间的不均匀沉降引起的负摩擦阻力等作用^[1⇓-3].虽然防渗墙位于覆盖层内部,受周围约束较大,地震作用下,和上部坝体相比加速度和动位移反应较小;但其本身具有体积较大,自重较重,弹性模量较大,变形协调能力差,抗拉强度低等缺点;当地震发生时,动应力变化比较复杂,叠加自身产生的较大静应力,其承受的荷载有可能超出抗拉强度,出现墙体裂缝^[4].目前,对于覆盖层上坝基混凝土防渗墙应力变形情况的研究大多采用的是静动力条件下基于材料参数确定性的传统有限元建模研究,这种研究通常采用的都是唯一的、确定性的材料参数.然而,由于坝基覆盖层钻孔取样时不可能无限密集取样,取样时也不可避免的会对土体存在扰动;另外,由于碾压施工的不均匀性、雨水、温度等条件的影响,造成坝体及覆盖层材料参数具有变异性,仅采用实验室得到的确定性材料参数进行传统的有限元建模分析可能无法准确地反映工程的实际情况^[5⇓-7].

近年来,对岩土材料参数空间变异性的研究取得了一定成果.文献[8-9]表明忽略土体参数空间变异性将高估边坡失效概率.王建娥等^[10]通过非侵入式随机有限元分析得出材料参数空间变异性对堆石坝变形计算结果影响明显.杨鸽等^[11]进行土石坝地震反应随机有限元分析,表明如果忽略材料参数的不确定性可能导致低估大坝的地震反应.郭晴等^[12]考虑坝体与覆盖层材料参数空间变异性的情况下,表明考虑材料参数不确定性的心墙坝的应力变形存在一定的概率会超过确定性有限元分析.朱晟等^[13]考虑级配空间随机特性对堆石坝应力变形进行分析,得出若不考虑随机性则有一定概率会低估坝体的应力变形.文献[14]研究了地震激励和材料参数的随机性对大坝动响应的影响,表明考虑相关性的随机切变模量引起的响应值要大于确定参数下的结果.

覆盖层具有土体类型及分布的复杂性、松散性等特点,有时还存在透镜体,因此相对于坝体材料参数覆盖层材料参数空间变异性更强,同时覆盖层材料和防渗墙直接接触,会对坝基防渗墙的地震反应分析结果产生影响.因此,本文将随机场理论与有限元分析结合,在蒙特卡罗框架下,采用基于Cholesky分解的高斯空间随机场离散方法,利用Python语言对Abaqus进行二次开发,考虑在覆盖层材料参数的不确定性下,建立了地震作用下堆石坝坝基混凝土防渗墙应力变形分析的非侵入式随机有限元计算分析方法,并对堆石坝坝基混凝土防渗墙的应力变形结果进行统计分析,结果可以为地震作用下堆石坝的安全分析和坝基混凝土防渗墙的设计优化提供参考.

1 堆石坝坝基覆盖层材料参数的随机性及其模拟

1.1 随机有限元计算流程

岩土体材料由于其本身的特点,力学参数具有一定的离散性,分布特征在空间上也具有一定不确定性.对于坝基覆盖层有时还存在透镜体,其材料参数的空间差异性更为显著.因此,为了研究在地震作用下材料参数的空间差异性对堆石坝坝基混凝土防渗墙应力变形的影响,有必要考虑采用随机有限元法进行分析.

随机有限元法是将随机场模拟技术与确定性有限元方法结合起来的一种方法.首先利用Python程序在蒙特卡罗法的框架下基于Cholesky分解法对系统中的变量进行自动抽样,然后利用Python程序直接修改Abaqus的全命令流inp文件,并提交计算求解结构的随机反应,最后对求得的多组结构随机反应样本进行统计分析.随机有限元计算流程如图1所示,该方法实现了随机场模拟的自动化以及随机参数的自动赋值,提高了计算效率.

图1

新窗口打开| 下载原图ZIP| 生成PPT

图1 随机有限元计算流程

Fig.1 Stochastic finite element calculation process

为了判断蒙特卡罗模拟次数是否足够,通过选取防渗墙上特征点,每次模拟结束后输出防渗墙特征点的沉降值,根据下面的公式计算其沉降平均值x_ave,并绘制它随模拟次数的变化曲线,当变化曲线趋于平稳时,认为蒙特卡罗模拟次数足够,其计算结果达到稳定.

(1)x_ave(N')=

\frac{1}{N'} \sum_{i}^{N'}

x_i

式中:x_ave(N')为第N'次模拟时防渗墙特征点的沉降平均值;x_i为第i次模拟的防渗墙特征点的沉降值;N'为模拟次数.

1.2 随机场模拟方法

根据随机场理论,覆盖层中任意两个不同空间位置处土体材料参数的自相关性可以采用自相关函数来描述^[10,15].对于随机场材料参数的统计特性可通过均值μ和标准差σ表征,其变异系数:

(2)C_v=μ/σ

对于随机场中任意两点间的自相关系数定义如下^[16]:

(3)ρ[H(x_i, y_i), H(x_j, y_j)]=

\frac{C o v [H (x_{i}, y_{i}), H (x_{j}, y_{j})]}{\sqrt{V a r [H (x_{i}, y_{i})]} \sqrt{V a r [H (x_{j}, y_{j})]}}

式中:H(x_i, y_i)、H(x_j, y_j)分别为随机场中任意两个单元网格的中心坐标;Cov(*)、Var(*)分别为中心点空间坐标的协方差和方差函数.

目前有多种描述随机场自相关性的相关函数,其中由于高斯型自相关函数平稳性好,所以本文选用高斯型自相关函数来描述各单元材料参数之间的空间相关性^[11].

(4)ρ[(x_i-x_j), (y_i-y_j)]=exp

[- π (\frac{{|x_{i} - x_{j}|}^{2}}{l_{h}^{2}} + \frac{{|y_{i} - y_{j}|}^{2}}{l_{v}^{2}})]

式中: $|x_{i} - x_{j}|$ 为空间两点水平方向的距离; $|y_{i} - y_{j}|$ 为空间两点竖直方向的距离;l_h和l_v分别为水平和竖直自相关距离.

一般而言,覆盖层材料参数均大于0且材料参数之间存在一定的相关性,因此,可以假定材料参数为服从对数正态分布的随机变量^[17].根据单元形心坐标和自相关函数确定随机场,离散后得到的一组相关随机变量的相关系数矩阵如下:

(5)ρ_n_×_n=

[ρ_{11} ρ_{12} \cdot\cdot\cdot ρ_{1 i} \cdot\cdot\cdot ρ_{1 n}]

将ρ_n_×_n进行Cholesky分解成一个上三角矩阵A₁和一个下三角矩阵A₂, 将下三角矩阵A₂左乘一个服从标准正态分布的随机样本B,即可得到相关标准正态分布随机场样本 $H_{i}^{D}$ (x, y),对随机场样本 $H_{i}^{D}$ (x, y)取指数^[18],即可得到所需的材料参数相关对数正态随机场:

(6)H_i(x, y)=exp(

μ_{l n x_{i}}

E+

σ_{l n x_{i}} H_{i}^{D}

(x, y))

(7)

μ_{l n x_{i}}

=ln

μ_{x_{i}}

σ_{l n x_{i}}^{2}

(8)

σ_{l n x_{i}}

\sqrt{l n [1 + {(σ_{x_{i}} / μ_{x_{i}})}^{2}]}

式中: $μ_{l n x_{i}}$ 、 $σ_{l n x_{i}}$ 分别为第i个单元变量取对数后的均值和标准差.

2 实例分析

某黏土心墙堆石坝,坝顶长度247.00 m,宽10.00 m,最大坝高74.50 m.心墙顶宽4.00 m,上下游坡度比为1:0.25,心墙最大底宽39.25 m,心墙与两岸坝肩混凝土垫层接触处采用厚度 2.00 m 的高塑性黏土过渡.覆盖层厚度最深达44.76 m,防渗采用混凝土防渗墙,墙厚1.20 m.防渗墙顶设置廊道,上部为黏土心墙,底部嵌入基岩.坝体分区如图2所示.

图2

新窗口打开| 下载原图ZIP| 生成PPT

图2 大坝典型断面材料分区图(m)

Fig.2 Material zoning map of a typical section of dam (m)

2.1 计算模型

单元网格划分时,大坝共划分 59 480 个节点、56 848 个单元.对土石料、覆盖层和混凝土等单元,大部分为8结点6面体单元,少数用6结点5面体单元和4结点4面体等单元过渡.由于需研究防渗墙和廊道的应力变形情况,所以对它们进行了精细化划分,过渡单元主要设置在防渗墙和廊道网格的附近以及材料分区的变化边界等位置,便于各区域过渡相连.静力边界条件设置为底面三向约束,河谷方向对y向约束;动力边界条件采用黏弹性边界,同时采用接触面单元模拟不同材料边界接触时的相互错动.大坝的静力边界条件和动力边界条件如图3和图4所示.地震波考虑同时从顺河水平向、横河水平向和竖向三向输入.根据《水工建筑物抗震设计规范》^[19],竖向加速度分量取水平向加速度分量的 2/3,即水平方向最大加速度为0.15g,竖直方向最大加速度为0.10g,地震波从基岩顶部输入.为了节约计算时间,本文截取0~15.54 s的地震波,输入间隔0.02 s.地震波加速度时程曲线(a-t)如图5所示.图中gal用来度量地震加速度,1 gal=1 cm/s².

图3

新窗口打开| 下载原图ZIP| 生成PPT

图3 大坝静力边界条件

Fig.3 Dam static boundary conditions

图4

新窗口打开| 下载原图ZIP| 生成PPT

图4 大坝动力边界条件

Fig.4 Dam dynamic boundary conditions

图5

新窗口打开| 下载原图ZIP| 生成PPT

图5 地震波加速度时程曲线

Fig.5 Seismic acceleration time history curve

2.2 随机材料参数的统计特性

对于静力计算,本文采用邓肯-张E-B非线性弹性模型.对于邓肯-张E-B非线性弹性模型参数敏感性的研究较为充分,尽管众多学者关心的敏感性指标不同,但根据文献[20 ⇓-22]对邓肯-张E-B模型参数的敏感性分析发现,土石坝静力参数密度ρ、弹性模量K、内摩擦角φ₀、弹性模量指数n、破坏比R_f对结果的影响较为显著,这5个参数中ρ、K和φ₀影响要比n和R_f的影响大很多.对于算例中的坝体,坝基混合卵石层较厚,与防渗墙直接接触的仅这一层覆盖层土体,故综合考虑选取混合土卵石层覆盖层参数ρ、K和φ₀为随机参数.一般情况下常用的随机变量的分布方式主要有正态分布、对数正态分布、极 III 型分布等,根据陈辉等^[23]对实际堆石坝邓肯-张E-B模型参数概率分布的统计发现,对邓肯-张E-B模型参数的概率分布模型选用对数正态分布是合适的.对于变异系数的选取,岩土体物理力学参数指标的变异系数取值通常为0.1~0.5,本文偏保守的选取材料参数的变异系数.材料随机参数的选取及分布方式如表1所示.表中:K_b为体积模量;m为体积模量指数;c为黏聚力;Δφ为内摩擦角变化量.

表1 覆盖层邓肯-张E-B模型随机参数的选取及分布方式

Tab.1 Selection and distribution of stochastic parameters in Duncan-Chang E-B model of overburden

材料参数	混合卵石层	均值	标准差	变异系数	分布方式
R_f	0.74	—	—	—	—
K	1062.8	1062.8	158.53	0.15	对数正态
n	0.33	—	—	—	—
K_b	681.3	—	—	—	—
m	0.08	—	—	—	—
φ₀/(°)	51.04	51.04	7.84	0.15	对数正态
Δφ/(°)	9.1	—	—	—	—
c/kPa	0	—	—	—	—
ρ/(g·cm^-3)	1.30	1.30	0.13	0.10	对数正态

新窗口打开| 下载CSV

目前,土石坝的地震反应可分为等效线性和非线性两大类.等效线性方法通过迭代的手段来近似地反映土体的非线性,这种方法计算效率较高,且对大多数情况,尤其是应变较小(小于1%~2%),地震加速度小于0.3g~0.4g的情况下都能提供较合理的结果^[24].非线性方法从理论上讲更有优势,但在大型土石坝有限元模型动力计算时,往往耗时较长.由于本文研究进行的是大型土石坝模型的随机有限元分析,需要进行很多次的重复计算来统计分析,所以对于动力计算,本文采用等效线性模型,模型参数为确定性模型参数.

2.3 地震作用下防渗墙随机有限元结果分析

由动力有限元计算结果可知,防渗墙位于覆盖层内部,受周围约束较大,地震作用下,防渗墙整体的动位移反应较小.防渗墙大主应力分布呈现出从底部到顶部逐渐增加趋势,最大值在防渗墙顶部,且为压应力,说明防渗墙整体受力状态良好.在初始分析时,选取防渗墙边界,中间和随机3种方式进行分析发现,边界位置和随机位置的均值及变异系数均小于防渗墙中部的均值和变异系数,所以中间位置更具有代表性,故本文在防渗墙中间部位选择A、B、C这3个特征点,对其应力变形规律进行分析.特征点位置的选取如图6所示.

图6

新窗口打开| 下载原图ZIP| 生成PPT

图6 防渗墙特征点位置选取

Fig.6 Selection of characteristic points of cut-off wall

根据式(1)绘制出防渗墙各特征点沉降平均值随蒙特卡罗模拟次数的变化曲线发现,当蒙特卡罗模拟次数达到500次后最终结果都已平稳,故本文选出特征点A绘制沉降平均值d_A随蒙特卡罗模拟次数N'的变化曲线,如图7所示.由图7可知,当蒙特卡罗模拟次数达到300次左右时计算结果已趋于平稳,因此下文对500次蒙特卡罗的模拟结果进行统计分析.

图7

新窗口打开| 下载原图ZIP| 生成PPT

图7 防渗墙沉降平均值随蒙特卡罗模拟次数的变化曲线

Fig.7 Variation curve of average settlement of cut-off wall with Monte Carlo simulation times

将随机参数计算结果与确定性有限元的计算结果进行比较,若随机参数计算结果大于确定性有限元的计算结果,将其定义为超标,由此确定在材料参数不确定性下的超标概率为

(9)P=

\frac{n_{A}}{N}

式中:n_A表示超过确定性有限元计算结果的次数;N表示总的随机有限元的模拟次数.

图8为特征点A的防渗墙横河向、顺河向和竖直向的动位移概率分布情况.图中:f为频数;U₁、U₂、U₃分别为横河向动位移、顺河向动位移、竖直向动位移;C为累计比例.实线为随机有限元计算结果的均值;虚线为确定有限元计算结果.表2给出了防渗墙3个特征点随机动位移统计结果.表中:U_ave为动位移均值.从图8中可以看出,该特征点的横河向、顺河向和竖直向动位移均值分别为2.464、0.875和1.777 cm;其中横河向和竖直向动位移均值均大于确定参数下计算结果,且超标概率基本都超过95%(见表2),可见参数随机情况下的计算结果与确定参数的计算结果有明显差异,且随机参数下的横河向和竖直向动位移的最大值超过确定参数下计算结果的3.2%和5.8%.顺河向动位移均值虽然小于确定参数下的计算结果,但仍然有超过43.2%的概率大于确定性有限元计算结果(见表3),并且顺河向动位移最大值超过确定参数下计算结果的87.6%.对其结果进行游程检验,如表3中所示.由表3可见,其显著性均大于0.1,可认为随机有限元计算结果呈现出一定的随机性,因此在实际工程中应关注这种随机情况的发生.

图8

新窗口打开| 下载原图ZIP| 生成PPT

图8 特征点A防渗墙变形概率分布

Fig.8 Deformation probability distribution of cut-off wall at characteristic Point A

表2 防渗墙特征点随机动位移计算结果

Tab.2 Calculation results of characteristic point displacement of diaphragm wall

方向	特征点A			特征点B			特征点C
方向	U_ave/cm	C_v/%	P/%	U_ave/cm	C_v/%	P/%	U_ave/cm	C_v/%	P/%
横河向	2.464	0.48	98.0	0.875	11.39	43.2	1.777	1.26	95.8
顺河向	2.491	0.32	98.0	1.815	3.66	4.2	1.774	0.79	97.0
竖直向	2.543	0.09	97.6	2.522	0.26	4.4	1.710	0.26	96.2

新窗口打开| 下载CSV

表3 防渗墙特征点随机动位移游程检验

Tab.3 Run test of characteristic points of cut-off wall with dynamic displacement

检验		特征点A			特征点B			特征点C
检验		U₁	U₂	U₃	U₁	U₂	U₃	U₁	U₂	U₃
渐近显著性(双尾)		0.179	0.128	0.371	0.325	0.325	0.152	0.474	0.929	0.592
MC显著性(双尾)		0.212	0.142	0.418	0.372	0.372	0.178	0.486	0.976	0.652
95%置信区间	上限	0.176	0.111	0.375	0.330	0.330	0.144	0.442	0.963	0.610
95%置信区间	下限	0.248	0.173	0.461	0.414	0.414	0.212	0.530	0.989	0.694

新窗口打开| 下载CSV

图9为特征点A的防渗墙大主应力S_max和小主应力S_min概率分布情况.图中:拉应力为负,压应力为正.表4给出了防渗墙3个特征点随机应力统计结果,表中S_ave为应力均值.从图9中可以看出,该特征点的大主应力和小主应力均值分别为1.228 MPa和 -1.053 MPa;其中大主应力的均值小于确定参数下的计算结果1.251 MPa,超标概率达到20%;且大主应力的最大值为1.722 MPa,超出确定参数下计算结果的37.6%.小主应力的均值大于确定有限元计算结果 -0.915 MPa,各特征点超标概率均超过95%;且小主应力的最大值为 -1.36 MPa,该结果已经超出C25混凝土的设计抗拉强度值1.27 MPa,超出确定参数下计算结果的48.6%.表4还显示防渗墙动应力的最大变异系数为14.9%,大于表2中动位移的最大变异系数11.39%,说明相比动位移,覆盖层参数的随机性对防渗墙的应力影响更明显.对3个特征点的计算结果进行游程检验,结果如表5所示.表5显示,其显著性均大于0.1,也呈现出一定的随机性.

图9

新窗口打开| 下载原图ZIP| 生成PPT

图9 特征点A防渗墙应力概率分布

Fig.9 Stress probability distribution of feature Point A cut-off wall

表4 防渗墙特征点随机动应力计算结果

Tab.4 Calculation results of dynamic stress at characteristic points of cut-off wall

应力	特征点A			特征点B			特征点C
应力	S_ave/MPa	C_v/%	P/%	S_ave/MPa	C_v/%	P/%	S_ave/MPa	C_v/%	P/%
大主应力	1.228	4.9	20	1.003	4.4	38	1.219	14.9	29
小主应力	-1.053	7.7	96	-0.994	10.1	95	-0.575	6.9	97

新窗口打开| 下载CSV

表5 防渗墙特征点随机动应力游程检验

Tab.5 Dynamic stress run test of characteristic points of cut-off wall

检验		特征点A		特征点B		特征点C
检验		S_min	S_max	S_min	S_max	S_min	S_max
渐近显著性(双尾)		0.858	0.210	0.371	0.858	0.244	0.858
MC显著性(双尾)		0.892	0.202	0.364	0.892	0.236	0.892
95%置信区间	上限	0.865	0.167	0.322	0.865	0.199	0.865
95%置信区间	下限	0.919	0.237	0.406	0.919	0.273	0.919

新窗口打开| 下载CSV

3 结论

本文将随机场理论与有限元方法结合起来,在蒙特卡罗框架下,基于Cholesky分解法的随机场模拟方法,讨论了考虑覆盖层材料参数空间变异性下的坝基防渗墙在地震作用下的应力变形情况,并结合工程实例对计算结果进行了统计分析.主要结论如下:

(1) 基于Cholesky分解法的随机有限元并考虑材料参数的不确定性,将有限元分析软件Abaqus与Python语言结合起来,实现了随机场模拟的自动化以及随机参数的自动赋值,提高了计算效率.

(2) 地震作用下随机有限元计算结果的离散程度显示,覆盖层参数的随机性对防渗墙的应力影响更明显.对累计500次随机有限元结果进行游程检验中发现,防渗墙动位移的最大变异系数为11.39%,小于动应力的最大变异系数14.9%.

(3) 考虑覆盖层参数的变异性,坝基防渗墙的应力变形均不同程度的超过确定参数下的计算结果.累计500次的随机有限元计算结果显示,防渗墙特征点A和特征点C的动位移超标概率均超过了95%;动应力小主应力的超标概率超过90%,且小主应力拉应力的最大值为1.36 MPa,超出确定参数下计算结果的48.6%.由此可见,忽略材料参数的不确定性将会低估坝基防渗墙的地震反应,可能会对坝基防渗墙的运行产生一定的威胁,在设计阶段应考虑到这种情况的发生.

参考文献

原文顺序

文献年度倒序

文中引用次数倒序

被引期刊影响因子

[1]

任翔, 高大水, 高全, 等.

土石坝超深混凝土防渗墙变形与受力分析

[J]. 长江科学院院报, 2018, 35(5): 120-124.