土石壩心墻裂縫發(fā)展的擴(kuò)展有限元模擬

趙曉龍,卞漢兵,鄭 威,孫兆輝,章賽澤,邱秀梅*

1.河海大學(xué) 巖土力學(xué)與堤壩工程教育部重點實驗室,江蘇 南京 210098

2.河海大學(xué) 巖土工程科學(xué)研究所,江蘇 南京 210098

3.山東農(nóng)業(yè)大學(xué) 水利土木工程學(xué)院,山東 泰安 271018

4.LEM3,CNRS 7239，洛林大學(xué),梅茲 法國 57045

5.徐州市水利建筑設(shè)計研究院,江蘇 徐州 221112

土石壩水力劈裂是亟待解決的問題之一[1,2]。目前關(guān)于水力劈裂的數(shù)值模擬，學(xué)者們做了很多工作，如有限元、邊界元和無單元法等，其中有限元應(yīng)用最為廣泛和成熟。李全明等[3]引入裂縫彌散理論來描述裂縫發(fā)展過程，從而實現(xiàn)水力劈裂的有限元模擬。王俊杰[4]等引入四節(jié)點等參單元建立裂縫及影響區(qū)的有限元模型，實現(xiàn)對裂縫的模擬。上述研究成果對于水力劈裂數(shù)值研究有很好的啟發(fā)意義，但由于有限元法本身的缺陷，使其模擬裂縫開展遇到很多問題。

有限元在處理裂紋問題時，需要根據(jù)裂縫的位置分布網(wǎng)格，在裂尖區(qū)域需要分布密集的網(wǎng)格來反映應(yīng)力的變化，當(dāng)裂縫擴(kuò)展時網(wǎng)格需要重新劃分，因此計算效率很低，對于復(fù)雜的裂縫擴(kuò)展問題甚至無法解答[5]。

在此背景下，Belytschko教授于1999年提出了擴(kuò)展有限元法XFEM[6,7]，用以解決有限元難以解決的不連續(xù)問題，如裂紋擴(kuò)展等。與常規(guī)有限元不同，擴(kuò)展有限元分析裂縫時，網(wǎng)格無需重新劃分，而是在包含裂縫的單元內(nèi)增加自由度來描述不連續(xù)體的位移間斷，位移間斷通過富集函數(shù)來描述，可以顯示裂紋的擴(kuò)展路徑，反映土體在流固耦合作用下變形特性，因而具有優(yōu)越性[5,8]。

近年來，國內(nèi)一些學(xué)者[9-12]開始將擴(kuò)展有限元應(yīng)用到混凝土壩體材料的開裂分析中，但土石壩粘土心墻水力劈裂的擴(kuò)展有限元模擬鮮有報道。

本文基于擴(kuò)展有限元，以損傷力學(xué)為理論基礎(chǔ)，參照混凝土等材料的建模思路，考慮到基質(zhì)吸力對黏土力學(xué)性能和裂縫發(fā)育的影響，建立了壓實黏土的孔隙擴(kuò)展模型和裂縫張開模型，對裂縫的發(fā)展進(jìn)行了模擬。該思路對于土石壩的水力劈裂數(shù)值研究有一定的借鑒性。

1 心墻黏土水力劈裂模型

1.1 黏土孔隙擴(kuò)展模型

1.1.1 模型本構(gòu)關(guān)系 根據(jù)損傷力學(xué)理論，將土體內(nèi)部的孔隙看作微裂紋（細(xì)觀裂紋），土體由于微裂紋發(fā)展而產(chǎn)生類似損傷的演化，將土體模擬為連續(xù)介質(zhì)，引入損傷變量，以熱力學(xué)為基礎(chǔ)，將孔隙的擴(kuò)展視為能量轉(zhuǎn)化的過程，由自由能和耗散勢導(dǎo)出孔隙擴(kuò)展的本構(gòu)關(guān)系和演化方程，描述材料從孔隙擴(kuò)展到出現(xiàn)宏觀裂紋的過程。

