高土石圍堰復(fù)合土工膜防滲三維有限元分析

徐 晗,陳 云,饒錫保,潘家軍

(長江科學(xué)院 水利部巖土力學(xué)與工程重點(diǎn)實(shí)驗(yàn)室，武漢 430010)

高土石圍堰復(fù)合土工膜防滲三維有限元分析

徐 晗,陳 云,饒錫保,潘家軍

(長江科學(xué)院 水利部巖土力學(xué)與工程重點(diǎn)實(shí)驗(yàn)室，武漢 430010)

如何在三維有限元計(jì)算中反映復(fù)合土工膜的力學(xué)特性是亟待解決的技術(shù)難題。總結(jié)了現(xiàn)有的土工膜數(shù)值模擬方法，在三維有限元計(jì)算中采用薄膜單元來模擬土工膜，建立了某高65 m的高土石圍堰土工膜防滲三維有限元模型，土工膜與上、下游壩料接觸部位均設(shè)置接觸。為模擬施工填筑與蓄水全過程，分別計(jì)算了土工膜邊界柔性約束與剛性約束2種工況。結(jié)果表明:薄膜單元較好地反映了土工膜柔性抗拉材料的力學(xué)特性，能適應(yīng)三維有限元數(shù)值模擬中常用的分步填筑施工模擬法，且后處理中能直接展示土工膜的位移、應(yīng)力云圖，具有方便高效的特點(diǎn)；土工膜邊界柔性約束位移分布均與剛性約束相同，只是最大主應(yīng)力及拉應(yīng)力分布區(qū)域略小于剛性約束，計(jì)算結(jié)果符合土石壩工程的一般規(guī)律。

復(fù)合土工膜；薄膜單元；三維有限元；高土石圍堰；變形特性

1 研究背景

復(fù)合土工膜(簡稱土工膜)具有防滲效果好、鋪設(shè)方便等優(yōu)點(diǎn)，可鋪設(shè)在壩體上游作斜墻或在土石壩中間作心墻，已廣泛應(yīng)用于水利工程建設(shè)中[1-2]。防滲體系是土石壩工程的生命線，如果土工膜拉裂或拉斷，將給整個(gè)工程的防滲帶來安全隱患?，F(xiàn)有規(guī)范[3]規(guī)定，高水頭(>50 m)擋水建筑物采用土工膜防滲應(yīng)經(jīng)過計(jì)算論證，由于土工膜厚度薄，不具備抗壓和抗彎特性，與壩料的力學(xué)特性差異較大，需要在土工膜與壩體之間設(shè)置非連續(xù)接觸面，因而在計(jì)算中如何考慮土工膜的實(shí)際作用存在較大的技術(shù)難度。

雖然學(xué)者們發(fā)展了解析法[4]、有限差分法[5]等多種方法來計(jì)算土工膜的應(yīng)力變形，但更多的是理論探討，較少運(yùn)用于實(shí)際土石壩工程。土石壩應(yīng)力變形計(jì)算多采用有限元法，在二維計(jì)算中可以將土工膜簡化為柔索單元[6]或不設(shè)置土工膜單元[7]，但土工膜三維有限元數(shù)值模擬研究甚少，所采用的力學(xué)模型仍需要進(jìn)一步研究論證，因而在三維有限元計(jì)算中考慮土工膜的作用是亟待解決的技術(shù)難題。

2 三維有限元土工膜模擬方法

在三維有限元計(jì)算中，主要采用如下幾種方法來考慮土工膜的作用。

(1) 土工膜內(nèi)嵌入實(shí)體單元：該方法將土工膜內(nèi)嵌入實(shí)體單元，并采用切向彈簧模擬土工膜與壩體之間的切向摩擦作用;在法向方向假定土工膜隨所嵌入的壩體單元一起運(yùn)動(dòng)。由于實(shí)際工程中蓄水時(shí)水壓力作用下土工膜與上、下游壩體之間有可能脫開，而采用該方法無法模擬此脫開性狀，同時(shí)水壓力也只能施加在所嵌入的單元面上，因而無法反映土工膜的真實(shí)受力性狀。

(2) 采用實(shí)體單元建立土工膜：該方法需要對(duì)土工膜建立實(shí)體單元，并劃分網(wǎng)格，土工膜上、下游側(cè)均與上、下游壩體之間設(shè)置接觸面，來模擬壩體與土工膜之間的接觸。雖然該方法可以在土工膜上加水壓力，并模擬土工膜與壩體之間的張拉脫開、剪切滑移等接觸力學(xué)特性，但是實(shí)際土工膜厚度只有1～2 mm，在三維模型中建模難度極大，且實(shí)體單元無法反映土工膜單元不具備抗壓和抗彎特性的特點(diǎn)。

