馮龍龍 蘇曉麗 李 星

（1.河海大學(xué) 水文水資源與水利工程科學(xué)國(guó)家重點(diǎn)實(shí)驗(yàn)室，南京 210098；2.河海大學(xué) 水資源高效利用與工程安全國(guó)家工程研究中心，南京 210098；3.河海大學(xué) 水利水電學(xué)院，南京 210098）

混凝土面板堆石壩是由防滲面板，防滲接地結(jié)構(gòu)以及壩體堆石體組成的混合結(jié)構(gòu)，由于混凝土面板堆石壩具有工程量小、安全、經(jīng)濟(jì)、施工方便和適應(yīng)性好的特點(diǎn)［1］，所以它在實(shí)際工程中得到廣泛的應(yīng)用.隨著技術(shù)的發(fā)展和設(shè)計(jì)經(jīng)驗(yàn)的積累，面板堆石壩的設(shè)計(jì)高度已經(jīng)突破了200m，如水布埡水電站的壩高超過(guò)了230m.由于壩高的增加，使得壩體的應(yīng)力和變形變得比較復(fù)雜.因此，在設(shè)計(jì)混凝土面板堆石壩時(shí)，對(duì)其進(jìn)行應(yīng)力和變形分析是十分重要的.本文計(jì)算時(shí)以某高面板堆石壩為例，利用Marc軟件的非線性求解功能進(jìn)行應(yīng)力和變形分析，給出其相應(yīng)的變化規(guī)律.

1 有限元計(jì)算原理

非線性有限元按位移模式求解的平衡方程為

式中，［K（u）］為勁度矩陣，｛u｝為節(jié)點(diǎn)位移列陣，｛R｝為節(jié)點(diǎn)荷載列陣.

本次計(jì)算采用中點(diǎn)增量法［2］求解非線性方程組.壩的填筑順序把自重荷載分為若干級(jí)，配合荷載分級(jí)進(jìn)行有限元網(wǎng)格劃分.中點(diǎn)增量法把每級(jí)自重荷載稱為荷載增量，記為｛ΔR｝，如果第i級(jí)為｛ΔR｝i，先根據(jù)上一級(jí)荷載計(jì)算末的｛σ｝i－1確定彈性常數(shù) Ei－1和νi－1，并組成勁度矩陣［K］i－1，在結(jié)構(gòu)中施加本級(jí)荷載增量的一半｛ΔR｝i／2，用式（2）求得位移增量：

再計(jì)算應(yīng)力與應(yīng)變?cè)隽?，疊加到上次計(jì)算的應(yīng)變與應(yīng)力上，由此得到本級(jí)計(jì)算的中點(diǎn)應(yīng)變｛ε｝i－1／2，然后可以得到彈性常數(shù)Ei－1／2和νi－1／2.再用式（3）求本級(jí)全荷載發(fā)生的位移增量：

最后計(jì)算本級(jí)的應(yīng)變與應(yīng)力增量，累加到上一級(jí)終了的應(yīng)變與應(yīng)力上，即為本級(jí)的應(yīng)變與應(yīng)力，壩體上游面的水壓力按照蓄水順序分級(jí)加荷.

2 堆石體的本構(gòu)關(guān)系以及在Marc中的實(shí)現(xiàn)

2.1 堆石體的本構(gòu)模型

堆石體的本構(gòu)關(guān)系表達(dá)了堆石體應(yīng)力應(yīng)變之間的關(guān)系，它無(wú)疑對(duì)應(yīng)力應(yīng)變的計(jì)算結(jié)果起決定性的作用［3］，因而采用恰當(dāng)?shù)谋緲?gòu)模型對(duì)面板堆石壩應(yīng)力變形分析是至關(guān)重要的.堆石料是非線性材料，根據(jù)三軸實(shí)驗(yàn)，鄧肯－張（E－B）模型能夠較好地反映土體材料應(yīng)力應(yīng)變的非線性，因此選擇鄧肯－張（E－B）模型作為堆石料、墊層和過(guò)渡料等的本構(gòu)模型.該本構(gòu)模型的基本原理［4］為

式中，Rf為破壞比；Pa為大氣壓力；K、Kb、m、n為試驗(yàn)常數(shù)；S為剪應(yīng)力水平，反映材料強(qiáng)度的發(fā)揮程度，S＝，（σ1－σ3）f為破壞時(shí)的偏應(yīng)力，由摩爾－庫(kù)倫破壞準(zhǔn)則得：

