基于Wolffersdorff模型面板堆石壩應(yīng)力應(yīng)變數(shù)值計(jì)算

黃小華，吳樂海

(福建水利電力職業(yè)技術(shù)學(xué)院，福建 永安 366000)

在全球能源緊缺，碳中和政策的背景下，新能源的開發(fā)已變得越來越受關(guān)注。利用流量和落差發(fā)電的水電能源能實(shí)現(xiàn)零碳排放，符合當(dāng)下的能源要求。作為當(dāng)前壩工界主要壩型之一的面板堆石壩因其具有施工方便，投資相對(duì)較低，在高度提升情況下占地較小等優(yōu)點(diǎn)，在水利工程中得到廣泛應(yīng)用[1]。由于面板堆石壩防滲的主體是混凝土面板，受力主體為堆石料，若在設(shè)計(jì)階段對(duì)堆石料變形預(yù)測(cè)不當(dāng)，后期可能會(huì)導(dǎo)致面板脫空、拉裂等情況的發(fā)生，從而造成嚴(yán)重滲漏，還有可能造成潰壩的風(fēng)險(xiǎn)，因此對(duì)堆石料應(yīng)力應(yīng)變的研究顯得尤為重要[2- 3]。由于堆石料具有剪脹剪縮、非線性等性質(zhì)，因此在堆石料的應(yīng)力應(yīng)變分析中選擇一個(gè)合理的本構(gòu)模型至關(guān)重要。目前，彈塑性本構(gòu)模型和非線性彈性模型是國(guó)內(nèi)外壩工界對(duì)混凝土面板堆石壩力學(xué)分析的2大主要本構(gòu)模型。但是，這2類本構(gòu)模型都不能很好地反映堆石料剪脹剪縮、非線性等力學(xué)特征，在考慮顆粒破損方面也較為欠缺，且忽略了屈服面與應(yīng)力路徑的關(guān)系[4- 6]。近年來，國(guó)內(nèi)外有學(xué)者將亞塑性理論應(yīng)用到無黏性土的本構(gòu)模型研究中，提成各種亞塑性本構(gòu)模型，該類模型能彌補(bǔ)傳統(tǒng)本構(gòu)模型的不足，能更好地反映堆石料特殊的力學(xué)特征。

本文簡(jiǎn)單介紹Wolffersdorff亞塑性本構(gòu)模型及其有限元計(jì)算方程。并將該亞塑性本構(gòu)模型應(yīng)用于福建仙游抽水蓄能電站下庫(kù)混凝土面板堆石壩力學(xué)數(shù)值計(jì)算中，同時(shí)將該模型計(jì)算的數(shù)值成果與沈珠江雙屈服面彈塑性本構(gòu)模型計(jì)算的數(shù)值成果進(jìn)行對(duì)比。旨在驗(yàn)證該本構(gòu)模型應(yīng)用于面板堆石壩應(yīng)力應(yīng)變分析的可行性。

1 Wolffersdorff亞塑性本構(gòu)模型

亞塑性本構(gòu)模型理論是在尋求顆粒骨架之間的級(jí)配性質(zhì)的基礎(chǔ)上[7]，以連續(xù)介質(zhì)力學(xué)理論，通過建立應(yīng)力率與應(yīng)變率的張量函數(shù)，在函數(shù)中加入其他本構(gòu)理論，如非線性理論、考慮顆粒破損[8]、濕化效應(yīng)[9]、剪脹性[10]等情況，對(duì)亞塑性本構(gòu)模型進(jìn)行不斷完善，從而得到適用于堆石料的本構(gòu)模型。Wolffersdorff[11]對(duì)砂礫料進(jìn)行研究，在前人研究壓縮性模型的基礎(chǔ)上引入Druker-Prager臨界準(zhǔn)則，加入了臨界狀態(tài)函數(shù)“F”，使得模型在模擬顆粒破壞狀態(tài)下的應(yīng)力應(yīng)變關(guān)系更加真實(shí)。Wolffersdorff模型表達(dá)式見式(1)所示。

(1)

(2)

模型中參數(shù)α和β的求取，吳長(zhǎng)彬[14]在常規(guī)三軸壓縮試驗(yàn)下得到與平均壓力ps有關(guān)，為簡(jiǎn)化其關(guān)系，引入?yún)?shù)k1、b1、k2、b2，如式(3)—(4)所示：

α=k1ln(T33/pα)+b1

(3)

β=k2ln(T33/pα)+b2

(4)

簡(jiǎn)化后的參數(shù)共11個(gè)：λ、hs、n、φc、ei0、ed0、ec0、k1、b1、k2、b2，使得參數(shù)獲取更加容易，根據(jù)福建仙游抽水蓄能電站的筑壩材料試驗(yàn)報(bào)告數(shù)據(jù)[15]計(jì)算得到主堆石料和次堆石料的模型參數(shù)見表1：