壓實黏土微裂紋演化的動力學(xué)特征在張拉和壓縮應(yīng)力作用下是不同的，為說明其演化對加載路徑的依賴性，使用標(biāo)量ωc、ωt表示微裂紋在受壓和受拉情況下的演變。構(gòu)造總體裂紋影響系數(shù)ω來表征兩種裂紋對黏土力學(xué)性能的影響：

式中αt和(1-αt)用來描述張拉和壓縮產(chǎn)生的微裂紋對材料力學(xué)性能總體影響的貢獻(xiàn)率。αt的定義見下式：

于是，壓實黏土的本構(gòu)關(guān)系可以直接通過自由能函數(shù)，即方程導(dǎo)出：

通過自由能函數(shù)，可推出兩個損傷變量，即在受壓和受拉狀態(tài)的關(guān)聯(lián)共軛熱力學(xué)力：

1.1.2 孔隙擴(kuò)展準(zhǔn)則的建立 基于不可逆轉(zhuǎn)的熱力學(xué)過程，孔隙擴(kuò)展需要一個適當(dāng)?shù)臏?zhǔn)則來決定。準(zhǔn)則可以是損傷共軛力(見式)的函數(shù)。但這種準(zhǔn)則難以通過試驗測定，因而現(xiàn)實多采納以物理學(xué)為依據(jù)的方法。由于微裂紋的形成和張拉應(yīng)變之間存在固然聯(lián)系，引入等效張拉應(yīng)變的定義，有其中，iε代表i方向的主應(yīng)變。括號表示只有x為正值才會被考慮，若x為負(fù)值或零，則〈x〉=0。將孔隙擴(kuò)展驅(qū)動力定義為等效張拉應(yīng)變加載歷史中達(dá)到的最大值：Yω=max(εeq，Yhis)。Yω取等效張拉應(yīng)變和加載歷史中達(dá)到的最高值Yhis之間的最大值。這表明微裂紋的形成和擴(kuò)展是一個不可逆的過程。

類似于混凝土材料中的Mazars模型，采用指數(shù)形式描述受拉和受壓狀態(tài)下微裂紋的演化過程：

式中，參數(shù)Bc和Bt用來控制在壓縮和張拉荷載條件下壓實黏土微裂紋的演化速度，很明顯Bt＞Bc，即在壓實黏土中，受拉條件下微裂紋更容易形成和擴(kuò)展。Yc0和Yt0定義為微裂紋的初始閾值。此外，ωc和ωt的演化法則應(yīng)當(dāng)驗證中的能量耗散條件。

1.1.3 基質(zhì)吸力對壓實黏土孔隙演化的影響 在非飽和土中，還要考慮到負(fù)空壓引起的基質(zhì)吸力對材料力學(xué)性能及微裂紋發(fā)展的影響。在非飽和條件下，線彈性孔隙介質(zhì)的本構(gòu)關(guān)系可以表達(dá)為：

式中π是等效孔隙壓力，b是比奧系數(shù)，由多孔介質(zhì)理論可以給出：

式中，Kb是材料的不排水體積模量，Km是矩陣體積模量，E是楊氏模量，v是泊松比。楊氏模量E通過函數(shù)E(ω)給出，說明比奧系數(shù)也受到孔隙擴(kuò)展的影響。對于壓實黏土，b可以直接取1。非飽和土中，基質(zhì)吸力對壓實粘土的力學(xué)性質(zhì)及微裂紋的擴(kuò)展有顯著影響，因此需要將其影響考慮進(jìn)來。如在土樣失水干燥過程中，可以觀察到干燥裂縫。這主要是由于土樣固有的非均質(zhì)性和結(jié)構(gòu)上的影響。目前工作中只認(rèn)為干燥裂縫是由于結(jié)構(gòu)的影響。為此，在損傷演化中，要考慮到由于受限制的收縮變形而引起的局部張拉應(yīng)力。提出了等效張拉應(yīng)變的一個新形式?εeq，該形式將力學(xué)荷載和基質(zhì)荷載引起的應(yīng)變結(jié)合起來。由于基質(zhì)荷載而引起的應(yīng)變計算公式如下：