2.2 在Marc軟件中的實(shí)現(xiàn)

Marc軟件作為大型有限元商業(yè)計(jì)算軟件之一，與其它軟件相比，具有能夠快速求解高度非線性問(wèn)題的優(yōu)點(diǎn).但是，它不包含工程中常用的鄧肯－張（E－B）非線性彈性模型，這嚴(yán)重影響了它在壩工領(lǐng)域中的使用.本文利用 Marc軟件提供二次開(kāi)發(fā)功能［5］，把通過(guò)實(shí)驗(yàn)得到的參數(shù)利用Fortran語(yǔ)言編寫子程序，利用Marc子程序來(lái)實(shí)現(xiàn)鄧肯－張非線性彈性模型的應(yīng)用，然后利用生死單元技術(shù)進(jìn)行分級(jí)加載，以進(jìn)行壩體施工和蓄水的模擬.

3 面板堆石壩實(shí)例分析

3.1 工程概況

某水庫(kù)工程位于河南省境內(nèi)，它的開(kāi)發(fā)任務(wù)主要是以防洪、供水為主，同時(shí)考慮發(fā)電和改善生態(tài)環(huán)境，并進(jìn)一步完善黃河下游調(diào)水調(diào)沙運(yùn)行條件.該工程為混凝土面板堆石壩，設(shè)計(jì)壩高124.0m，壩頂高程288.5m，頂部設(shè)防浪墻，其高度為1.2m，壩頂總長(zhǎng)為481.0m，寬為10.0m，設(shè)計(jì)時(shí)上游采用1∶1.5的坡率，下游的為1∶1.6.壩體主要由主堆石區(qū)、次堆石區(qū)、混凝土面板和防滲墻等幾部分組成，主堆石區(qū)與次堆石區(qū)分界線為壩軸線向下游1∶0.6.根據(jù)面板壩規(guī)范要求，100m以上的壩面板上游面下部設(shè)有上游鋪蓋區(qū)及蓋重區(qū)，上游鋪蓋的頂高程，取為壩體高度的35%，即頂部高程在200.00m，最大高度為35 m；在周邊縫下設(shè)小料區(qū)（特殊墊層）；下游坡面為大塊石護(hù)坡.

3.2 有限元模型及其計(jì)算參數(shù)

考慮壩體分區(qū)、施工程序及加載過(guò)程，并考慮到防滲墻的連接型式，對(duì)壩體及壩基進(jìn)行剖分，建立三維有限元模型如圖1所示，總共剖分9 862個(gè)單元，11 489個(gè)結(jié)點(diǎn).

圖1 有限元計(jì)算模型

混凝土面板、趾板、連接板等在達(dá)到破壞強(qiáng)度之前線性關(guān)系一般較好，故按線彈性材料處理，參數(shù)見(jiàn)表1.

表1 面板堆石壩的彈性材料參數(shù)

面板堆石壩中主堆石料、次堆石料、粘土夾層和過(guò)渡石料等按照非線性彈性材料進(jìn)行模擬，計(jì)算使用的參數(shù)見(jiàn)表2.

荷載加載時(shí)，利用Marc軟件提供的生死單元功能實(shí)現(xiàn)分級(jí)加載來(lái)模擬面板堆石壩的澆筑過(guò)程以及蓄水期的情況.按照面板堆石壩施工進(jìn)度和蓄水計(jì)劃，先后將荷載分為24級(jí)模擬.

表2 面板堆石壩的鄧肯－張模型（E－B）參數(shù)

3.3 計(jì)算結(jié)果分析

分析計(jì)算結(jié)果時(shí)，取3個(gè)剖面1、2、3（如圖2所示）進(jìn)行分析，由于篇幅有限，文中只給出2號(hào)剖面的等值線圖，其它剖面的極值列入相應(yīng)的表中.計(jì)算時(shí)選取的坐標(biāo)方向：x向?yàn)轫樅酉?，y方向?yàn)閴胃叻较?，向上為正，z軸為壩軸線方向，從左向右；應(yīng)力按照土力學(xué)的習(xí)慣，以壓應(yīng)力為正，拉應(yīng)力為負(fù).

圖2 剖面示意圖

3.3.1 壩體應(yīng)力和位移分析

