前坪水庫砂礫石黏土心墻壩三維流固耦合分析

歷從實，武穎利，皇甫澤華，張兆省，韓 艷

（1.河南省前坪水庫建設管理局，河南 鄭州450003；2.河南省水利第一工程局，河南 鄭州450000；3.南京水利科學研究院，江蘇南京210098）

流固耦合作為自然界普遍存在的一種現象,其蘊含著復雜的物理現象和深刻的數學原理［1］。近年來,國內外學者在很多領域對流固耦合數值模擬和算法進行了探索研究,并取得了一定成果［2-6］。在水利工程領域,王曉玲等［7］對考慮流固耦合作用的拱壩抗滑穩定可靠度進行了研究；李冰冰等［8-9］對考慮流固耦合作用的尾礦壩、土石壩滲流穩定進行了研究；劉占濤等［10-11］對基于流固耦合作用的深厚覆蓋層和復雜壩基土石壩穩定進行了研究；陳文元等［12］對考慮流固耦合作用的壩體地震動力特性進行了研究；唐克東等［13］對不同開度弧形閘門流固耦合進行了數值模擬,取得了一定的成果。砂礫石黏土心墻壩作為土石壩的一種,相關流固耦合模擬研究文獻記載不多。本研究針對前坪水庫工程砂礫石黏土心墻壩的特點,結合地質、設計和相關試驗資料,采用南水雙屈服面本構模型,結合三維流固耦合計算方法,模擬壩體的填筑過程、施工后至蓄水前、蓄水過程以及運行期的壩體變形規律、應力分布規律和工作性態。

1 前坪水庫工程壩體結構特點

前坪水庫工程是國務院確定的172項節水供水重大水利工程之一。壩型為砂礫石黏土心墻壩,設計壩高90.3 m,壩頂長810 m,壩頂寬10 m,上游坡比1∶2.0～1∶2.5,下游坡比1∶2.0。黏土心墻頂寬為4 m,上下游坡比均為1∶0.3,緊貼心墻外側為2 m寬粗砂反濾層,粗砂反濾層外側為2 m寬碎石反濾層,碎石反濾層外側為砂礫石壩殼料。右岸壩肩坡比為1∶0.75,左岸壩肩坡比為1∶1。

2 三維有限元模型及計算參數

采用基于ABAQUS有限元軟件二次開發工具UMAT數據接口編寫的南水雙屈服面本構模型進行計算。根據設計圖紙,建立三維有限元模型。模型深度自建基面向下取80 m,兩側自壩底邊線向外延伸350 m作為截斷邊界。模型底部施加三向約束,四周施加法向約束。壩體三維有限元模型如圖1所示。模型共326 648個單元、345 836個節點,單元類型為C3D8P,在反濾料-心墻、防滲墻-壩體間設立薄層單元替代接觸面單元。

計算模擬分79個計算步。第1步為地應力平衡步,建立初始地應力場,消除覆蓋層初始沉降位移影響；第2～72步為土石壩施工填筑過程；第73～78步為土石壩蓄水至校核水位的過程；第79步模擬校核水位下持續運行10 000 d,研究壩體的工后變形。采用流固耦合方法對填筑施工期、竣工60 d后、蓄水期及蓄水運行20 a四種工況土石壩的受力變形特性進行計算分析。

圖1 壩體三維有限元模型

南水雙屈服面模型兼顧了鄧肯模型和劍橋模型的優點,考慮了土體的剪脹特性,且對各類土體適應性強,砂礫石壩殼料、心墻料、高塑性黏土、砂礫石覆蓋層均采用南水雙屈服面彈塑性本構模型［14］進行計算,對于基巖、混凝土防滲墻、防滲帷幕采用線彈性模型進行計算。

通過試驗和地質勘查資料獲取計算參數,見表1。

表1 壩體各區域主要計算參數

3 壩體應力變形計算分析

3.1 施工期壩體應力變形特性

圖2為壩體0+280斷面施工期應力應變特性圖。

由圖2可知,施工期0+280斷面最大沉降量為118.8 cm,發生在約1/2壩高位置,最大水平位移為16.3 cm,發生在心墻與上游側反濾料交接面靠近壩頂附近,方向指向下游側,最小水平位移為-17.3 cm,發生在最大水平位移的對稱位置,方向指向上游側。從大小主應力分布規律來看,施工期壩體大小主應力均為壓應力；大小主應力在心墻與反濾層交界位置有明顯跌落,心墻內部應力拱比較明顯。施工期壩體整體應力水平不高,僅在混凝土防滲墻與高塑性黏土交界位置附近因變形不協調導致應力劣化而出現較高的應力水平,壩體絕大部分位置應力水平較低。

3.2 竣工后60 d壩體應力變形特性

