【www.zhangdahai.com--其他范文】
邹德高，姜秋婷，刘京茂，金 伟，朱先文
(1.大连理工大学 海岸与近海工程国家重点实验室，辽宁 大连 116024；
2.大连理工大学 水利工程学院，辽宁 大连 116024；
3.中国电建集团成都勘测设计研究院有限公司，四川 成都 610072)
实测资料显示，高土石坝施工过程中心墙内会产生较高的超静孔隙水压力(简称孔压)。例如，小浪底斜心墙坝(坝高154 m)竣工时最大孔压为1.37 MPa[1]；
糯扎渡心墙坝(坝高261.5 m)满蓄时最大孔压为2.30 MPa[2]。这些监测成果为发展和验证高心墙坝安全评价方法、推动科学认识大坝工作性态提供了宝贵支撑，但已有工程普遍存在有效渗压计数目相对较少的问题。在建的高295 m两河口心墙坝建立了完备的监测系统，获得了丰富的孔压监测资料，目前最大监测孔压已达3.74 MPa[3]，同时发现孔压呈显著不均匀分布的特征，超出了以往的认识，引起了工程界对两河口大坝心墙安全的担忧。
为了分析孔压对工程安全的影响，众多学者提出并发展了基于有效应力原理的流固耦合分析方法，主要包括Biot流固耦合法[4-7]和不排水流固耦合法[8-9]两种，两者各有优缺点。Biot流固耦合法在理论上更严密，但面临着方程组系数矩阵病态问题，病态程度与时间步长、单元尺寸、渗透性等多种因素有关，并且求解效率低、易失败[10]，相关研究成果主要集中在二维分析。近年来，许多学者[4-5]开始对土石坝进行大型三维Biot流固耦合分析，但相关研究仍较少。不排水流固耦合法虽仅适用于土体饱和度较高、孔隙流体渗透性极低的情况，但该方法可在传统有限元计算软件中植入，实现方便，并且求解效率高、稳定性好。同时，大量研究[11-12]表明孔隙流体的压缩性与饱和度密切相关。为反映饱和度对孔压的影响，曹雪山等[13]在Biot流固耦合法中考虑了孔隙流体的压缩性，但其使用时仍存在前述问题。此外，心墙孔压模拟的准确性还取决于采用的本构模型。目前心墙力学特性模拟大都采用邓肯E-B本构模型，该模型计算的水平位移较实测偏大[14-15]，这会导致心墙孔压模拟失真。相比之下，弹塑性本构模型在模拟土体力学特性时结果更合理[16]。
综上，为了真实反映土石坝心墙孔压与饱和度的耦合效应及其三维分布特征，本文首先将饱和度作为初始材料参数引入到广义塑性本构模型中，建立了不排水条件下简化、高效的流固耦合分析方法，并验证了该方法在模拟心墙坝填筑孔压和变形时的适用性。基于心墙实测饱和度随机分布规律和三维流固耦合数值模拟，系统研究了两河口大坝心墙孔压分布特征及其形成机理。本文研究可为分析超高心墙坝填筑期孔压发展规律及其影响提供可靠的支撑手段。
两河口心墙坝位于四川省甘孜州雅江县，最大坝高295 m，坝顶宽16 m，河床部位心墙底高程2580 m，坝顶高程2875 m，正常蓄水位2865 m，上下游坝坡坡度分别为1∶2.0和1∶1.9(考虑马道的综合坝坡)。大坝施工过程中以左岸坝顶为原点，沿坝轴线共布设5个监测横断面。根据已有现场监测资料[3]，选择0+340 m横断面为典型断面，心墙中渗压计(PDB)和电磁式沉降环(DC4)埋设位置如图1所示。截至2021年8月，大坝填筑至高程2860 m，此时心墙最大沉降为2985 mm，位于0.5H(H为坝高)，如图2所示。
图1 典型断面心墙渗压计和电磁式沉降环分布
图2 电磁式沉降环测值沿高程分布
图3为心墙底部渗压计实测孔压时程以及填筑和水位的变化过程。可以看到：(1)不同位置的渗压计孔压发展是同步的，并且与填筑高程变化呈正相关关系；
(2)2021年5月至8月水位增长108 m，实测孔压仅增长0.1 MPa左右；
(3)实测最大孔压为3.74 MPa，明显高于同类工程(见研究背景部分)。说明填筑期上游库水还未来得及渗入心墙，高孔压主要由填筑重力引起的，因此本文分析中并未考虑上游水头入渗的影响。
图3 心墙渗压计实测孔压时程
图4为心墙渗压计测值沿高程的分布。如图4(a)所示，2019年8月不同横断面同一高程孔压的相对差值达52%；
同一横断面上下游对称位置的渗压计测值相对差值达91%。如图4(b)(c)，孔压随深度呈非线性递增规律；
不同高程孔压最大值所在位置沿高程是往复变化的，部分高程位于心墙边缘处。说明两河口大坝心墙孔压呈显著不均匀分布特征。
图4 心墙渗压计孔压测值
3.1 方法介绍为了反映两河口大坝心墙孔压与饱和度的耦合效应及其三维分布特征，首先将饱和度作为材料参数引入到广义塑性本构模型中，建立了不排水条件下简化、高效的流固耦合分析方法。具体实现过程如下。
(1)心墙有效应力-应变关系模拟。土石坝中筑坝料表现出明显的弹塑性变形特性，为了更合理地描述筑坝料的力学特性，邹德高等[17]考虑筑坝料的压力相关性和循环滞回特性，对原始广义塑性本构模型进行了改进，并成功应用于多座土石坝的静、动力反应分析[18-19]。因此，本文采用改进后广义塑性本构模型模拟心墙有效应力与应变的关系，两者关系表达为：
dσ′=Dep：dε
(1)
dε=dεe+dεp
(2)
(3)
式中：dεe为弹性应变增量；
dεp为塑性应变增量；
Dep为弹塑性矩阵；
De为弹性矩阵；
ng为塑性流动方向；
n为加载方向矢量；
H为塑性模量；
下标“L”和“U”分别表示加载和卸载。该模型的详细介绍本文不再赘述，模型参数及其物理意义详细简介见文献[20]。
(2)心墙孔压-应变关系模拟。土石坝心墙土填筑饱和度通常在80%～93%[13]，气体以封闭气泡的形式存在于水中，孔隙流体的压缩性由气体和水共同决定。因此，将孔隙水气作为一种可压缩混合流体进行模拟时，孔压和应变关系可描述为[8]：
dpw=Dw：dε
(4)
(5)
(6)
式中：Dw为孔隙流体的刚度矩阵；
dpw为孔压增量；
dε为应变增量；
Kaw为孔隙流体的体积压缩模量；
I为元素都为1的3×3矩阵；
O为3×3的零矩阵。
已有研究[21-22]表明，对于土体中液相相连通、气相不连通的气封闭状态(饱和度大于80%左右)，可忽略基质吸力的影响，则pw=pair，pair为孔隙气体压力。在该状态下，土体中孔隙水气的体积压缩模量可以表示为[23]：
(7)
式中：n为孔隙率；
Kf为纯水的压缩模量，Pa；
pw为孔压，Pa；
pa为1个大气压，约等于1×105Pa；
S为当前饱和度；
h为空气在水中的体积溶解系数，一般取值0.02。
随着pair(等于pw)的增加，水中的气泡会进一步溶解到水中，土体饱和度会逐级增大。依据Boyle定律和Henry定律，孔隙气体压力变化条件下当前饱和度S与初始饱和度S0关系可以表示为[24]：
(8)
(3)流固耦合方法的实现。土石坝心墙的渗透系数极小[13，25]，在填筑施工过程中流体渗流速度十分缓慢，此时土与孔隙流体呈同步变形的状态[8]。根据有效应力原理，将式(1)和(4)相加，可得心墙总应力增量与应变增量关系：
dσ=(Dep+Dw)：dε
(9)
从上式可以看到，只需在已有的本构模型中增加Dw并记录高斯点孔压即可实现心墙孔压的分析，并且当孔隙流体模量取0时，可自动退化为式(1)。由于有限元分析时外力是与总应力平衡的，因此上式很容易在单相介质中实现，并且无需修改应力积分和迭代求解流程。更为重要的是，因本文流固耦合方法本质是单相问题求解，因此避免了Biot流固耦合分析矩阵病态、求解速度慢的问题。
3.2 方法验证基于两河口心墙坝二维有限元模型，采用Biot流固耦合法对上述不排水流固耦合法的适用性进行了验证。两种方法中材料参数均相同，坝体填筑过程与实际工程保持一致(共历时6 a)，材料分区、填筑模拟等详细资料可见第4部分三维分析介绍。不排水流固耦合法中仅心墙采用不排水方法模拟，Biot流固耦合法中心墙两侧设定为排水边界。已有研究[1，26-28]表明，试验得到渗透系数会明显高估现场孔压的消散。为了与实际孔压累积和消散趋势保持一致，开展了不同心墙渗透系数的计算分析，结果表明，当渗透系数取1×10-10m/s时，计算的孔压消散发生的时间与实测时间是基本一致的(见图5渗透系数1×10-9m/s和1×10-10m/s的Biot流固耦合计算结果对比)，因此可认为心墙的渗透系数是1×10-10m/s。其值明显小于试验值，这与现场心墙非饱和、封闭气泡对孔隙水以及反滤料对细颗粒的阻滞作用有助于提高心墙的抗渗能力等因素有关。
图5也给出了不排水流固耦合法与Biot流固耦合法(心墙渗透系数为1×10-10m/s)计算结果的对比。结果表明，Biot流固耦合法心墙计算沉降略大于不排水流固耦合法计算沉降，最大差值为9%；
孔压则略小于不排水流固耦合法计算结果，最大差值为7.5%左右。综上，两种分析方法的计算结果差别较小，不排水流固耦合法可用于两河口心墙坝填筑期孔压和变形模拟。
图5 心墙沉降和孔压时程
4.1 有限元模型与材料参数
4.1.1 有限元模型 本次研究采用如图6所示的有限元模型。为模拟两河口心墙坝施工过程及水位变化，大坝分为49层填筑，其中围堰9层，坝体40层。采用实际建设过程的蓄水方案[3]，边填筑边蓄水，库水位分47步蓄至坝顶以下10 m。有限元模型采用了六面体等参单元和比例边界多面体单元模拟，单元总数为57 754个，结点总数为61 230个，自由度总数为183 690个。计算采用大连理工大学抗震研究所开发的岩土工程非线性有限元分析程序GEODYNA[29]。
图6 大坝有限元网格
4.1.2 材料参数 大坝计算模型包括围堰、堆石料、过渡料和心墙料等分区(见图6)。其中，心墙料参数根据室内试验[30]确定，如表1所示。
表1 广义塑性本构模型参数
4.2 计算工况两河口心墙坝施工过程中对心墙料含水量进行了监测，本文依据现场检测提供的160个数据，对心墙初始饱和度进行了分析。如图7所示，心墙实测初始饱和度S0基本分布在80%～100%，满足均值88%，方差0.001的正态分布N(88%，0.001)。
图7 心墙实测初始饱和度分布
为全面分析心墙饱和度对孔压和变形的影响，开展了如表2所示工况的比较研究。其中工况1：S0为100%；
工况2：S0为88%；
工况3：S0按N(88%，0.001)随机取值。考虑到监测条件有限，无法监测每个位置的饱和度，为了充分考虑初始饱和度的随机性，工况3进行了7组数值模拟。
表2 计算工况
4.3 结果分析
4.3.1 工况1和工况2 图8对比了工况1和工况2中渗压计PDB-60位置计算与实测孔压时程，两种工况中计算孔压与实测值均差别较大，最大差值分别为43%和-46%。图9对比了两个工况下心墙计算和实测孔压分布，可以得到：工况1较实测偏大，工况2较实测偏小，并且两种工况计算孔压离散性都较小，而实测值离散性较大。说明假定初始饱和度为常数难以反映实际孔压量值大小和孔压分布的显著不均匀性。
图8 心墙孔压时程(渗压计PDB-60)
图9 心墙孔压分布
图10给出了心墙沉降沿高程的分布图。工况1心墙最大沉降为3235 mm，位于0.63H，与实测最大沉降量值(2985 mm)和位置(0.5H)均存在较大的差别；
并且0.4H以下计算沉降较实测值偏小，0.4H以上较实测值偏大。工况2心墙最大沉降为2729 mm，位于0.53H，与实测吻合较好；
并且计算沉降与实测结果沿坝高分布规律是基本一致的，说明反映现场平均饱和度影响的变形计算结果与实际更吻合。同时对比3.2节二维计算沉降可知，不考虑大坝的三维空间效应会显著高估两河口大坝的计算变形。
图10 坝轴线沉降沿高程的分布
4.3.2 工况3 (1)孔压时程和分布。图11给出了渗压计PDB-60位置实测和工况3中7组计算孔压时程的对比。可以得到：计算孔压时程与初始饱和度密切相关，总体上初始饱和度越大计算孔压越大，但两者关系呈非线性变化关系；
当S0在92%～93%时，计算孔压时程与实测基本吻合，因此PDB-60位置的初始饱和度约为92%～93%。
图11 工况3中7组数值模拟的孔压时程(渗压计PDB-60)
图12对比了心墙实测孔压与工况3计算结果沿坝高的分布，限于篇幅仅给出了2组计算结果。如图12所示，3和7号数值模拟计算孔压分布与监测值基本一致(图12(a)(b))，7组数值模拟孔压基本能够覆盖实测孔压分布区域(图12(c))，说明考虑心墙饱和度随机性后，可数值再现心墙孔压不均匀分布规律。此外，不同数值模拟中孔压最大值在3.5～4.7 MPa，均位于心墙底部，但计算最大值大于监测最大值(3.74 MPa)，说明增加渗压计数量可能会监测到更大的孔压值。但高孔压也间接说明了心墙防渗性能良好，如果心墙渗透系数很大(见3.2节)，孔压会迅速渗透转移，不会出现局部高孔压的现象。
图12 心墙孔压分布散点图
综上，结合工况1和工况2计算结果可知，两河口大坝的心墙孔压的不均匀分布特征是心墙饱和度随机性导致的。
(2)沉降变形分布。图13给出了工况3心墙沉降沿高程的分布。7组数值模拟沉降完全重合，并与工况2的结果基本一致，最大差值仅为3%左右。说明保持平均饱和度不变，考虑饱和度随机性分布对心墙沉降变形影响不大，这是因为大坝沉降是土体变形的累积效果，与土体平均饱和度相关，受材料的局部特性影响较小。
图13 工况3中坝轴线沉降沿高程分布(未竣工)
(3)主应力分布。图14为工况3心墙上游侧单元的有效小主应力沿高程的分布，为简洁，图中仅给出了1号、3号和7号数值模拟的结果。可以看到：不同数值模拟中同一高程心墙有效小主应力具有明显的差别，但有效小主应力均大于0，未出现拉应力。
图14 工况3心墙上游侧有效小主应力(未竣工)
(4)应力水平分布。土体的受力破坏程度主要与其应力水平相关(或与破坏面的距离有关)。图15为工况3心墙中部一列单元的应力水平(应力水平=当前应力比/破坏应力比，等于1时表示破坏，等于0时表示等压状态)。可以看到：考虑饱和度的随机性对心墙应力水平分布也具有较大影响，不同数值模拟中同一高程心墙应力水平具有明显的差别，但心墙应力水平均低于0.85；
此外，虽然心墙底部中间孔压最大，但该位置应力水平低于0.45。说明两河口大坝心墙高孔压状态下的力学性态是良好的。
图15 工况3心墙应力水平(未竣工)
(1)将饱和度作为初始材料参数引入广义塑性本构模型，并考虑了饱和度随着围压变化的过程，建立了不排水条件下简化、高效的三维流固耦合分析方法，为分析高土石坝施工期心墙孔压分布规律及其影响提供了有效的工具。
(2)心墙饱和度初始值显著影响其孔压，饱和度越大孔压越大。当心墙初始饱和度以实测正态分布随机取值时，数值分析再现了两河口大坝心墙实测的高孔压以及孔压明显不均匀分布的特征，阐明了孔压分布规律与饱和度随机性是直接相关的。
(3)饱和度随机性对心墙沉降变形的影响低于3%，对心墙的局部应力和应力水平影响较为明显，但心墙有效小主应力均大于0，应力水平均低于0.85，且底部中间高孔压处的应力水平低于0.45。因此，饱和度随机性导致的高孔压以及孔压不均匀分布不会影响两河口大坝心墙的安全。
猜你喜欢 孔压心墙本构 时间平方根法评价隔离墙t50及固结系数广东土木与建筑(2022年9期)2022-10-15花扬子江诗刊(2021年3期)2021-11-11300 m级超高直心墙和斜心墙土石坝应力变形分析丽水学院学报(2021年5期)2021-10-31地铁行车荷载作用下粉质黏土累积孔压特性研究水利与建筑工程学报(2019年2期)2019-05-13竹节桩复合地基沉桩施工超孔隙水压力研究水文地质工程地质(2019年1期)2019-02-18锯齿形结构面剪切流变及非线性本构模型分析西南交通大学学报(2018年5期)2018-11-08高密度聚乙烯单轴拉伸力学性能及本构关系研究中国塑料(2016年4期)2016-06-27Therapeutic effect of okra extract on gestational diabetes mellitus rats induced by streptozotocinAsian Pacific Journal of Tropical Medicine(2015年12期)2015-10-31ABH沥青混凝土心墙坝应力应变分析中国工程咨询(2015年6期)2015-02-16常见的钛合金热变形本构模型分析及未来发展钛工业进展(2014年5期)2014-02-28本文来源：http://www.zhangdahai.com/shiyongfanwen/qitafanwen/2023/0823/643665.html