式第二項為基質(zhì)吸力產(chǎn)生的等效拉應(yīng)變。很明顯，基質(zhì)吸力引起的等效拉應(yīng)變與等效孔隙壓力π、比奧系數(shù)b、土體的體積模量K0有直接關(guān)系，取決于飽和度，因而?εeq受土樣吸水飽和或失水干燥過程的控制。在等效應(yīng)變的新形式中，式中的損傷驅(qū)動力是?εeq而不是εeq的函數(shù)。

1.1.4 模型參數(shù)的確定 在提出的壓實黏土的力學(xué)模型中，有六個參數(shù)需要確定：兩個未損傷材料的彈性參數(shù)(K0,G0)，及四個損傷特征參數(shù)(Bc，Bt，Yt0，Yc0)。參數(shù)確定方法如下：初始彈性常數(shù)，即初始體積模量K0和剪切模量G0，與初始楊氏模量E0和泊松比v0相關(guān)的，可從初始損傷前單軸張拉(或單軸壓縮)試驗應(yīng)力應(yīng)變曲線的線性部分獲得。與張拉損傷有關(guān)的參數(shù)Bt和Yt0和與壓縮損傷有關(guān)的參數(shù)Bc和Yc0，可分別通過拉伸和壓縮試驗確定。

1.2 壓實黏土的裂縫本構(gòu)模型

在壓實黏土未形成宏觀裂紋時，假定為均質(zhì)，此時采用1.1中建立的黏土孔隙擴(kuò)展模型。在圖1中，應(yīng)力達(dá)到峰值前，采用上述模型，這是典型的非線性過程。當(dāng)孔隙發(fā)展到臨界值，即圖中標(biāo)注的ω0時，應(yīng)力也即將達(dá)到峰值，宏觀裂縫出現(xiàn)。此時連續(xù)介質(zhì)力學(xué)不再適用，需要將不連續(xù)體和宏觀裂紋單獨考慮。

在出現(xiàn)宏觀裂紋的應(yīng)力軟化階段需要引入軟化曲線。以試驗測量結(jié)果為依據(jù)，給出了一個經(jīng)驗公式如下：

式中，ft是材料的單軸抗拉強(qiáng)度，[un]m是對應(yīng)應(yīng)力為0的最大裂縫寬度，c1是模型參數(shù)，用來控制應(yīng)力—裂紋寬度曲線的形狀。實際中，由于裂紋展開而耗散的斷裂能Gf可以直接通過峰后曲線來獲得，如圖2所示，其大小代表了曲線下的面積。上述三個參數(shù)表述了裂紋出現(xiàn)后的力學(xué)行為特征，可通過單軸張拉試驗確定。

圖1 巖土材料單軸拉伸應(yīng)力應(yīng)變?nèi)€及孔隙演化過程Fig.1 Complete stress-strain curve of geotechnical material under uniaxial tension and its pore evolution process

圖2 巖土材料宏觀裂縫的應(yīng)力—裂縫寬度關(guān)系Fig.2 Relationship between stress and crack width of macroscopic crack of geotechnical material

連續(xù)介質(zhì)力學(xué)和宏觀裂隙的力學(xué)行為通過臨界損傷值緊密地聯(lián)系在一起，形成一個完整的分析，得到了連續(xù)的應(yīng)力應(yīng)變?nèi)€。如圖3所示，在上升曲線部分，為連續(xù)介質(zhì)力學(xué)的分析結(jié)果，下降曲線則是運(yùn)用了裂縫張開模型分析的結(jié)果。

圖3 巖土材料單軸拉伸試驗曲線及其數(shù)值模擬Fig.3 Uniaxial tensile test curve of geotechnical material and its numerical simulation

圖4 帶有裂紋的不連續(xù)體Fig.4 Discontinuity with crack

2 裂縫擴(kuò)展數(shù)值方法

2.1XEFM方法