圖3為竣工期堆石體2號(hào)斷面的應(yīng)力分布，圖4為滿蓄期堆石體2號(hào)剖面的應(yīng)力等值線圖，表3分別列出了3個(gè)剖面的應(yīng)力極值.

表3 剖面的應(yīng)力極值統(tǒng)計(jì)表 （單位：MPa）

圖3 竣工期2號(hào)剖面應(yīng)力等值線圖（MPa）

圖4 滿蓄期2號(hào)剖面應(yīng)力等值線（MPa）

由圖3～4和表3可以得到壩體應(yīng)力分布規(guī)律：竣工期壩體主應(yīng)力等值線與壩坡基本平行，從壩頂向壩基呈逐漸加大的趨勢(shì)，堆石體應(yīng)力極大值出現(xiàn)在2號(hào)剖面，第一主應(yīng)力為2MPa，第三主應(yīng)力為0.80 MPa；蓄水后，受水荷載作用，堆石應(yīng)力極值增大，所處的位置進(jìn)一步向上游主堆石區(qū)靠近，主堆石區(qū)應(yīng)力較次堆石區(qū)應(yīng)力有所增加，第一主應(yīng)力最大值為2 MPa，第三主應(yīng)力最大值為0.85MPa，極值仍然出現(xiàn)在2號(hào)剖面.

圖5～6分別為2號(hào)斷面竣工期和滿蓄期堆石體的位移等值線圖，3個(gè)斷面的位移極大值列于表4.

圖5 竣工期2號(hào)剖面位移等值線圖（cm）

圖6 滿蓄期2號(hào)剖面位移等值線圖（cm）

表4 堆石體的位移統(tǒng)計(jì)表 （單位：cm）

由圖5～6以及表4可以看出：由于堆石體的泊松效應(yīng)，使得橫向剖面上水平位移分布規(guī)律基本上是上游堆石區(qū)位移指向上游，下游堆石區(qū)位移指向下游，這符合竣工期面板堆石壩上下游方向位移分布的一般規(guī)律.

1號(hào)剖面向上游最大位移為3cm，位于上游側(cè)剖面1／3壩高處；向下游最大位移為7cm，位于下游側(cè)剖面1／3壩高處；最大豎向位移為40cm，位于剖面的1／2壩高處.2號(hào)剖面向上游最大位移為9cm，位于上游側(cè)基礎(chǔ)覆蓋層以上壩體1／4壩高位置，向下游最大位移為23cm，位于基礎(chǔ)覆蓋層以上壩體1／3壩高位置；最大豎向位移為96cm，位于剖面的1／3壩高處.3號(hào)剖面向上游最大位移為3cm，位于上游側(cè)剖面1／3壩高處；向下游最大位移為9cm，位于下游側(cè)剖面1／3壩高處；豎直向位移極值為42cm，位于斷面的1／2壩高處.

水庫(kù)蓄水后，在水荷載的作用下，壩體3個(gè)剖面的順河向位移均有較大變化，但仍然有一定的規(guī)律.1號(hào)剖面向上游最大位移為1.4cm，位于壩底附近；向下游最大位移為9.3cm，位于下游側(cè)剖面1／3壩高處；豎直向位移極值為40cm，位于剖面的1／2壩高處.2號(hào)剖面向上游最大位移為2cm，位于上游堆石區(qū)壩底部位置，向下游最大位移為27cm，位于下游次堆石區(qū)壩底部附近；最大豎向位移為98cm，位于剖面的1／3壩高處.3號(hào)剖面向上游最大位移為1.4cm，位于壩底部附近；向下游最大位移為11cm，位于下游側(cè)剖面1／3壩高處，豎直向位移為43cm，位于剖面的1／2壩高處.

3.3.2 面板應(yīng)力與變形分析

圖7～8分別是竣工期和滿蓄期面板的應(yīng)力圖.

圖7 竣工期面板應(yīng)力分布（MPa）

圖8 滿蓄期面板應(yīng)力分布（MPa）

由圖可知：竣工期，面板順坡向主要表現(xiàn)為壓縮變形，最大壓應(yīng)力為3.8MPa，出現(xiàn)在2號(hào)剖面182m高程位置.面板軸向應(yīng)力與順坡向應(yīng)力相比數(shù)值較小，基本呈兩岸受拉、中間受壓分布，最大壓應(yīng)力為0.8MPa，出現(xiàn)在2號(hào)剖面180m高程位置，最大拉應(yīng)力值為1MPa，出現(xiàn)在3號(hào)剖面190m高程位置.