2 亞塑性材料矩證模量[D]的確定

有限元的理論是將整體進(jìn)行分割成單元，將單元利用虛功原理外力功與應(yīng)力功平衡，列出平衡方程，再對(duì)單元進(jìn)行疊加求得整體平衡方程。單元的平衡方程如下：

(5)

式中，{δ*}e—節(jié)點(diǎn)虛位移；{ε*}—虛應(yīng)變列陣，可表示為{ε*}=[B]{δ*}e；{σ}—單元應(yīng)力列陣{σ}=[D][B]{δ}e。

[K(D)]{δ}={R}

(6)

式中，{δ}—整體節(jié)點(diǎn)位移列陣；{R}—整體結(jié)點(diǎn)荷載列陣；[K]—整體勁度矩陣；[D]—材料矩陣模量。因此上式有限元平衡方程的求解取決于[D]，由率型本構(gòu)方程[16]：

(7)

(8)

因此，可得到矩陣[D]的表達(dá)式為：

(9)

對(duì)式(8)進(jìn)行展開：

(10)

表1 模型參數(shù)表

提出式(2)中的Δε，結(jié)合式(10)，可得到矩陣模量[D]分量：

(11)

式中，下標(biāo)i=1，3;j=1，3，從而可得整體的勁度矩陣[K]求解出有限元平衡方程(6)。

3 亞塑性本構(gòu)模型在工程中的應(yīng)用

3.1 工程概況及計(jì)算方法

為了驗(yàn)證改進(jìn)的Wolffersdorff亞塑性本構(gòu)模型應(yīng)用于面板堆石壩應(yīng)力應(yīng)變分析的可行性，本次選用的工程案例為福建仙游抽水蓄能電站下庫(kù)混凝土面板堆石壩，該大壩的基本參數(shù)如下：壩高75.1m，壩頂寬度8.0m，壩長(zhǎng)263m，前坡坡比：1∶1.4，后坡坡比：1∶1.35，混凝土面板的材料特性參數(shù)取γ=24.0kN/m3、μ=0.167、E=2.8×107kPa。大壩的斷面圖如圖1所示，對(duì)大壩有限元網(wǎng)格剖分如圖2所示。

圖1 大壩斷面圖(單位：m)

圖2 大壩有限元網(wǎng)格圖

對(duì)式(6)式非線性有限元方程的求解，采用中點(diǎn)增量法，對(duì)仙游面板堆石壩的下庫(kù)大壩進(jìn)行分28級(jí)加載，施加本級(jí)荷載的一半力{ΔR}i于結(jié)構(gòu)，得到本級(jí)的平均值，使得結(jié)構(gòu)的精度得到提高，中點(diǎn)增量法示意圖如圖3所示。

圖3 中點(diǎn)增量法

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

3.2.1壩體應(yīng)力應(yīng)變結(jié)果分析

在對(duì)福建仙游抽水蓄能電站下庫(kù)面板堆石壩的三維有限元分析中，堆石料的本構(gòu)模型分別采用改進(jìn)的Wolffersdorff亞塑性本構(gòu)模型和沈珠江彈塑性模型，計(jì)算工況為竣工期和正常蓄水期。2種本構(gòu)模型計(jì)算的壩體應(yīng)力及變形結(jié)果見表2，各種工況下的水平位移、垂直位移及大小主應(yīng)力等值線圖如圖4—7所示。其中，水平方向的位移以向下游為正，垂直方向的位移以向下為正，應(yīng)力以壓為正。

表2 不同模型壩體應(yīng)力應(yīng)變結(jié)果對(duì)比

圖4 大壩河床斷面水平位移等值線分布圖

圖5 大壩河床斷面垂直位移等值線分布圖

圖6 大壩河床斷面第一主應(yīng)力等值線分布圖

圖7 大壩河床斷面第三主應(yīng)力等值線分布圖

由表2可以看出，竣工期Wolffersdorff模型計(jì)算出的最大水平位移為4.3cm，最大垂直位移為24.8cm，沈珠江模型最大水平位移6.9cm，最大垂直位移24.0cm?？⒐て赪olffersdorff模型計(jì)算出的最大水平位移比沈珠江模型更大，增幅達(dá)到60%，2種本構(gòu)模型計(jì)算出的最大垂直位移較為接近，約占?jí)胃叩?.33%；蓄水期Wolffersdorff模型最大水平位移5.1cm，最大垂直位移28.8cm，沈珠江模型最大水平位移7.9cm，最大垂直位移26.5cm。蓄水期Wolffersdorff模型計(jì)算出的最大水平位移比沈珠江模型增大55%，最大垂直位移的比較中，Wolffersdorff模型也比沈珠江模型更大2.3cm。從計(jì)算結(jié)果上看，2種模型位移的計(jì)算結(jié)果偏差不大，在水壓力作用下，下游向的水平位移比上游向的水平位移來的更大。從位移等值線分布云圖分析，兩者模型分布特征基本相似，在堆石的自重作用下，壩體垂直位移最大值出現(xiàn)在1/2壩高處，并向四周逐漸減小。水平位移大致成對(duì)稱性，呈下凸上凹形狀，略有不同在于蓄水期Wolffersdorff模型壩體上游側(cè)中上部位指向下游側(cè)，沈珠江模型指向上游側(cè)。

