薛志平 ，孫一清 ，張 震

（1.山西省柏葉口水庫建設(shè)管理局，山西 呂梁 033000；2.河海大學(xué)水利水電學(xué)院，江蘇 南京 210098）

0 引言

柏葉口水庫位于山西省呂梁市交城縣會立鄉(xiāng)文峪河干流上，控制流域面積 875 km2。壩址區(qū)河谷底寬約為 120 m，河谷形狀呈 U 形，主要建筑物為面板堆石壩，壩軸線走向?yàn)?N 63.23°E。壩頂高程為1 138.30 m，最大壩高為 88.30 m，是山西省最高的大壩，其安全對于呂梁市和汾河中下游的防洪安全具有十分重要的意義。

現(xiàn)場檢查發(fā)現(xiàn)大壩壩頂及上游面板存在多條裂縫，由于該壩施工速度比較快，前期監(jiān)測系統(tǒng)不夠完備，更新后的監(jiān)測系統(tǒng)運(yùn)行時間不長，監(jiān)測資料匱乏。為清楚地分析大壩變形規(guī)律，全面地分析快速施工對大壩安全的影響，有必要通過數(shù)值手段全面了解大壩各種不利工況下的變形情況[1-4]。為此結(jié)合柏葉口水庫大壩實(shí)際工程資料，仿真模擬大壩堆石料的填筑、混凝土面板澆筑及蓄水過程，建立考慮具體施工過程分級加載的大壩三維有限元分析模型，運(yùn)用 COMSOL Multiphysics 數(shù)值分析軟件獲取壩體在竣工期、正常蓄水位、設(shè)計及校核洪水位工況下的變形特性。

1 有限元分析模型

1.1 建模過程分析

要建立有限元分析模型，先要在 Auto CAD 中建立準(zhǔn)確的混凝土面板堆石壩幾何模型。運(yùn)用 Auto CAD 繪圖軟件繪制大壩三維幾何模型，坐標(biāo)系的選取如下：X 軸沿壩軸線指向右岸方向?yàn)檎?，零點(diǎn)取在壩軸線左岸處；Y 軸指向水流方向?yàn)檎泓c(diǎn)取在壩橫右 0 + 000.00 m 斷面處；Z 軸垂直指向上方為正，與高程一致。為提高計算結(jié)果的準(zhǔn)確性，計算模型在 X，Y，Z 各軸的邊界范圍如下：在豎直方向上，底部邊界為 Z = 913.5 m 的截面；在順河向上，上游的邊界為 Y = -306.718 m 的截面，下游邊界為Y = -298.293 m 的截面；在壩軸線方向上，左邊界為 X = 486.6 m 的截面，右邊界為 X = -176.6 m 的截面。為避免導(dǎo)入 COMSOL Multiphysics 仿真計算時畸變單元的產(chǎn)生，對大壩幾何模型進(jìn)行優(yōu)化，避免扭曲面的產(chǎn)生[5-7]。優(yōu)化過程中不斷使用 Live Link功能，連接 Auto CAD 與 COMSOL Multiphysics，實(shí)現(xiàn)邊修改邊導(dǎo)入的數(shù)據(jù)同步。柏葉口大壩壩體三維有限元網(wǎng)格模型如圖 1 所示，剖分后得到 336 632 個單元。

圖 1 柏葉口大壩三維有限元網(wǎng)格模型（單位：m）

1.2 荷載分級

真實(shí)模擬施工荷載是混凝土面板堆石壩有限元分析中至關(guān)重要的一步，可以體現(xiàn)結(jié)構(gòu)體本身隨施工過程的變化情況，從而更好地反映材料的非線性，使計算結(jié)果更符合實(shí)際。通過工程資料可知，壩頂高程為 1 138.30 m，正常蓄水位為 1 133.00 m。在混凝土面板堆石壩的施工中，壩體的填筑過程是分層碾壓，具體施工程序如下：1）混凝土趾板澆筑；2）壩體次堆石區(qū)高程 1 098.00 m 以下填筑；3）壩體主堆石區(qū)高程 1 098.00 m 以下及過渡層、墊層填筑；4）壩體堆石區(qū)高程 1 098.00 m 至壩頂及過渡層、墊層填筑；5）面板澆筑；6）鋪蓋、防浪墻及下游護(hù)坡填筑。