擴(kuò)展有限元的基本觀點是“擴(kuò)展”被裂縫分割的單元節(jié)點自由度?！皵U(kuò)展”的自由度用來描述不連續(xù)體的“位移間斷”。

2.1.1 XEFM位移模式 擴(kuò)展有限元通過對裂縫分割單元結(jié)點增加額外自由度來表征不連續(xù)體，位移模式的一般形式：

式中Ω為帶有裂紋的面積域(見圖4)，Γ為裂紋，Ni為單元插值函數(shù)，ai為單元節(jié)點位移，Hj(x)為Heaviside函數(shù)，在裂紋面的上下兩側(cè)分別取值+1和-1，作為富集函數(shù)用以表示裂紋處位移的間斷，

bi為“擴(kuò)展”自由度的節(jié)點位移。式中第一項即為有限元的公式，而第二項為擴(kuò)展自由度的貢獻(xiàn)。

2.1.2 XEFM支配方程 利用虛功原理，建立擴(kuò)展有限元法的支配方程：

式中，da是常規(guī)有限元節(jié)點位移增量，db是擴(kuò)展自由度的節(jié)點位移增量，為總體勁度矩陣，fext為外界施加在節(jié)點上的節(jié)點荷載，為由內(nèi)力引起的常規(guī)自由度和擴(kuò)展自由度的節(jié)點荷載列陣。實際上，上述支配方程可以簡寫為KU=F，這一表達(dá)和有限元法完全一致。

總體勁度矩陣由Kaa、Kab、Kba和Kbb四個單元勁度矩陣構(gòu)成，其中與常規(guī)有限元單元勁度矩陣一致；考慮了常規(guī)和擴(kuò)展自由度之間的耦合；描述了不連續(xù)體對勁度矩陣的貢獻(xiàn)，其中T為裂縫的勁度矩陣，為下標(biāo)n和t代表裂縫的法線和切線方向。

2.1.3 裂紋積分方法 被裂縫分割的單元，由于跳躍方程在整個方程區(qū)域內(nèi)是不連續(xù)的，故在計算單元剛度矩陣和單元節(jié)點力時不像有限元法中那樣進(jìn)行Gauss積分，而是需要進(jìn)行分區(qū)處理。被裂縫分割以后，分區(qū)則是不規(guī)則的四邊形或三角形。為了解決這一問題，需將裂縫兩邊的區(qū)域全部離散為三角形。值得注意的是，這些三角形并不是新生成的單元，而是積分過程中使用的積分域，不形成額外的節(jié)點和單元。為了考慮裂縫對單元剛度矩陣和節(jié)點荷載的貢獻(xiàn)，裂縫表面也進(jìn)行了線積分。

2.2 裂縫擴(kuò)展準(zhǔn)則和擴(kuò)展方向

關(guān)于黏土開裂，目前常用的準(zhǔn)則主要有拉應(yīng)力開裂準(zhǔn)則和剪切開裂準(zhǔn)則。拉應(yīng)力開裂準(zhǔn)則，即當(dāng)土體中拉主應(yīng)力達(dá)到土體抗拉強(qiáng)度時，土體開裂，裂縫的方向垂直于該拉主應(yīng)力方向。該準(zhǔn)則物理意義明確，但土石壩心墻內(nèi)部通常都處于受壓狀態(tài)，并不存在拉應(yīng)力。很多學(xué)者都認(rèn)為心墻發(fā)生開裂是由于其滿足拉應(yīng)力開裂準(zhǔn)則[13-15]，但即便存在“拱效應(yīng)”，也很難使心墻的應(yīng)力為零或負(fù)值[16]。剪切開裂準(zhǔn)則即莫爾—庫倫強(qiáng)度理論。在這里，還可以采用拉應(yīng)變開裂準(zhǔn)則，即當(dāng)土中拉應(yīng)變達(dá)到土體抗拉強(qiáng)度對應(yīng)的拉應(yīng)變時，土體開裂。相對于拉應(yīng)力準(zhǔn)則，對于心墻土體，拉應(yīng)變準(zhǔn)則更加適合，因為即便土體為壓應(yīng)力的條件下，拉應(yīng)變也可能存在。