由表2可知，竣工期Wolffersdorff模型第一主應(yīng)力最大值為1.02MPa，第三主應(yīng)力最大值為0.37MPa，沈珠江模型相應(yīng)的應(yīng)力為1.04，0.48MPa，從計(jì)算結(jié)果可以看出，在竣工期，2種本構(gòu)模型第一主應(yīng)力計(jì)算數(shù)值非常接近，第三主應(yīng)力二者差值為0.11MPa，從應(yīng)力的等值線分布圖來看，兩者模型的應(yīng)力分布規(guī)律相似，大小主應(yīng)力均出現(xiàn)在大壩的底部。蓄水期Wolffersdorff模型第一主應(yīng)力最大值為1.11MPa，第三主應(yīng)力最大值為0.42MPa，對(duì)比沈珠江模型為1.12、0.52MPa，2種模型計(jì)算結(jié)果非常接近。從圖7還可以看出，2種模型的第三主應(yīng)力最大值較竣工期均有往上游移動(dòng)的趨勢(shì)。

3.2.2面板應(yīng)力應(yīng)變結(jié)果分析

蓄水期2種本構(gòu)模型計(jì)算的面板應(yīng)力變形計(jì)算成果見表3、如圖8—9所示，規(guī)定應(yīng)力以壓為正，變形以下游向?yàn)檎?/p>

表3 不同模型蓄水期面板應(yīng)力應(yīng)變結(jié)果對(duì)比

如圖8蓄水期面板位移云圖所示，Wolffersdorff模型和沈珠江模型撓度位移和軸線向位移分布規(guī)律相同，最大撓度部分均出現(xiàn)在壩體1/2～1/3區(qū)域，并從壩體兩側(cè)向中間逐漸變大。從云圖可以看出，混凝土面板最大擾度的位置與壩體堆石料最大垂直位置非常接近。從表3的變形計(jì)算結(jié)果可知，沈珠江模型最大撓度為7.8cm，軸線向位移最大為1.1cm，最小為-1.1cm。Wolffersdorff模型最大撓度為11.3cm，軸線向位移最大為1.7cm，最小為-1.6cm。Wolffersdorff模型計(jì)算的最大擾度比沈珠江模型增加45%，相比下Wolffersdorff模型的計(jì)算結(jié)果偏于安全。

圖8 蓄水期面板位移云圖

如圖9蓄水期面板應(yīng)力云圖所示，2種模型的順坡向和軸線向的應(yīng)力分布規(guī)律相似，都能夠體現(xiàn)順坡向應(yīng)力中部受壓，頂部和底部受拉的狀態(tài)，軸線向應(yīng)力中部受壓，兩邊受拉狀態(tài)，應(yīng)力狀態(tài)分布合理。沈珠江模型最大順坡向應(yīng)力為2.11MPa，最小為-0.22MPa，軸線向最大應(yīng)力為1.1MPa，最小為-0.4MPa。Wolffersdorff模型最大順坡向應(yīng)力為5.06MPa,最小為-0.1MPa，軸線向最大應(yīng)力為2.68MPa，最小為-0.4MPa。從計(jì)算結(jié)果可以看出，Wolffersdorff模型計(jì)算的面板應(yīng)力分布比沈珠江模型更為合理。沈珠江模型計(jì)算的混凝土面板底部順坡向拉應(yīng)力存在應(yīng)力過大的問題，與實(shí)際情況不符。

圖9 蓄水期面板應(yīng)力云圖

4 結(jié)論

(1)Wolffersdorff模型的建立是以連續(xù)介質(zhì)力學(xué)為基礎(chǔ)，沈珠江模型是以彈塑性模型為基礎(chǔ)，兩者理論不同，但兩模型通過對(duì)仙游面板堆石壩的壩體及面板應(yīng)力應(yīng)變計(jì)算對(duì)比分析，兩者結(jié)果相差不大，從實(shí)際驗(yàn)證了改進(jìn)后的Wolffersdorff模型可以應(yīng)用于實(shí)際工程的應(yīng)力應(yīng)變分析。

(2)在混凝土面板應(yīng)力應(yīng)變計(jì)算方面，相比于沈珠江彈塑性模型，Wolffersdorff模型計(jì)算參數(shù)簡(jiǎn)單，容易確定，使用簡(jiǎn)便，且能夠反映堆石料剪脹性，應(yīng)力應(yīng)變非線性等力學(xué)特征。但是該模型也存在一些不足，如不能考慮堆石料流變特征等，有待進(jìn)一步研究。