初次蓄水時，根據(jù)蓄水安全鑒定報告與壩體結(jié)構(gòu)設(shè)計和度汛要求，參考已建工程經(jīng)驗(yàn)，水荷載是分級施加的，結(jié)合壩體網(wǎng)格單元的剖分情況，具體加載方式如下：1）水位 1 061.9～1 080.0 m 蓄水；2）水位 1 080.0～1 100.0 m 蓄水；3）水位 1 100.0～1 128.0 m 蓄水；4）水位 1 128.0～1 133.0 m 蓄水到正常蓄水位。

計算工況中還應(yīng)考慮運(yùn)行期水位升至設(shè)計和校核洪水位的工況，選取竣工期未蓄水、正常蓄水位、設(shè)計和校核洪水位 4 種計算工況依次分析面板堆石壩變形分布。

1.3 參數(shù)選取及計算工況確定

1.3.1 材料參數(shù)的選取

壩體填筑料采用鄧肯-張（E-B）模型，相關(guān)參數(shù)根據(jù)現(xiàn)場試驗(yàn)參照類似工程確定，如表 1 所示。基巖地層和混凝土采用線彈性模型，計算參數(shù)依據(jù)表 2 所示資料。

1.3.2 沉降結(jié)果的分析

大壩外部沉降從 2011 年 6 月開始觀測（1 078，1 098，1 118 高程觀測點(diǎn)），2012 年 6 月壩頂高程開始觀測，測點(diǎn)分布如圖 2 所示。通過觀測資料分析，最大沉降發(fā)生在 W12觀測點(diǎn)，沉降量為 70.32 mm，壩頂最大沉降發(fā)生在 W9觀測點(diǎn)，沉降量為 45.63 mm，沉降遵循逐漸下沉規(guī)律，沒有明顯的突變點(diǎn)。觀測統(tǒng)計資料[8]如表 3 所示。

可見，通過 Auto CAD 三維模型導(dǎo)入 COMSOL Multiphysics 進(jìn)行網(wǎng)格剖分，選取合適的材料及材料參數(shù)、本構(gòu)模型進(jìn)行模擬，通過數(shù)值模擬出的沉降計算值與實(shí)際觀測的沉降值對比，發(fā)現(xiàn)計算值偏小。原因分析如下：1）未將溫度及時效對壩體沉降的影響考慮進(jìn)數(shù)值模擬中；2）面板縫與周邊縫建模考慮不夠全面；3）水滲入面板造成的沉降改變等諸多因素造成。

表 1 非線性材料本構(gòu)關(guān)系的計算參數(shù)表

表 2 線彈性材料應(yīng)力應(yīng)變關(guān)系計算參數(shù)表

圖 2 柏葉口大壩外部測點(diǎn)分布圖

表 3 大壩外部沉降觀測資料統(tǒng)計表 mm

2 壩體變形分析

分析過程中，最先模擬基巖的初始地應(yīng)力場。隨著自重和水壓力等的施加，堆石料和面板的位移及應(yīng)力不斷變化。到最后加載完成，會較穩(wěn)定地呈現(xiàn)一定的規(guī)律[9-11]。

計算結(jié)果中位移的正負(fù)與坐標(biāo)正負(fù)一致，在分析該混凝土面板堆石壩的變形時，選取大壩典型斷面河床處 1-1 斷面（壩 0 + 180.00 斷面），左岸 2-2斷面（壩 0 + 060.00 斷面），右岸 3-3 斷面（壩 0 +245.00 斷面），分析這 3 個典型斷面在竣工期、正常蓄水位、設(shè)計及校核洪水位 4 種工況下的變形。

2.1 竣工期變形分析

2.1.1 河床 1-1 斷面

壩體河床處 1-1 斷面竣工期位移分布如圖 3 所示。由于斷面上下游相對呈現(xiàn)大致對稱結(jié)構(gòu)，計算結(jié)果顯示，竣工期的位移也具有較好的對稱性。由圖 3 a 可得，竣工期未蓄水時，壩體上下游最大位移約為 4.50 mm；由圖 3 c 可得，壩體的最大沉降約為55.90 mm，由于存在鋪蓋層（1 A，1 B），沉降最大值出現(xiàn)在壩體的中部偏上 2/3 處。