在本文計算中，極限拉應(yīng)變、拉應(yīng)力、抗剪強(qiáng)度、臨界損傷同時作為裂縫開裂和擴(kuò)展的準(zhǔn)則。在每一個加載步后，那些沒有被宏觀裂紋分割的單元內(nèi)的每個節(jié)點的損傷值都會被校核并和開裂準(zhǔn)則比較。若開裂準(zhǔn)則滿足，那么有兩種被區(qū)分的情況。第一種情況，若該單元在裂紋尖端范圍內(nèi)，這屬于裂紋擴(kuò)展情況。另一種情況，若該單元周圍是普通的單元，那么在普通單元內(nèi)會產(chǎn)生新的裂紋尖端。為了便于計算，允許裂紋擴(kuò)展以確保尖端總在單元的邊界上。同時，將富集函數(shù)添加到普通單元節(jié)點上以描述單元的位移間斷。隨著新裂紋的形成，在當(dāng)前加載步下，整個結(jié)構(gòu)的平衡會重新進(jìn)行計算和校核，直到達(dá)到新的平衡。

裂縫擴(kuò)展方向的確定是宏觀裂紋擴(kuò)展模型中最難解決的問題之一。在本文，采用了以平均張拉應(yīng)力為基礎(chǔ)的準(zhǔn)則。宏觀裂紋從它的尖端以垂直最大張拉應(yīng)力的方向進(jìn)行擴(kuò)展?？紤]到局部相互作用的影響，后期通過非局部平均應(yīng)力張量來確定。

式中，w是Gauss型權(quán)函數(shù)，定義如下：

式中，r代表當(dāng)前點和裂紋尖端之間的距離，lc是材料的特征長度，該特征長度決定了裂紋尖端周圍相互作用域的范圍。lc的數(shù)值取決于特征元的尺寸，在二維問題中，可以表示為：其中，Ae是二維網(wǎng)格中的平均單元面積。

3 數(shù)值模擬算例

工程背景為蒙陰縣岸堤水庫大壩，壩高29.80 m，頂寬7 m，為混凝土路面。主河槽壩段為粘土心墻砂殼壩。該壩修建于上世紀(jì)50年代末，當(dāng)時建設(shè)標(biāo)準(zhǔn)偏低，施工質(zhì)量較差。心墻局部土料含砂粒較高，滲透系數(shù)偏大，且土料整體壓實不均勻，因而存在水力劈裂破壞的風(fēng)險。

模擬試樣為長20 cm高5 cm的長方形，其邊界條件如圖5所示。為了模擬土石壩心墻的邊界條件，限制水平面前后方向的側(cè)向變形，實際上就是一個平面應(yīng)變問題，在試件的上部施加恒定的初始壓力。在試件的左右兩邊水平位移全部限定，試件下部的豎向位移限定。試件的左邊施加隨時間變化的水壓力，而試件右邊水壓力一直設(shè)定為0，試件的上部和下部皆為不透水條件。在試件中央設(shè)定一條預(yù)設(shè)的裂縫，其寬度為0.5 cm，長度為10 cm。庫水壓力的大小和變化可通過調(diào)整水壓力的大小來模擬。

整個試件被劃分1600個單元，1701個節(jié)點（圖6），預(yù)設(shè)裂縫為藍(lán)色區(qū)域部分，用于模擬土石壩中的軟弱滲透面，豎向應(yīng)力為0.5 MPa，左邊的水壓力在0.05 s內(nèi)升高到0.2 MPa，然后保持這一水壓力直到試驗結(jié)束。土樣力學(xué)參數(shù)通過三軸試驗，并參考大壩試驗報告相關(guān)成果綜合確定，其中粘土的滲透系數(shù)為3.99×10-9m/s，砂土部分的滲透系數(shù)為1×10-3m/s。因為土樣一直處于恒定壓力狀態(tài)下，采用開裂準(zhǔn)則為拉應(yīng)變準(zhǔn)則。在裂縫擴(kuò)展過程中，裂縫穿越單元的滲透系數(shù)也隨之增加，采用和砂土一致的滲透系數(shù)。這是一個流固耦合加裂縫擴(kuò)展的復(fù)雜過程。