(3) 土工膜依附在單元表面：該方法假定土工膜依附在單元表面上，不單獨(dú)建立土工膜單元，而根據(jù)加載前后土工膜所在位置單元表面的節(jié)點(diǎn)位移變化，來換算土工膜應(yīng)變，此計(jì)算方法非常粗略，僅能反映土工膜大致的受力性狀。

上述的幾種三維有限元計(jì)算中模擬土工膜的方法均存在明顯的缺點(diǎn)，本文采用ABAQUS提供的薄膜單元來模擬土工膜[8]。薄膜單元是一種面單元，只能傳遞面內(nèi)的力，不承受彎矩，即它沒有彎曲剛度，在ABAQUS中常用的薄膜單元為M3D4R (三維4節(jié)點(diǎn)減縮積分薄膜單元)，單元的正、負(fù)法向面均可以設(shè)置非連續(xù)接觸面。雖然有學(xué)者對(duì)薄膜單元在土石壩中的應(yīng)用進(jìn)行了初步的探討[9]，但并沒有模擬土工膜隨堰體施工填筑而逐步增高的過程，因而需要對(duì)薄膜單元是否適用于土石壩數(shù)值模擬中常見的填筑過程進(jìn)行進(jìn)一步的論證。

上述薄膜單元屬于柔性抗拉材料，可以直接設(shè)置土工膜的厚度，但不具備抗壓和抗彎特性；同時(shí)在土工膜與上、下游壩料之間均可設(shè)置接觸，接觸面的力學(xué)模型又可分為法向模型和切向模型，可以模擬土工膜與壩料之間的法向張拉脫開、切向剪切滑移等接觸力學(xué)特性；在此方法中，水壓力也直接加在土工膜上，因而在三維有限元計(jì)算中較好地反映了土工膜自身的基本力學(xué)特性。

表1 堰體材料的鄧肯-張E-μ模型參數(shù)Table 1 Parameters of Duncan-Zhang E-μ model of cofferdam material

3 有限元模型及模擬的填筑過程

某圍堰堰頂寬15 m，最大堰高65 m，堰體由石渣和砂礫石組成，采用混凝土防滲墻上接土工膜防滲，建立的三維有限元模型如圖1所示，防滲墻+土工膜防滲縱斷面如圖2所示，其中土工膜最大高度為45 m。

圖2 防滲墻+土工膜縱斷面Fig.2 Longitudinal section of impervious wall and geotechnical membrane

圖3 土工膜接觸面 設(shè)置示意Fig.3 Schematic diagram of the contact interface between geotechnical membrane and dam

三維有限元計(jì)算中采用薄膜單元來模擬土工膜，對(duì)于土工膜與壩料接觸部位也設(shè)置接觸單元，反映土工膜與壩料之間的接觸特性。土工膜與壩體之間的相互接觸如圖3所示。根據(jù)室內(nèi)試驗(yàn)成果，土工膜與堰體填料間的界面摩擦系數(shù)為0.29，防滲墻與覆蓋層之間的摩擦系數(shù)參照三峽二期圍堰實(shí)測(cè)值取0.19。復(fù)合土工膜采用線彈性模型，其厚度取2.5 mm，彈性模量為100 MPa。

模擬施工的全過程，荷載分級(jí)按堰體填筑次序進(jìn)行，模擬的施工程序?yàn)椋簢叩鼗鯌?yīng)力分析→分級(jí)填筑截流堤→水下分級(jí)填筑上游堆石→水下分級(jí)填筑下游堆石→施工防滲墻→水上分級(jí)填筑碾壓堆石→鋪設(shè)過渡料、壩料及土工膜→圍堰分級(jí)蓄水，下游降水到圍堰底。

堰體材料采用的鄧肯-張E-μ模型參數(shù)見表1，其中：c為內(nèi)聚力；φ為摩擦角；Rf為破壞比；K為初始變形模量的基數(shù)；n反映變形模量隨圍壓增長的關(guān)系；G為初始泊松比的基數(shù)；F反映初始泊松比隨圍壓增長而降低的速率；D反映初始泊松比隨軸向應(yīng)變?cè)鲩L的關(guān)系。

4 土工膜應(yīng)力變形分析

4.1 土工膜變形分析