圖 3 壩體河床 1-1 斷面在竣工期的位移（單位：mm）

2.1.2 左岸 2-2 斷面

竣工期壩體左岸 2-2 斷面位移分布如圖 4 所示。由圖 4 a 可得，竣工期未蓄水時，左岸壩體上下游最大位移約為 2.80 mm；由圖 4 c 可得，壩體的最大沉降約為 34.10 mm，由于存在 1 A，1 B 鋪蓋層，沉降最大值出現(xiàn)在壩體的中部偏上 2/3 處。

圖 4 壩體左岸 2-2 斷面在竣工期的位移（單位：mm）

2.1.3 右岸 3-3 斷面

竣工期壩體右岸 3-3 斷面位移分布如圖 5 所示。由圖 5 a 可得，竣工期未蓄水時，右岸壩體上下游最大位移約為 2.88 mm；由圖 5 c 可得，壩體的最大沉降約為 37.50 mm，由于存在 1 A，1 B 鋪蓋層，沉降最大值出現(xiàn)在壩體的中部偏上 2/3 處。

2.2 正常蓄水位變形分析

正常蓄水位工況下，壩體河床處 1-1、左岸 2-2及右岸 3-3 等斷面中，壩體上下游最大位移出現(xiàn)在河床 1-1 斷面，且大于左右岸斷面，約為 7.30 mm。由于存在 1 A 和 1 B 鋪蓋層，3 個斷面的沉降最大值均出現(xiàn)在壩體的中部偏上 2/3 處，且河床 1-1 斷面的壩體沉降最大，約為 60.30 mm。

圖 5 壩體右岸 3-3 斷面在竣工期的位移（單位：mm）

2.3 設(shè)計洪水位變形分析

設(shè)計洪水位工況下，壩體河床處 1-1、左岸 2-2及右岸 3-3 等斷面中，壩體上下游最大位移出現(xiàn)在河床處 1-1 斷面，且大于左右岸斷面，約為 7.63 mm。由于存在 1 A 和 1 B 鋪蓋層，3 個斷面的沉降最大值均出現(xiàn)在壩體的中部偏上 2/3 處，且河床 1-1 斷面的壩體沉降最大，約為 61.30 mm。

2.4 校核洪水位變形分析

校核洪水位工況下，壩體河床處 1-1、左岸 2-2及右岸 3-3 等斷面中，壩體上下游最大位移同樣在河床處 1-1 斷面，且大于左右岸斷面，約為 7.75 mm。由于存在 1 A 和 1 B 鋪蓋層，3 個斷面的沉降最大值均出現(xiàn)在壩體的中部偏上 2/3 處，且河床 1-1 斷面的壩體沉降最大，約為 61.60 mm。

3 結(jié)語

本研究用非線性有限元法的理論及計算方法、堆石料的本構(gòu)模型、施工逐級加載方式對柏葉口混凝土面板堆石壩變形進(jìn)行了研究，解決了以下問題：1）導(dǎo)入 COMSOL Multiphysics 商用軟件，使用自由四面體網(wǎng)格，進(jìn)行了模型的有限元網(wǎng)格剖分；2）竣工期垂直方向壩體的最大位移（沉降）為 55.60 mm，出現(xiàn)在壩體的中上部分，約占最大壩高的 0.33%，向下游的水平位移大于向上游的水平位移。正常蓄水位下，壩體在豎直方向上，壩體的最大沉降為 60.30 mm。隨著水位升高壩體沉降不斷增大，設(shè)計和校核洪水位相對應(yīng)的沉降量為 61.30和 61.60 mm。由此可知柏葉口混凝土面板堆石壩整體安全，本工程施工速度快對整體變形規(guī)律改變不大。下一步將引入考慮蠕變的本構(gòu)模型，進(jìn)行深入研究，從而為動態(tài)變形監(jiān)控指標(biāo)奠定基礎(chǔ)。