圖5 水力劈裂模擬試件邊界條件Fig.5 Boundary condition of hydraulic fracture simulation specimen

圖6 試件網(wǎng)格劃分Fig.6 Specimen mesh generation

圖7 試件不同階段水壓力分布及裂縫擴(kuò)展示意圖Fig.7 Distribution of hydraulic pressure and crack extension of specimen in different stages

圖7 給出了對應(yīng)于不同階段的試件的水壓力分布和變形網(wǎng)格豎向位移的分布?？梢钥闯?，在0.5 MPa的豎向壓力下，在0.2 MPa的水壓力沖擊下，土樣被劈裂。在0.05 s時，由于土樣的滲透率較低，試件內(nèi)部的水壓力還基本保留為初始水壓力。預(yù)設(shè)滲透弱面的滲透系數(shù)遠(yuǎn)遠(yuǎn)大于試件的滲透系數(shù)，其水壓力在加壓的同時達(dá)到了0.2 MPa。在這一水壓力沖擊下，預(yù)設(shè)的滲透弱面被水壓力沖開，呈張開趨勢，如圖7(b.1)示意。使得裂縫前端的單元內(nèi)部產(chǎn)生較大的拉應(yīng)變，雖然在該時刻，試件仍然處于受壓狀態(tài)，但是預(yù)設(shè)裂縫的前部單元的拉應(yīng)變已經(jīng)達(dá)到了裂縫開裂準(zhǔn)則規(guī)定的拉應(yīng)變，裂縫開始擴(kuò)展。

在裂縫擴(kuò)展過程中，由于裂縫的張開，水壓沿著裂縫向試件深部，即試件右邊擴(kuò)展，在水壓力作用下，裂縫被撐開，使得裂縫前部單元一直處于拉伸狀態(tài)，因此裂縫一直快速向前擴(kuò)展，如圖7(a.2)和(b.2)示意。在裂縫即將貫通時，裂縫被撐開達(dá)到最大程度。在裂縫貫通的一瞬間，由于試件右側(cè)處于排水狀態(tài)，裂縫中的水壓力，特別是試件右邊，水壓力立即降到0。對應(yīng)的裂縫同時閉合，隨著時間的增加，裂縫中的水壓力基本保持不變，相應(yīng)的裂縫開度變化也不大。但是由于裂縫的擴(kuò)展，土樣與水的接觸面加大，土樣中的水壓力很快達(dá)到穩(wěn)定狀態(tài)。

此算例雖然簡單，但卻驗證了水力劈裂發(fā)生的兩個基本條件，即“缺陷”和“快速蓄水的初期”[17,18]。缺陷即指心墻中應(yīng)該存在薄弱面，算例中用滲透系數(shù)較大的單元來表示；而快速蓄水指需要一定的水力梯度，心墻內(nèi)部的水壓力本來很小，幾乎為零，0.2 MPa的水頭很快地進(jìn)入到薄弱面中，才導(dǎo)致了裂縫的擴(kuò)展，從而發(fā)生了水力劈裂。從圖7中還可以看到，隨著水力劈裂的發(fā)展，裂縫周圍土樣內(nèi)部的孔隙水應(yīng)力不斷升高，即紅色區(qū)域不斷增大，這在一定程度上抵消裂縫和心墻的水頭差，從而減緩裂縫擴(kuò)展程度。當(dāng)裂縫貫穿后，裂縫內(nèi)的高水壓力突然消失，裂縫閉合，土樣內(nèi)的水壓力重新分布并達(dá)到穩(wěn)定。正是由于水力劈裂后裂縫的快速閉合，才導(dǎo)致實際中的水力劈裂現(xiàn)象難以觀察。算例不僅可以反映裂縫連續(xù)的擴(kuò)展?fàn)顟B(tài)，還能反映裂縫周圍孔隙水應(yīng)力的連續(xù)變化，這是擴(kuò)展有限元的優(yōu)越性。