圖4為土工膜沉降等值線，可知土工膜變形在堰體最大橫斷面處位移值最大。表2為蓄水期堰體與土工膜變形極值，可知土工膜與壩體豎向位移量值基本一致，沉降最大值發(fā)生在堰體中部，壩軸向位移量值也一致，只是順河向位移極值有所差異，這是由于堰體水平位移極值發(fā)生在下游堰體內(nèi)，如圖5所示。通過分析圖5可知，在土工膜與上、下游堰體交接處，豎向位移與水平位移等值線均有明顯的突變，表明該模型較好地反映了土工膜與上、下游堰體之間的非連續(xù)接觸變形特性。

圖4 土工膜變形等值線Fig.4 Deformation contours of geomembrane表2 蓄水期堰體與土工膜變形極值Table 2 Extreme values of the deformation of cofferdam and geomembrane impoundment stage

cm

圖5 最大橫斷面位移等值線Fig.5 Displacement contours at the maximum cross section

4.2 土工膜應(yīng)力分析

土工膜為柔性構(gòu)件，只能承受拉應(yīng)力，不能承受壓應(yīng)力，因而計(jì)算結(jié)果顯示了土工膜所受到的拉應(yīng)力數(shù)值及區(qū)域。

表3為剛性約束條件下土工膜拉應(yīng)力、應(yīng)變極值，圖6為土工膜最大主應(yīng)力、最大主應(yīng)變等值線。由表3可知，土工膜的拉應(yīng)力可以直觀地顯現(xiàn)出來，剛性約束下的最大主應(yīng)力0.67 MPa，最大主應(yīng)變0.71%。土工膜剛性約束后的拉應(yīng)力主要集中在土工膜的上部及周邊約束地區(qū)，這是由于該地區(qū)變形較大或者受約束所致；而由于土工膜與防滲墻交接地帶處于壩體中下部，該部位埋深較大，下游堰體的抗力較強(qiáng)，土工膜變形較小，因而沒有拉應(yīng)力產(chǎn)生。

表3 土工膜應(yīng)力、應(yīng)變極值Table 3 Extreme values of the stress and strain of geomembrane

圖6 土工膜最大主應(yīng)力、最大主應(yīng)變等值線Fig.6 Contours of maximum principal stress and maximum principal strain of geomembrane

綜上所述，土工膜不具備抗壓和抗彎特性，因而主要呈現(xiàn)拉應(yīng)力，沒有壓應(yīng)力產(chǎn)生，計(jì)算結(jié)果較好地體現(xiàn)了土工膜柔性抗拉材料的力學(xué)特性。

5 邊界約束對(duì)土工膜應(yīng)力變形的影響

工程實(shí)踐表明，在土工膜與周邊結(jié)構(gòu)剛性約束時(shí)，在接頭處容易產(chǎn)生應(yīng)力集中，致使土工膜拉裂[10]。為了驗(yàn)證土工膜邊界約束施加方法的正確性，在薄膜單元與結(jié)構(gòu)連接處，分別施加了剛性約束與柔性約束2種邊界條件，剛性約束為土工膜邊界各位移分量全限制，柔性約束則釋放部分位移分量，土工膜的應(yīng)力變形結(jié)果見表4。

表4 邊界約束對(duì)土工膜應(yīng)力變形的影響結(jié)果Table 4 Influence of boundary constraint on the stress and deformation of geomembrane

由表4可知柔性約束條件下，蓄水期土工膜水平位移最大值為11.6 cm(向下游)，豎向位移最大值為43.2 cm，壩軸向位移最大值為5.9 cm，與剛性約束邊界位移值相比，蓄水期土工膜各向位移均與剛性約束下相同；只是最大主應(yīng)力為0.59 MPa，最大主應(yīng)變?yōu)?.57%，略小于剛性約束條件。

圖7為蓄水期不同邊界約束條件土工膜拉應(yīng)力區(qū)域云圖，可知在土工膜與兩側(cè)山體搭接部位容易產(chǎn)生應(yīng)力集中，柔性約束條件下土工膜的拉應(yīng)力區(qū)域及數(shù)值明顯小于剛性約束條件，因而在土石壩工程中，土工膜邊界應(yīng)盡量采用柔性約束。

6 結(jié) 論

(1) 采用ABAQUS提供的薄膜單元來模擬土工膜，可在土工膜與上、下游壩料接觸部位設(shè)置接觸，薄膜單元能較好地反映土工膜柔性抗拉材料的力學(xué)特性，在土石壩三維計(jì)算中能適應(yīng)常用的分步填筑施工模擬法，且后處理中能直接展示土工膜的位移、拉應(yīng)力與應(yīng)變等值線，具有方便高效的特點(diǎn)。

