吴之宝

(甘肃省武威市水利综合事务中心，甘肃 武威 733000)

虽然我国面板坝筑坝技术引入较晚，但随着西部大开发一大批混凝土面板堆石坝的兴建，现如今我国混凝土面板堆石坝筑坝技术已达到世界先进水平行列，积淀了较为丰富的设计、施工和运行管理经验。随着计算机科学技术的发展，我国在混凝土堆石坝数值模拟方面，取得了一大批科研成果，计算方法的研究也得到了迅速发展[1-3]，对于混凝土堆石坝的筑坝技术和数值仿真研究起到了积极的推进作用。堆石坝蓄水过程会引起大坝一系列应力和变形响应，如何合理准确模拟堆石坝的蓄水过程以及其对大坝应力变形规律的影响仍是大坝仿真分析的重点课题[4-5]。能否实现使用有限元方法对堆石坝应力场和位移场进行精准的计算，选择一种合乎实际情况地堆石体本构模型至关重要。本文采用邓肯张E-B模型[6-7]来模拟堆石体材料非线性特征，基于有限元软件ABAQUS平台，通过对三维有限元模型施加不同高度的静水压力来模拟堆石坝蓄水过程，对完建期大坝蓄水过程中引起的坝体应力和位移变化规律进行详细的分析研究。

1 计算理论与方法

土体的非线性应力应变关系可采用双曲线方程按式(1)作近似描述[8]。

(1)

式中：σ1为第一主应力；σ3为第三主应力；εa为轴向应变；Ei为初始切线弹模；(σ1-σ3)f为破坏时的第一、第三主应力之差；Rf为破坏比，见式(2)：

(2)

式中：(σ1-σ3)ult表示双曲线渐近线相映射的第一和第三主应力差；b为双曲线函数参数。其中，黏土的Rf=0.70～0.90，砂的Rf=0.60～0.85，砂卵石的Rf=0.65～0.85。

切线体积模量Bt与围压pa的关系如式(3)。

(3)

式中：Kb为体积模量数；m为无量纲系数。由于粗粒料的莫尔包线通常情况下并非为直线，故而要使用下列内摩擦角计算如式(4)：

(4)

式中：φ1为σ3等于大气压pa时的φ值；Δφ为σ3增大10倍时φ的减小值。

该模型用卸荷-再加荷模量Eur来代替Et，Eu r定义如下式(5)：

(5)

式中：Ku r为试验常数。一般情况下Ku r>K。

2 面板堆石坝分析模型的建立

由工程资料可知，大坝建基面高程2210.00 m，最大坝高133.00 m，某面板堆石坝正常蓄水位高程2271 m，上游回填区顶面高程2222 m，计算得上游水深49 m。坝体材料分区设计，共分为6个区，分别为垫层区、过渡料区、主堆石料区、下游堆石料区和下游回填区。大坝基岩主要材料组成为不透水的砂岩、砾岩。坝体-地基几何模型及材料分区如图1所示。图2显示出了堆石坝有限元模型。有限元网格采用六面体8节点solid单元划分，模型共划分314 646个单元和331 456个节点。

图1 堆石坝几何模型及材料分区

图2 堆石坝有限元模型

在模拟完建期后大坝蓄水过程时，堆石体采用邓肯E-B模型来模拟，坝体及上下游回填区力学参数如表1所示。

表1 坝体及上下游回填区邓肯张模型力学参数

3 堆石坝施工全过程应力及变形规律研究

3.1 堆石坝蓄水过程模拟

基于有限元软件ABAQUS平台，通过对三维有限元模型施加不同高度的静水压力来模拟堆石坝蓄水过程。为详细了解堆石坝蓄水过程中对坝体应力位移产生的影响，本文取坝前水深为15 m、25 m、35 m、49 m,分别对应工况一、工况二、工况三和工况四进行分析比较，其中工况一和工况四对应的坝前水位分别为死水位和正常蓄水位。

3.2 蓄水过程坝体应力变化

混凝土面板堆石坝坝前水深为15 m(工况一)时的计算结果如图3所示，为便于分析，文中显示的云图为坝体及上下游回填区1/2模型(下同)。从图中可以看出，当坝前水深为15 m时(此时为死水位时坝前水深)，坝体及上下游回填区最大拉应力为0.27 MPa，出现在下游回填区与坝肩接触部位中部。压应力峰值为2.00 MPa，位于坝体底部。与完建期相比，最大拉应力和主压应力均没有发生明显变化，说明在该堆石坝蓄水至死水位工况时，水荷载作用对坝体整体应力场改变较小。

图3 坝体应力云图(坝前水位15 m)

