網(wǎng)格劃分對(duì)邊坡穩(wěn)定有限元分析的影響

李永亮，余成華,馬佐成,張殿義

(1.長(zhǎng)安大學(xué) 公路學(xué)院，陜西 西安 710064；2.深圳市勘察研究院有限公司，廣東 深圳 518026)

?

網(wǎng)格劃分對(duì)邊坡穩(wěn)定有限元分析的影響

李永亮1，余成華2,馬佐成1,張殿義1

(1.長(zhǎng)安大學(xué) 公路學(xué)院，陜西 西安 710064；2.深圳市勘察研究院有限公司，廣東 深圳 518026)

摘要：網(wǎng)格劃分對(duì)邊坡穩(wěn)定性有限元分析的準(zhǔn)確性有很大影響，而目前對(duì)此影響的研究相對(duì)不夠。以有限元強(qiáng)度折減法在邊坡穩(wěn)定性分析中的一標(biāo)準(zhǔn)算例為依據(jù)，探討ABAQUS中網(wǎng)格劃分技術(shù)、網(wǎng)格形狀、網(wǎng)格密度、局部加密網(wǎng)格、網(wǎng)格生成算法和單元積分方法對(duì)邊坡穩(wěn)定性計(jì)算的影響。最后用一復(fù)雜土坡案例分析以上諸因素對(duì)邊坡穩(wěn)定性分析的敏感性，得出在適當(dāng)范圍內(nèi)加密網(wǎng)格可節(jié)約計(jì)算時(shí)間，二次積分單元計(jì)算精度高但更耗成本等結(jié)論，對(duì)邊坡穩(wěn)定分析的網(wǎng)格劃分具有一定參考價(jià)值。

關(guān)鍵詞：ABAQUS；強(qiáng)度折減法；網(wǎng)格劃分；邊坡穩(wěn)定性

強(qiáng)度折減理論概念由Zienkiwicz等(1975年)最早提出并運(yùn)用于邊坡的穩(wěn)定性分析，之后Griffith，Dwason等對(duì)其發(fā)展做了大量工作。我國(guó)鄭穎人[1-3]、趙尚毅[1-2]、宋二祥[4]、鄭宏[5]等擴(kuò)大了有限元強(qiáng)度折減法的應(yīng)用范圍。不同流動(dòng)法則、屈服準(zhǔn)則、計(jì)算模型的范圍和計(jì)算參數(shù)對(duì)強(qiáng)度折減法計(jì)算精度的影響已得到大量研究[5-7]，而對(duì)網(wǎng)格劃分影響的研究相對(duì)不夠。因此，探究網(wǎng)格劃分對(duì)邊坡穩(wěn)定性分析的影響極有必要。

1　有限元強(qiáng)度折減法及邊坡整體失穩(wěn)的判據(jù)

邊坡穩(wěn)定分析的有限元強(qiáng)度折減法是通過(guò)不斷降低巖土體抗剪切強(qiáng)度參數(shù)直至達(dá)到極限破壞狀態(tài)為止，程序自動(dòng)根據(jù)彈塑性有限元計(jì)算結(jié)果得到滑動(dòng)破壞面，同時(shí)得到邊坡的強(qiáng)度儲(chǔ)備安全系數(shù)[8]。折減后的抗剪強(qiáng)度參數(shù)表達(dá)為:

(1)

(2)

式中：c和φ是土體所能夠提供的抗剪強(qiáng)度；cm和φm是維持平衡所需要的或土體實(shí)際發(fā)揮的抗剪強(qiáng)度；Fr是強(qiáng)度折減系數(shù)。ABAQUS通過(guò)將Fr定義為場(chǎng)變量，初始分析時(shí)為避免模型破壞，可取Fr<1，在后續(xù)分析步中線性增加場(chǎng)變量Fr，在坡體達(dá)到臨界狀態(tài)以前，每個(gè)折減系數(shù)對(duì)應(yīng)著一個(gè)塑性區(qū)發(fā)展?fàn)顩r，計(jì)算中止(數(shù)值不收斂)后按照失穩(wěn)評(píng)價(jià)標(biāo)準(zhǔn)確定安全系數(shù)[9]。

采用有限元法計(jì)算邊坡穩(wěn)定時(shí)判斷邊坡失穩(wěn)的判據(jù)有以下幾種類型[10]：