4 結(jié)論

本文建立了孔隙擴(kuò)展模型和裂縫張開模型，并用擴(kuò)展有限元法模擬水力劈裂過程，其模型參數(shù)可通過試驗測得，且模擬過程嚴(yán)格服從熱力學(xué)的不可逆定律。相比文獻(xiàn)[4]，考慮了心墻中裂縫周圍水壓力隨裂縫擴(kuò)展過程的變化，也是一種進(jìn)步。本方法對水力劈裂的數(shù)值研究提供了一種思路，但研究成果還有待室內(nèi)試驗和工程實例的驗證。

[1]張丙印,于玉貞,張建民.高土石壩的若干關(guān)鍵技術(shù)問題[C].北京：中國土木工程學(xué)會第九屆土力學(xué)及巖土工程學(xué)術(shù)會議,2003：163-186

[2]郭沖,邱秀梅,趙曉龍,等.土石壩心墻填土恒壓吸濕變形試驗研究[J].山東農(nóng)業(yè)大學(xué)學(xué)報：自然科學(xué)版,2013,44(1)：76-80

[3]李全明,張丙印,于玉貞,等.土石壩水力劈裂發(fā)生過程的有限元數(shù)值模擬[J].巖土工程學(xué)報,2007,29(2)：212-217

[4]王俊杰,朱俊高.堆石壩心墻抗水力劈裂性能研究[J].巖石力學(xué)與工程學(xué)報,2007,26(增1)：2880-2886

[5]李錄賢,王鐵軍.擴(kuò)展有限元法XFEM及其應(yīng)用[J].力學(xué)進(jìn)展,2005,35(1)：5-20

[6]Belytschko T,Black T.Elastic crack growth in finite elements with minimal remeshing[J].International Journal for Numerical Methods in Engineering,1999,45：601-620

[7]Mo?s N,Dolbow J,Belytschko T.A finite element method for crack growth without remeshing[J].International Journal for Numerical Methods in Engineering,1999,46：131-150

[8]郭歷倫,陳忠富,羅景潤,等.擴(kuò)展有限元方法及應(yīng)用綜述[J].力學(xué)季刊,2011,32(4)：612-625

[9]杜效鵠,段云嶺,王光綸.重力壩斷裂數(shù)值分析研究[J].水利學(xué)報,2005,36(9)：1035-1042

[10]方修君,金 峰,王進(jìn)廷.基于擴(kuò)展有限元法的Koyna重力壩地震開裂過程模擬[J].清華大學(xué)學(xué)報：自然科學(xué)版,2008,48(12)：2065-2069

[11]方修君,金 峰.裂隙水流與混凝土開裂相互作用的耦合模型[J].水利學(xué)報,2007,38(12)：1466-1474

[12]董玉文,任青文.重力壩水力劈裂分析的擴(kuò)展有限元法[J].水利學(xué)報,2011,42(11)：1361-1367

[13]黃文熙.對土石壩科研工作的幾點看法[J].水利水電技術(shù),1982(4)：23-27

[14]孫亞平.水力劈裂機(jī)理研究[D].北京：清華大學(xué),1985

[15]馮曉瑩,徐澤平.心墻水力劈裂機(jī)理的離心模型試驗研究[J].水利學(xué)報,2009,40(10)：1259-1263

[16]王俊杰,朱俊高,張 輝.關(guān)于土石壩心墻水力劈裂研究的一些思考[J].巖石力學(xué)與工程學(xué)報,2005,24(增2)：5664-5668

[17]朱俊高,王俊杰,張 輝.土石壩心墻水力劈裂機(jī)制研究[J].巖土力學(xué),2007,28(3)：487-492

[18]趙曉龍,邱秀梅,韓慧敏,等.土石壩帶裂縫黏土心墻破壞機(jī)理試驗研究[J].中國農(nóng)村水利水電,2016(2)：134-138