蓄水后，面板順坡向應(yīng)力仍然以壓為主，最大壓應(yīng)力為4.1MPa，出現(xiàn)在2號(hào)剖面182m高程位置.但靠近河床兩端出現(xiàn)較大的拉應(yīng)力區(qū)，極值為1.8 MPa，出現(xiàn)在面板的最底部.相對(duì)于竣工期，軸向應(yīng)力有所增大，但總體呈兩岸面板受拉，中間受壓分布，符合軸向應(yīng)力的一般分布規(guī)律：2號(hào)剖面高程217m位置軸向壓應(yīng)力最大值為1.5MPa；相對(duì)于竣工時(shí)面板的拉應(yīng)力明顯增大，而且在靠近河床的兩端及兩岸周邊處會(huì)出現(xiàn)較大的拉應(yīng)力，最大值為2.5MPa.

圖9為竣工期面板位移分布等值線圖，圖10為滿蓄期面板位移分布等值線圖.

圖9 竣工期面板位移等值線圖（cm）

圖10 滿蓄期面板位移等值線圖（cm）

由圖9～10可以看出：竣工期，由于已經(jīng)歷了一期面板擋水，因此面板撓度基本指向壩內(nèi)，位于河床中央的面板底部撓度較大，撓度最大值為7cm；面板在壩軸向上由河谷向中央擠壓［6］，左岸軸向位移最大值為1.6cm（指向右岸）；軸向位移在右岸的最大值為1.8cm，方向指向河床中央.蓄水后，面板變形分布規(guī)律較好，面板撓度指向壩內(nèi)，面板中間區(qū)域數(shù)值較大，最大值為24cm；面板軸向有進(jìn)一步由河谷向中央擠壓的趨勢(shì)，軸向位移在左岸最大為3.4cm，右岸為3.7 cm，且分布具有一定的對(duì)稱性.

4 結(jié) 語(yǔ)

本文通過(guò)Marc軟件提供的子程序功能，利用Fortran語(yǔ)言編寫子程序以實(shí)現(xiàn)非線性彈性E－B模型在Marc中的應(yīng)用.根據(jù)三維非線性有限元分析可以得出：

1）三維情況下，壩體具有較強(qiáng)的空間效應(yīng)，河床中央剖面的最大主應(yīng)力值不在壩基面底部，而位于距離底部有一定距離的位置.

2）由于基礎(chǔ)有較深的覆蓋層，導(dǎo)致堆石體的最大沉降發(fā)生的位置以及面板最大撓度出現(xiàn)的位置相對(duì)于修建在基巖上的面板壩位置偏低一些，而且該工程由于次堆石下部參數(shù)較低、變形較大，故造成最大值位置偏向下游次堆石下部附近.

3）面板變形分布規(guī)律較好，面板撓度指向壩內(nèi)，壩軸向上由河谷向中央擠壓；蓄水后面板應(yīng)力仍然以壓為主，但是在靠近河床兩端出現(xiàn)較大的拉應(yīng)力，在實(shí)際工程中應(yīng)予以重視.

［1］ 孔憲京，張 宇，鄒德高.高面板堆石壩面板應(yīng)力分布特性及其規(guī)律［J］.水利學(xué)報(bào)，2013，44（6）：631－639.

［2］ 田 艷，王瑞駿，呂海東.積石峽面板堆石壩應(yīng)力變形非線性有限元分析［J］.水資源與水工程學(xué)報(bào)，2007，18（5）：71－74.

［3］ 高 輝，李守義.面板堆石壩非線性有限元應(yīng)力變形分析［J］.西北水資源與水工程，2000，11（4）：33－37.

［4］ 付國(guó)棟，王瑞駿，趙一新.施工期壩面過(guò)水面板堆石壩非線性有限元分析［J］.水資源與水工程學(xué)報(bào)，2009，20（6）：116－119.

［5］ 陳火紅，尹偉奇，薛小香.MSC.Marc二次開(kāi)發(fā)指南［M］.北京：科學(xué)出版社，2004.

［6］ 張 峰.深覆蓋層面板堆石壩應(yīng)力變形分析［J］.工程地質(zhì)計(jì)算機(jī)應(yīng)用，2012，68（4）：30－35.