(1)以有限元計(jì)算不收斂為標(biāo)準(zhǔn)。其認(rèn)為當(dāng)邊坡處于極限平衡狀態(tài)時(shí)，非線性有限元方程組的迭代過(guò)程將不收斂(內(nèi)含迭代次數(shù)或殘差應(yīng)力與外力比值)。

(2)以廣義剪應(yīng)變標(biāo)準(zhǔn)或廣義塑性應(yīng)變標(biāo)準(zhǔn)為標(biāo)準(zhǔn)。其認(rèn)為當(dāng)邊坡達(dá)到極限平衡狀態(tài)時(shí)，廣義剪應(yīng)變或廣義塑性應(yīng)變的等值線圖中，必有一條等值線(通常是數(shù)值較小的那條等值線)由坡底貫通到坡頂。

(3)以位移為標(biāo)準(zhǔn)。如坡頂位移增量與折減系數(shù)增量之比、殘差位移、無(wú)量綱位移、特征部位位移或位移的突變等。

(4)以塑性區(qū)的貫通為標(biāo)準(zhǔn)，可由塑性應(yīng)變及塑性功表示。

(5)聯(lián)合使用(1)～(4)中的某幾條。

2　計(jì)算模型與參數(shù)

2.1力學(xué)模型和屈服準(zhǔn)則

采用理想彈塑性本構(gòu)模型，屈服準(zhǔn)則采用Mohr-Coulomb強(qiáng)度準(zhǔn)則，屈服條件為：

(3)

式中：θα為洛得應(yīng)力角；I1為應(yīng)力張量的第一不變量；J2為偏應(yīng)力張量的第二不變量[11]。

2.2邊坡模型

采用參考文獻(xiàn)[2]中鄭穎人等的一個(gè)算例，邊坡模型如圖1所示。

圖1邊坡模型圖

2.3邊坡模型網(wǎng)格劃分

根據(jù)ABAQUS的網(wǎng)格劃分技術(shù)，四邊形結(jié)構(gòu)化網(wǎng)格如圖2所示。

2.4土體參數(shù)

土體的本構(gòu)關(guān)系采用Mohr-Coulomb彈塑性模型，其參數(shù)見(jiàn)表1。

圖2邊坡網(wǎng)格劃分圖

表1　土體物理力學(xué)參數(shù)

3　網(wǎng)格劃分計(jì)算及結(jié)果分析

本文所采用的邊坡算例，鄭穎人等得出其安全系數(shù)為1.21，Spencer法算得為1.20。以下分析網(wǎng)格劃分對(duì)該邊坡穩(wěn)定性計(jì)算的影響。

3.1網(wǎng)格劃分技術(shù)的影響

3.1.1四邊形網(wǎng)格下不同劃分技術(shù)

計(jì)算結(jié)果見(jiàn)表2。

表2　四邊形網(wǎng)格下的計(jì)算結(jié)果

注：表中網(wǎng)格數(shù)均是在相同的網(wǎng)格種子布置下而得。

Structured技術(shù)是將一些標(biāo)準(zhǔn)的單元(如四邊形、正方體)等應(yīng)用于一些簡(jiǎn)單的沒(méi)有獨(dú)立的點(diǎn)、線、面和洞的幾何區(qū)域。Sweep技術(shù)是先在源邊或面上生成網(wǎng)格，然后沿著掃略路徑復(fù)制節(jié)點(diǎn)，直到目標(biāo)邊或面，得到網(wǎng)格。Free技術(shù)是最為靈活的網(wǎng)格劃分技術(shù)，幾乎可以用于任何形狀[9]。

由上表可知Structured技術(shù)和Free技術(shù)劃得的網(wǎng)格數(shù)少于Sweep技術(shù)，并且其安全系數(shù)、總位移和剪應(yīng)力值均相同而較Sweep技術(shù)的結(jié)果精確一些。

3.1.2三角形網(wǎng)格下不同劃分技術(shù)

計(jì)算結(jié)果見(jiàn)表3。

表3　三角形網(wǎng)格下的計(jì)算結(jié)果

由表3可知三角形網(wǎng)格下Structured技術(shù)和Free技術(shù)所得網(wǎng)格數(shù)、安全系數(shù)等值相同，此時(shí)兩種劃分技術(shù)所得計(jì)算精度和耗時(shí)一樣。

3.2網(wǎng)格密度的影響