(2) 當(dāng)土工膜邊界為剛性約束時(shí)，在蓄水期水壓力的作用下，在土工膜與兩側(cè)山體搭接部位易產(chǎn)生應(yīng)力集中；土工膜邊界柔性約束下各向位移分布均與剛性約束相同，只是最大主應(yīng)力及拉應(yīng)力分布區(qū)域小于剛性約束，計(jì)算結(jié)果符合土石壩的一般規(guī)律。

[1] 顧淦臣.復(fù)合土工膜或土工膜堤壩的實(shí)例述評(píng)(續(xù))[J].水利水電技術(shù), 2003, 34(1): 56-59.

[2] 包承綱. 三峽二期圍堰創(chuàng)新及其對(duì)我國高圍堰建設(shè)的意義[J]. 長江科學(xué)院院報(bào), 2014, 31(9): 33-42.

[3] SL/T 225—98, 水利水電工程土工合成材料應(yīng)用技術(shù)規(guī)范[S]. 北京: 中國水利水電出版社, 1998.

[4] 高 登, 朱 斌, 陳云敏. 圓形構(gòu)筑物周邊土工膜沉降時(shí)的褶皺產(chǎn)生機(jī)理及分析模型[J]. 巖土工程學(xué)報(bào), 2008, 30(8): 1126-1132.

[5] 吳兆和, 沈長松, 瞿忠烈, 等. 復(fù)合土工膜堆石壩膜受力變形特性分析[J]. 水電能源科學(xué), 2010, 28(3): 91-94.

[6] 花加鳳, 束一鳴, 張貴科, 等. 土石壩壩面防滲膜中的夾具效應(yīng)[J]. 水利水電科技進(jìn)展，2007, 27(2): 66-68.

[7] 姜海波. 復(fù)合土工膜心墻與斜墻高土石壩應(yīng)力應(yīng)變研究[J]. 長江科學(xué)院院報(bào), 2014, 31(1): 53-57.

[8] 莊 茁. ABAQUS有限元軟件6.4版入門指南[M]. 北京: 清華大學(xué)出版社, 2004.

[9] 尚 層.土石壩復(fù)合土工膜防滲斜墻應(yīng)力變形分析[D].烏魯木齊: 新疆農(nóng)業(yè)大學(xué), 2012.

[10]李 波, 程永輝, 程展林. 圍堰防滲墻與復(fù)合土工膜聯(lián)接型式離心模型試驗(yàn)研究[J]. 巖土工程學(xué)報(bào), 2012, 34(11): 2081-2086.

(編輯：黃 玲)

Three-dimensional Finite Element Analysis on Composite Geomembraneof High Earth Rockfill Cofferdam

XU Han，CHEN Yun，RAO Xi-bao，PAN Jia-jun

(Key Laboratory of Geotechnical Mechanics and Engineering of Ministry of Water Resources,Yangtze River Scientific Research Institute, Wuhan 430010, China)

To reflect the mechanical properties of composite geomembrane in 3-D finite element calculation is a technical problem to be solved. On the basis of summarizing current numerical methods of geomembrane simulation, we propose to simulate the geomembrane using membrane element in 3-D finite element calculation. By building a finite element model of high earth-rock cofferdam of 65 m height with geomembranes contacted with dam in the upstream and downstream, we simulated the flexible and rigid constraint conditions of geomembrane boundary respectively in filling construction and impounding stages. Results suggest that the membrane element could well reflect the flexible tensile properties of geomembrane. The proposed method could also adapt to the stepped filling in 3-D finite element simulation. Moreover, it could conveniently and effectively display the displacements and stresses of geomembranes in post-processing. The displacement distribution of geomembrane under flexible boundary restraint was the same with that under rigid boundary constraint; while the principal tensile stress and its distributed region under flexible boundary constraint were slightly smaller than those under rigid constraint. The calculation results are consistent with the general rules of earth-rockfill dam.

composite geomembrane; membrane element; 3D finite element; high earth rockfill cofferdam；deformation characteristics

2016-06-23;

2016-08-23

國家自然科學(xué)基金項(xiàng)目(51309029)；長江科學(xué)院創(chuàng)新團(tuán)隊(duì)項(xiàng)目(CKSF2015051/YT)；長江科學(xué)院院所基金項(xiàng)目(CKSF2015036/YT)

徐 晗(1978-)，男，湖北公安人，教授級(jí)高級(jí)工程師，博士，主要從事巖土工程數(shù)值分析研究，(電話)027-82829743(電子信箱)mechanics007@aliyun.com。

10.11988/ckyyb.20161056

TV641.1

A

1001-5485(2017)02-0110-04

2017,34(2):110-113