壩體0+280斷面竣工結束60 d的應力變形特性為：壩體最大沉降量為123.1 cm,較施工期增大了4.3 cm,最大值發生的位置與施工期一致。最大水平位移為16.3 cm,最小水平位移為-17.4 cm,水平位移分布規律與施工期一致,最小水平位移略有減小?？⒐そY束60 d壩體大小主應力分布規律與施工期基本一致,心墻內部大小主應力和應力水平隨著孔壓的消散略有調整。

4.3 蓄水期壩體應力變形特性計算分析

圖3為壩體0+280斷面蓄水期的應力、變形特性圖。

由圖3可知,蓄水期0+280斷面最大沉降量為126.7 cm,較竣工期沉降量增大了3.6 cm,沉降分布規律與施工期基本一致。最大水平位移為34.6 cm,發生在上游壩體壩殼料上部區域,方向指向下游；最小水平位移發生在下游壩坡中上部位置,最小值為-6.0 cm,方向指向上游。蓄水期大小主應力分布有所調整,心墻內部大小主應力變化比較明顯。蓄水期壩體應力水平整體不高。蓄水期壩體滲流場分布合理,符合一般心墻壩的孔壓分布規律。

圖2 壩體0+280斷面施工期不考慮流固耦合作用下的應力變形特性

圖3 壩體0+280斷面蓄水期應力變形特性

3.4 滿蓄狀態運行20 a壩體應力變形特性

圖4 為壩體0+280斷面蓄水運行20 a后的應力應變特性圖。

由圖4可知,壩體滿蓄狀態運行20 a后0+280斷面壩體沉降量為128.9 cm,較蓄水期增大了2.2 cm,沉降分布規律基本不變。壩體應力分布和應力水平分布與蓄水期基本保持一致。

3.5 壩體耦合（蓄水）前后最大應力變形特性對比

選取壩體典型斷面0+280、0+550、0+650,對不同工況計算數據進行統計分析,見表2。

壩體大部分沉降在施工期完成,施工期黏土心墻殘余超孔壓未完全消散,在竣工60 d后,壩體各典型斷面的最大沉降量均略有增大；蓄水后,受水荷載影響以及心墻內部孔壓分布的調整,壩體最大沉降量繼續增大；在大壩完成蓄水,保持滿庫運行20 a后,壩體最大沉降量較蓄水期有所增大,并趨于穩定。竣工期壩體最大斷面位置水平位移基本保持上下游側對稱分布,其他斷面受地形和壩體不對稱的影響水平位移極值有一定差異；蓄水期受水荷載影響,蓄水對壩體水平變形比較敏感；滿庫運行20 a后,壩體水平位移極值受壩體內部孔壓調整影響有微小變化,極值為35.2 cm,整體規律與蓄水期保持一致。黏土心墻在竣工期和蓄水期的大小主應力最大值受孔壓調整的影響逐漸減小,黏土心墻應力水平最大值也表現為逐漸減小。

圖4 壩體0+280蓄水運行20 a后的應力變形特性

表2 壩體耦合（蓄水）前后最大應力、變形比較

3.6 壩體施工期實際監測最大沉降量與理論計算最大沉降量對比

壩體施工期埋設多種安全監測儀器,施工期統計至2019年2月的沉降量監測數據見表3。

由表3可知壩體施工期黏土心墻最大沉降率為1.35%,由表2可知理論計算的壩體施工期黏土心墻最大沉降率為1.25%,二者吻合較好。

4 結 論

通過對前坪水庫砂礫石黏土心墻壩采用三維流固耦合計算分析,得出如下結論：

（1)通過對比監測資料和計算結果可知,壩體施工期實際監測最大沉降率與理論計算沉降率比較吻合,施工期完成了大部分沉降。

（2)黏土心墻大小主應力分布較為合理,不同時期心墻的應力水平最大值均出現在高塑性黏土與混凝土防滲墻交界區域,心墻區域絕大部分應力水平在0.6以下,表明黏土心墻有較高的強度安全儲備,心墻在不同時期不會發生剪切破壞。

（3)從心墻水平位移分布規律來看,竣工期壩體最大最小水平位移分別為17.7 cm和-17.3 cm,水平位移基本沿壩軸線左右對稱；竣工60 d水平位移較竣工期略有變化,量值很?。恍钏谛膲ψ畲?、最小水平位移分別為34.6 cm和-7.4 cm,發生在心墻頂部,心墻表現為向下游側變形；滿蓄運行20 a后,心墻水平位移的分布規律基本與蓄水期一致。

（4)采用流固耦合分析計算方法,通過考慮施工時間效應、蓄水時間效應以及運行期長期變形的模擬,分析了壩體在各個時期的變形特征和應力分布,體現了心墻固結特性對壩體變形和應力分布的影響。根據計算結果,工程建設管理部門在對壩頂道路、防浪墻以及圍護欄施工時,應選擇合適的施工時期和分縫方式,以消除壩體工后變形對建筑物的影響。

表3 壩體施工期實際監測沉降量統計