3.2.1四邊形網(wǎng)格Structured技術(shù)下不同網(wǎng)格密度

計(jì)算結(jié)果見(jiàn)表4。

表4　四邊形網(wǎng)格Structured技術(shù)時(shí)的計(jì)算結(jié)果

根據(jù)表4可知當(dāng)網(wǎng)格數(shù)為3328和2726時(shí)所得的安全系數(shù)和剪應(yīng)力都基本不再變化，說(shuō)明網(wǎng)格密度已經(jīng)足夠[12]。網(wǎng)格數(shù)為2726時(shí)邊坡總位移云圖如圖3所示。

圖3網(wǎng)格數(shù)為2726時(shí)的總位移云圖

3.2.2四邊形網(wǎng)格Sweep技術(shù)下不同網(wǎng)格密度

計(jì)算結(jié)果如表5。

表5　四邊形網(wǎng)格Sweep技術(shù)時(shí)的計(jì)算結(jié)果

類似于表4，當(dāng)網(wǎng)格數(shù)為2766和2532時(shí)Fr和剪應(yīng)力值變化很小，說(shuō)明此時(shí)網(wǎng)格數(shù)目已能滿足計(jì)算精度要求。將表5與表4比較可知，同樣條件下Sweep技術(shù)比Structured技術(shù)劃出的網(wǎng)格偏少而又能保證計(jì)算結(jié)果的精度且計(jì)算歷時(shí)相對(duì)較小。

3.3網(wǎng)格形狀的影響

3.3.1Structured技術(shù)下不同網(wǎng)格形狀

計(jì)算結(jié)果見(jiàn)表6。由表6可知，在Structured技術(shù)下四邊形網(wǎng)格數(shù)目最少，而計(jì)算得安全系數(shù)較其他兩者精確，且其總位移和剪應(yīng)力值介于后兩者值之間，所以此類情況應(yīng)優(yōu)先考慮采用四邊形形狀網(wǎng)格。

表6　Structured技術(shù)下不同網(wǎng)格形狀計(jì)算結(jié)果

3.3.2Sweep技術(shù)下不同網(wǎng)格形狀

計(jì)算結(jié)果見(jiàn)表7。

表7　Sweep技術(shù)下不同網(wǎng)格形狀計(jì)算結(jié)果

由表7可知，這兩種形狀的網(wǎng)格此時(shí)所得的網(wǎng)格數(shù)、安全系數(shù)、總位移和應(yīng)力值均相同。

3.4局部加密網(wǎng)格的影響

根據(jù)有限元計(jì)算的特點(diǎn)并結(jié)合該邊坡模型，在坡頂之后一定范圍內(nèi)加密網(wǎng)格，計(jì)算結(jié)果見(jiàn)表8。

表8　網(wǎng)格局部加密的計(jì)算結(jié)果

注：H表示坡高，本文H=20 m。

通過(guò)表8與表4的比較可知，在1.5H和1.0H范圍內(nèi)加密網(wǎng)格，在模型其他部位劃分相對(duì)稀疏的網(wǎng)格既能保證計(jì)算的精確性，又能減少有限元迭代計(jì)算的總次數(shù)而節(jié)約計(jì)算時(shí)間。而在0.5H范圍內(nèi)的加密計(jì)算結(jié)果并不理想，故建議在進(jìn)行分析時(shí)可適當(dāng)局部細(xì)化網(wǎng)格來(lái)提高計(jì)算效益。

3.5網(wǎng)格生成算法的影響

在Free技術(shù)下不同的網(wǎng)格生成方法計(jì)算結(jié)果如表9。

表9　不同網(wǎng)格生成算法的計(jì)算結(jié)果

中性軸算法生成的網(wǎng)格與種子吻合程度較差，進(jìn)階算法生成的網(wǎng)格和種子可以較好地吻合，更容易得到大小均勻的網(wǎng)格。由表9可知進(jìn)階算法所得安全系數(shù)和應(yīng)力值等較中性軸算法精確。

3.6單元積分方法的影響

不同單元積分方法的計(jì)算結(jié)果見(jiàn)表10。

表10　不同單元積分方法的計(jì)算結(jié)果

由表10可知Quadratic對(duì)安全系數(shù)和應(yīng)力等的計(jì)算結(jié)果較其他兩種積分方法準(zhǔn)確，但其計(jì)算耗時(shí)大。線性縮減積分限于其自身特點(diǎn)，所得數(shù)值不夠精確。