混凝土面板堆石坝坝前水深为25 m(工况二)时计算结果如图4所示。从图中可以看出，当坝前水深为25 m时，坝体及上下游回填区最大拉应力为0.27 MPa，同样出现在下游回填区与坝肩接触部位中部。值得注意的是，此时上游面板与坝肩接触部位也有较大拉应力出现，说明此时上游水荷载作用能够引起较大的局部应力场变化。与工况一相比，坝体及上下游回填区最大压应力略微增大，其值为2.01 MPa，同样出现在坝体底部位置。

图4 坝体应力云图(坝前水位25 m)

混凝土面板堆石坝坝前水深为35 m(工况三)时计算结果如图5所示。从图中可以看出，当坝前水深为35 m时，坝体及上下游回填区最大拉应力出现在上游面板与坝肩接触部位，与工况一和工况二相比拉应力有较大增长，其值为0.35 MPa。坝体及上下游回填区最大压应力则变化不大，其值为2.02 MPa，同样出现在坝体底部位置。

图5 坝体应力云图(坝前水位35 m)

混凝土面板堆石坝坝前水深为49 m(工况四)时计算结果如图6所示。从图中可以看出，当坝前水深为49 m时(正常蓄水位)，坝体及上下游回填区最大拉应力出现在上游面板与坝肩接触部位，相比工况三有较大增长，其值为0.51 MPa。坝体及上下游回填区最大压应力变化不大，其值为2.06 MPa，出现在坝体底部位置。

图6 坝体应力云图(坝前水位49 m)

3.3 蓄水过程坝体位移变化

图7显示给出了四种工况下坝体及上下游回填区的顺河向位移云图。从图7(a)可以看出，堆石坝坝前水深为15 m时，模型最大正向顺河向位移出现在坝体下游堆石区中部位置，其值为10.70 cm，坝体上游面靠近坝顶位置也出现较大正向顺河向位移。模型最大负向顺河向位移出现在坝体主堆石区中部位置，其值为4.96 cm。当坝前水深为25 m时(图7(b))，与工况一类似，模型最大正向顺河向位移同样出现在坝体下游堆石区中部位置，其值为11.00 cm。与前两种工况不同，当坝前水深为35 m时(图7(c))，大坝顺河向位移最大值出现在坝体上游坝面中心位置，其值为11.50 cm，模型最大负向顺河向位移出现在坝体主堆石区下部位置，其值为4.14 cm。随着坝前水位增加，大坝顺河向位移继续增大，在正常蓄水位工况下(图7(d))最终达到17.60 cm。

图7 蓄水期坝体位移云图

3.4 蓄水过程大坝应力及变形规律分析

混凝土面板堆石坝应力位移特征值随坝前水深变化情况如表2所示，从表中可以看出，在堆石坝蓄水过程中，随着上游水深增加，坝体最大主拉应力发生位置逐渐由下游回填区与坝肩接触部位转移至上游面板底部，正常蓄水位工况(工况四)下最大主拉应力值达到0.51 MPa。与主拉应力相比，大坝主压应力变化不大，正常蓄水位工况(工况四)下最大主压应力值为2.06 MPa。

表2 混凝土面板堆石坝应力位移随坝前水深变化情况

顺河向位移变化情况：与完建期相比，工况一引起的顺河向位移场没有发生明显变化，随着坝前水深的增加，水荷载引起的顺河向位移最大值逐渐增大，且顺河向位移最大值位置发生改变，出现在上游面板中部位置，在堆石坝蓄水过程中，坝体顺河向位移由完建期的10.60 cm增大至17.60 cm，增幅为66.0%；横河向位移变化情况：随着坝前水深的增加，水荷载引起的横河向位移最大值呈现出先减小后增大的趋势，在堆石坝蓄水过程中，坝体横河向位移由完建期的7.82 cm增大至8.64 cm，增幅为10.5%。从总体来看，堆石坝蓄水过程引起的坝体横河向位移变化量较小。竖向位移变化情况：与完建期相比，工况一引起的竖向位移场没有发生明显变化，随着坝前水深的增加，水荷载引起的坝体竖向位移最大值逐渐增大。在堆石坝蓄水过程中，坝体竖向位移由完建期的54.80 cm增大至57.70 cm，增幅为5.3%。

4 结 论

(1)在堆石坝蓄水过程中，坝前水深较小时，坝体最大主拉应力发生在下游回填区与坝肩接触部位中部位置。当坝前水深大于35 m时，最大主拉应力值明显增大，发生在上游面板底部位置，正常蓄水位工况下最大主拉应力值达到0.51 MPa，最大主压应力值为2.06 MPa。

(2)随着坝前水深的增加，与横河向位移和竖向位移相比，堆石坝蓄水过程能够引起较大的顺河向位移，正常蓄水位工况下顺河向位移最大值出现在上游面板中部位置。在堆石坝蓄水过程中，坝体顺河向位移由完建期的10.60 cm增大至17.60 cm，增幅为66.0%。由此可见，堆石坝蓄水过程会引起上游面板中部出现较大的顺河向位移，这在堆石坝设计时应引起注意。