4　案例分析

為了探究網(wǎng)格劃分諸因素對(duì)邊坡穩(wěn)定性分析的敏感性，采用一復(fù)雜邊坡案例對(duì)其進(jìn)行分析[8]。

4.1邊坡模型如圖4所示。

圖4邊坡模型圖

4.2土體物理力學(xué)參數(shù)

表11　土體物理力學(xué)參數(shù)

4.3不同網(wǎng)格劃分因素下的邊坡穩(wěn)定性分析

計(jì)算結(jié)果見(jiàn)表12。由表12可知當(dāng)網(wǎng)格形狀和密度種子一定，不同的劃分技術(shù)所得邊坡安全系數(shù)差異很小。當(dāng)網(wǎng)格形狀和劃分技術(shù)相同而網(wǎng)格密度不同時(shí)，分析得出的安全系數(shù)數(shù)值變化并不明顯。在劃分技術(shù)一定而網(wǎng)格形狀不同時(shí)，Structured技術(shù)下三角形網(wǎng)格在同樣的網(wǎng)格密度下所得結(jié)果不理想而且計(jì)算耗時(shí)多。在坡頂后適當(dāng)區(qū)域內(nèi)加密網(wǎng)格求得安全系數(shù)數(shù)值變動(dòng)不大，只是起到相對(duì)節(jié)約計(jì)算成本的作用，而對(duì)計(jì)算結(jié)果的精度影響不顯著。同樣的網(wǎng)格密度下進(jìn)階算法所得安全系數(shù)相對(duì)較中性軸算法的結(jié)果準(zhǔn)確，但若進(jìn)行調(diào)整會(huì)發(fā)現(xiàn)二者分析得出的結(jié)果精確程度相同。不同的單元積分方法對(duì)計(jì)算結(jié)果有很大影響，二次型積分方法求得的安全系數(shù)較線性積分精確，而相同情況下縮減積分的結(jié)果明顯偏大。因此，對(duì)于此邊坡，網(wǎng)格劃分諸因素中單元積分方法對(duì)邊坡穩(wěn)定性計(jì)算影響最大，網(wǎng)格形狀和網(wǎng)格生成方法次之，最后是網(wǎng)格劃分技術(shù)和局部區(qū)域加密網(wǎng)格。

表12　網(wǎng)格劃分因素影響程度分析結(jié)果

5　結(jié)　語(yǔ)

采用ABAQUS探討了不同的網(wǎng)格劃分因素對(duì)強(qiáng)度折減法分析邊坡穩(wěn)定性的計(jì)算結(jié)果及敏感性，得出如下結(jié)論：

(1)在四邊形和三角形網(wǎng)格下Structured技術(shù)劃分的網(wǎng)格數(shù)相對(duì)其他技術(shù)所得網(wǎng)格較少，求出安全系數(shù)和應(yīng)力值等的精確程度和計(jì)算時(shí)間優(yōu)于Sweep技術(shù)和Free技術(shù)。

(2)在Structured技術(shù)和Sweep技術(shù)下對(duì)模型劃分不同密度的網(wǎng)格時(shí)當(dāng)求出的安全系數(shù)、最大位移和應(yīng)力值等趨于穩(wěn)定時(shí)說(shuō)明計(jì)算結(jié)果已滿足要求。實(shí)際建模中，應(yīng)在保證分析精度的前提下，選擇適當(dāng)?shù)木W(wǎng)格密度，對(duì)模型進(jìn)行局部網(wǎng)格適當(dāng)細(xì)化并盡量保證單元形狀是規(guī)則的。

(3)只要網(wǎng)格種子布置得不是過(guò)于稀疏，使用進(jìn)階算法往往更容易得到形狀規(guī)則的網(wǎng)格，而且其結(jié)果也相對(duì)較為準(zhǔn)確，但這不是絕對(duì)的，保證單元尺寸均勻十分重要。另外二次完全積分所得應(yīng)力值和安全系數(shù)更精確，但其計(jì)算成本高，分析時(shí)應(yīng)據(jù)實(shí)際情況選用最優(yōu)的單元積分方法。

(4)通過(guò)復(fù)雜邊坡案例對(duì)網(wǎng)格劃分因素的敏感性分析可知，對(duì)邊坡穩(wěn)定性計(jì)算結(jié)果影響最大的是單元積分方法，其次為網(wǎng)格形狀和網(wǎng)格生成方法次之，最后為網(wǎng)格劃分技術(shù)和局部區(qū)域加密網(wǎng)格的影響。

上述結(jié)論僅是數(shù)值模擬結(jié)果,有待于在工程實(shí)踐中進(jìn)一步驗(yàn)證。

參考文獻(xiàn)：

[1]鄭穎人,趙尚毅,張魯渝.用有限元強(qiáng)度折減法進(jìn)行邊坡穩(wěn)定分析[J].中國(guó)工程科學(xué),2002,4(10):57-61,78.

[2]鄭穎人,趙尚毅.有限元強(qiáng)度折減法在土坡與巖坡中的應(yīng)用[J].巖石力學(xué)與工程學(xué)報(bào),2004,23(19):3381-3388.

[3]劉明維,鄭穎人.基于有限元強(qiáng)度折減法確定滑坡多滑動(dòng)面方法[J].巖石力學(xué)與工程學(xué)報(bào),2006,25(8):1544-1549.

[5]鄭宏,劉德富.彈塑性矩陣Dep的特性和有限元邊坡穩(wěn)定性分析中的極限狀態(tài)標(biāo)準(zhǔn)[J].巖石力學(xué)與工程學(xué)報(bào),2005,24(7):1099-1105.

[6]李宗坤,霍香麗,李艷.本構(gòu)模型對(duì)邊坡穩(wěn)定有限元分析的影響[J].人民黃河,2009,31(12):89,91.

[7]李垠,蘇凱,李杰.Mohr-Coulomb等面積圓屈服準(zhǔn)則在邊坡穩(wěn)定性分析中的應(yīng)用[J].大地測(cè)量與地球動(dòng)力學(xué),2009,29(2):135-139.

[8]鄭穎人，等.邊坡與滑坡工程治理(第2版)[M].北京:人民交通出版社,2010.

[9]費(fèi)康,張建偉.ABAQUS在巖土工程中的應(yīng)用[M].北京:中國(guó)水利水電出版社,2010.

[10]謝新宇,劉開(kāi)富,張繼發(fā).邊坡及基礎(chǔ)工程數(shù)值分析新進(jìn)展[M].北京:科學(xué)出版社,2010.

[11]張宏.應(yīng)用彈塑性理論[M].西安:西北工業(yè)大學(xué)出版社,2011.

[12]曹金鳳,石亦平.ABAQUS有限元分析常見(jiàn)問(wèn)題解答[M].北京:機(jī)械工業(yè)出版社,2009.

TheEffectoftheMeshontheFEMAnalysisofSlopeStability

LI Yong-liang1,YU Cheng-hua2,MA Zuo-cheng1,ZHANG Dian-yi1

(1.SchoolofHighways,Chang’anUniversity,Xi’an,Shaanxi710064,China;2.ShenzhenInvestigation&ResearchInstituteCo.,Ltd.,Shenzhen,Guangdong518026,China)

Abstract：Meshing has great impact on the accuracy of the FEM analysis of slope stability,but its research is currently insufficient.By employing ABAQUS and on the basis of an standard calculation example of slope stability which was analyzed by FEM Shear Strength Reduction Method,the effects of meshing techniques,mesh shapes,mesh density,partially-refined mesh,mesh algorithm and element integration on account of slope stability were investigated.At last,the effects of above factors on the sensitivity of slope stability were analyzed on a complex soil slope.the results indicate that the increase of the refinement of the mesh in proper regions can reduce computing time,quadratic integration element method has higher calculation accuracy but with larger computational expense and so on.These findings have a certain reference value for the meshing of the analysis of slope stability.

Keywords：ABAQUS;shear strength reduction method;meshing;slope stability

DOI:10.3969/j.issn.1672-1144.2014.06.026

中圖分類號(hào)：TU413.6+2

文獻(xiàn)標(biāo)識(shí)碼：A

文章編號(hào)：1672—1144(2014)06—0133—04

作者簡(jiǎn)介：李永亮(1992—)，男，甘肅天水人，碩士研究生，研究方向?yàn)檫吰路€(wěn)定分析、加筋邊坡和加筋擋墻。

收稿日期：2014-07-18修稿日期：2014-08-21