趙杰　徐劍　王桂萱　尹訓(xùn)強(qiáng)

摘要：從巖石材料細(xì)觀(guān)非均勻性特點(diǎn)出發(fā)，運(yùn)用連續(xù)介質(zhì)損傷力學(xué)理論，推導(dǎo)了巖石材料細(xì)觀(guān)損傷本構(gòu)關(guān)系式。以ANSYS二次開(kāi)發(fā)平臺(tái)為主體，基于細(xì)觀(guān)損傷力學(xué)，提出基于宏觀(guān)均質(zhì)假定的等效損傷單元來(lái)模擬隧道洞口段巖質(zhì)邊坡的損傷開(kāi)裂過(guò)程。計(jì)算結(jié)果表明：在細(xì)觀(guān)方面邊坡?lián)p傷破壞是個(gè)不斷損傷的過(guò)程，在地震荷載作用下，邊坡結(jié)構(gòu)首先會(huì)出現(xiàn)微裂紋，隨著地震動(dòng)時(shí)間增加，已出現(xiàn)的微裂紋逐漸擴(kuò)展形成可見(jiàn)的裂縫，最后形成連續(xù)的損傷區(qū)域；尤其在隧道拱頂上方形成連續(xù)貫通的損傷區(qū)域，可看出拱頂已完全破壞，給隧道安全運(yùn)行帶來(lái)不利影響，因此在設(shè)計(jì)施工時(shí)要對(duì)拱頂上方進(jìn)行加固處理。該計(jì)算方法和分析成果可用于同類(lèi)復(fù)雜結(jié)構(gòu)的損傷破壞分析。 關(guān)鍵詞：巖質(zhì)邊坡；非均質(zhì)性；損傷力學(xué)；等效損傷單元；損傷破壞；ANSYS

中圖分類(lèi)號(hào)：U45 文獻(xiàn)標(biāo)識(shí)碼：A 文章編號(hào)：1000-0666（2016）01-0053-07

0 引言

近年來(lái)發(fā)生的重大巖石動(dòng)力災(zāi)害事件表明：礦山巷道、地下洞室、隧道工程、邊坡工程垮塌和破壞，多數(shù)是由于巖體內(nèi)的節(jié)理、裂隙及其擴(kuò)展造成（羅國(guó)煌等，2004；李明田等，2003）。巖體內(nèi)的節(jié)理、裂隙就是其損傷的實(shí)質(zhì)性表現(xiàn)。所以，從細(xì)觀(guān)損傷理論出發(fā)，研究圍巖體細(xì)觀(guān)損傷擴(kuò)展機(jī)理和演化特征對(duì)巖體工程穩(wěn)定性的評(píng)價(jià)與支護(hù)設(shè)計(jì)是非常有意義的。

Kachanov在1958年提出損傷力學(xué)思想，至今已經(jīng)過(guò)50多年的發(fā)展，一些理論已日益成熟并在各學(xué)科和工程技術(shù)領(lǐng)域得以應(yīng)用。周維垣等（1998）應(yīng)用連續(xù)介質(zhì)損傷力學(xué)理論，從巖體內(nèi)部微裂紋產(chǎn)生和擴(kuò)展的損傷機(jī)理出發(fā)，建立了堅(jiān)硬巖體的彈性一損傷耦合的各向異性彈脆性損傷本構(gòu)模型。朱萬(wàn)成等（2003）從巖石的細(xì)觀(guān)非均勻性特點(diǎn)出發(fā)，應(yīng)用彈性損傷力學(xué)對(duì)巖石材料細(xì)觀(guān)的損傷演化過(guò)程進(jìn)行了描述。趙寶友等（2003）基于損傷力學(xué)原理，通過(guò)引入損傷變量的方法詳細(xì)推導(dǎo)了ABAQUS中的損傷塑性本構(gòu)模型屈服階段及其卸載后的損傷機(jī)制和力學(xué)行為，并以某水電站地下廠(chǎng)房洞室結(jié)構(gòu)為例，進(jìn)行二維動(dòng)力時(shí)程非線(xiàn)性分析，得出強(qiáng)震下采用此損傷塑性模型是合理和必要的。隋斌等（2008）利用FLAC3D軟件模擬了地震作用下某地下洞室群圍巖的動(dòng)態(tài)響應(yīng)，采用新的劈裂破壞判據(jù)對(duì)震后可能出現(xiàn)的劈裂損傷區(qū)范圍進(jìn)行了預(yù)測(cè)。左雙英和肖明（2009）基于深埋大型地下洞室群圍巖損傷機(jī)理，對(duì)映秀灣水電站洞室群在三向地震荷載耦合作用下的動(dòng)力響應(yīng)進(jìn)行數(shù)值模擬，得出地震作用下圍巖損傷分布的特征。張志國(guó)等（2010）采用動(dòng)力時(shí)程分析法研究地下洞室地震響應(yīng)特性，結(jié)合三維彈塑性損傷有限元理論，編寫(xiě)了動(dòng)力時(shí)程計(jì)算程序，分析了漁子溪水電站地下廠(chǎng)房洞室群的圍巖穩(wěn)定。

筆者在上述研究的基礎(chǔ)上，基于巖石材料的宏觀(guān)均質(zhì)假定，運(yùn)用連續(xù)介質(zhì)彈性損傷力學(xué)模型來(lái)描述單元的材料特性，采用最大拉應(yīng)力與摩爾一庫(kù)倫破壞準(zhǔn)則作為損傷閾值，以通用有限元軟件ANSYS開(kāi)發(fā)平臺(tái)為主體，運(yùn)用UPFs的二次開(kāi)發(fā)特點(diǎn)，基于用戶(hù)自定義單元的功能將損傷單元法引入到ANSYS中，用來(lái)分析地震作用下隧道洞口段巖質(zhì)邊坡從微裂紋產(chǎn)生、擴(kuò)展到形成可見(jiàn)裂縫，直至最終的損傷過(guò)程。

1 損傷演化本構(gòu)模型

1.1 巖石材料非均質(zhì)性特征

巖石作為一種自然地質(zhì)體，由于各種因素使巖石內(nèi)部存在大量不同的微缺陷，所以巖石介質(zhì)的構(gòu)成非常復(fù)雜，通常對(duì)其進(jìn)行數(shù)學(xué)描述也是非常困難的。唐春安（1993）認(rèn)為可以利用統(tǒng)計(jì)的方法對(duì)離散后的巖石介質(zhì)進(jìn)行近似描述。建立材料細(xì)觀(guān)非均質(zhì)特性隨機(jī)模型時(shí)，所謂非均質(zhì)性是指材料的力學(xué)參數(shù)，如強(qiáng)度、彈性模量、泊松比等。在給定的空間上對(duì)模型進(jìn)行網(wǎng)格劃分，為了描述材料參數(shù)的這種空間變化，需要解決兩個(gè)問(wèn)題：第一，如何使材料參數(shù)在數(shù)值上滿(mǎn)足非均勻分布規(guī)律；第二，如何將這些不同材料參數(shù)隨機(jī)地分布到劃分的單元中去。本文將巖質(zhì)邊坡看成一個(gè)大的樣本空間，采用等效損傷單元對(duì)巖石材料進(jìn)行離散，使得離散后的單元形式具有類(lèi)似于普通有限單元的特性，并且認(rèn)為其本身為均勻、連續(xù)的介質(zhì)，然后在樣本空間內(nèi)引入概率統(tǒng)計(jì)的方法，假定各單元的材料力學(xué)性質(zhì)存在差異且服從Weibull概率統(tǒng)計(jì)分布，從而可表現(xiàn)出一定的非均勻型和離散性。楊天鴻等（2004）認(rèn)為采用Weibull隨機(jī)分布可以較好地模擬材料性質(zhì)的非均勻性，該分布可以用分布密度函數(shù)來(lái)定義：式中，μ為巖石材料隨機(jī)分布力學(xué)性質(zhì)參數(shù)（強(qiáng)度、彈性模量、泊松比等）；μ⁰為巖石材料基元體力學(xué)性質(zhì)平均值；m為巖石力學(xué)性質(zhì)分布函數(shù)的密度形狀函數(shù)，其物理意義反應(yīng)了巖石材料的均質(zhì)性；f（μ）是巖石力學(xué)性質(zhì)μ的統(tǒng)計(jì)分布函數(shù)。

式（1）很好地反映了巖石材料的細(xì)觀(guān)力學(xué)性質(zhì)非均勻性分布情況。隨著均勻性系數(shù)m的增加，基元體的力學(xué)性質(zhì)將集中于一個(gè)狹窄的范圍之內(nèi)，表明單元強(qiáng)度趨于均勻；當(dāng)m值減小時(shí)，基元體的力學(xué)性質(zhì)分布范圍變寬，表明巖石介質(zhì)的性質(zhì)趨于均勻。可見(jiàn)m值反映了統(tǒng)計(jì)模型中材料結(jié)構(gòu)的離散程度，能夠較為真實(shí)地反映巖石材料的非均質(zhì)性。

1.2 巖體本構(gòu)關(guān)系與屈服準(zhǔn)則

余大慶（1993）認(rèn)為巖石應(yīng)力-應(yīng)變曲線(xiàn)是由于其受力后的不斷損傷引起微裂紋萌生和擴(kuò)展而造成的，而不是由于其塑性變形。尤其是在拉伸應(yīng)力作用下，其脆性更加明顯。因此，用彈性損傷力學(xué)本構(gòu)模型來(lái)描述巖石材料細(xì)觀(guān)單元力學(xué)性質(zhì)是合適的。在本文的研究方法中，等效單元應(yīng)力一應(yīng)變本構(gòu)關(guān)系亦采用了彈性損傷模型。在施加地震荷載前，用彈性模型和泊松比來(lái)定義每一個(gè)等效損傷單元，并且將其假定為線(xiàn)彈性、均值、無(wú)損傷的。當(dāng)單元進(jìn)入損傷階段后，依據(jù)Le-maitre應(yīng)變等價(jià)原理，認(rèn)為應(yīng)力作用在受損材料上引起的應(yīng)變與有效應(yīng)力作用在無(wú)損材料引起的應(yīng)變等價(jià)。據(jù)這一原理，受損材料單元的本構(gòu)關(guān)系可通過(guò)無(wú)損材料中的名義應(yīng)力得到：式中，E⁰為材料初始彈性模量；E為材料損傷后的彈性模量，D為損傷變量。D=0表示材料處于無(wú)損傷狀態(tài)，D=1表示處于完全損傷（斷裂或者破壞）狀態(tài)，00

在巖石材料損傷演化中，隨著等效損傷單元應(yīng)力會(huì)不斷增加，當(dāng)單元的應(yīng)力或應(yīng)變狀態(tài)滿(mǎn)足某個(gè)給定的損傷閾值（準(zhǔn)則）時(shí)，單元開(kāi)始損傷。這里選擇兩個(gè)損傷閾值準(zhǔn)則，第一個(gè)采用最大拉應(yīng)變準(zhǔn)則，當(dāng)細(xì)觀(guān)單元的最大拉伸主應(yīng)變達(dá)到其給定的應(yīng)變閥值時(shí)，該單元開(kāi)始發(fā)生拉伸損傷；第二個(gè)采用摩爾一庫(kù)侖準(zhǔn)則，當(dāng)細(xì)觀(guān)單元的應(yīng)力狀態(tài)滿(mǎn)足摩爾一庫(kù)侖準(zhǔn)則時(shí)，該單元發(fā)生剪切損傷。需要說(shuō)明的是，巖石材料在外部荷載作用下，其抗拉強(qiáng)度遠(yuǎn)低于抗壓強(qiáng)度的特點(diǎn)，且單元的拉伸損傷和剪切損傷不會(huì)同時(shí)發(fā)生于同一荷載步下。換句話(huà)說(shuō)拉伸準(zhǔn)則具有優(yōu)先權(quán)，若細(xì)觀(guān)單元滿(mǎn)足拉伸準(zhǔn)則，則不需再判斷其是否滿(mǎn)足摩爾一庫(kù)侖準(zhǔn)則。只有在損傷單元未滿(mǎn)足拉伸準(zhǔn)則之后才判定其是否滿(mǎn)足摩爾一庫(kù)侖則。圖1、2分別給出了單元材料拉伸本構(gòu)關(guān)系曲線(xiàn)與壓縮本構(gòu)關(guān)系曲線(xiàn)圖。單元滿(mǎn)足拉損傷模型以及摩爾一庫(kù)侖準(zhǔn)則對(duì)應(yīng)的損傷變量表達(dá)為式中，f_tr為單元參與強(qiáng)度；ε_t0為拉損傷應(yīng)變閾值；ε_tu為單元材料極限拉伸應(yīng)變；ε_c0為單元應(yīng)力達(dá)到單元材料單軸抗壓強(qiáng)度時(shí)的主要應(yīng)變；ε_cu為單元材料的極限應(yīng)變；ε_c為單元的單軸壓縮應(yīng)變；N為與下降段形狀有關(guān)的參數(shù)。

2 等效損傷單元在ANSYS中的實(shí)現(xiàn)

本文等效損傷單元采用的本構(gòu)關(guān)系與屈服準(zhǔn)則是其與普通有限單元最大的區(qū)別之處?；诘刃卧厥庑?，需要以ANSYS為開(kāi)發(fā)平臺(tái)，利用UPFs為開(kāi)發(fā)工具進(jìn)行二次開(kāi)發(fā)。圖3給出了基于ANSYS為二次開(kāi)發(fā)平臺(tái)的等效損傷單元建立流程，其中虛線(xiàn)方框部分由UPFs的功能實(shí)現(xiàn)。詳細(xì)的程序流程如下：

首先，建立隧道洞口段巖質(zhì)邊坡有限元模型，對(duì)用戶(hù)進(jìn)行單元編譯和鏈接，施加位移邊界條件，控制求解流程以及結(jié)構(gòu)的輸出（此部分利用GUI和APDL完成）；然后，進(jìn)入地震動(dòng)力分析求解階段，在每一荷載步施加地震荷載向量；最后，在ANSYS數(shù)據(jù)庫(kù)中將建立單元的所有信息輸入到接口子程序UserElem.f中，并完成以下步驟（此部分利用子程序接口UserElem.f完成）：

①檢驗(yàn)是否第一次進(jìn)入等效損傷單元子程序，若“是”，則利用蒙特卡羅方法生成服從Weibull分布的材料參數(shù)隨機(jī)數(shù)，然后生成最終的單元材料力學(xué)參數(shù)，將此數(shù)據(jù)存儲(chǔ)在數(shù)組中并使用FOR-TRAN語(yǔ)言中的COMMON將其定義為全局變量，以便以后的平衡迭代步和荷載步調(diào)用，進(jìn)入②；若“否”，則將直接獲得上一迭代步的材料力學(xué)參數(shù)，進(jìn)入②。

②通過(guò)節(jié)點(diǎn)真實(shí)位移{u}_t+{△u}^it計(jì)算單元的應(yīng)變s；由損傷變量Dⁱ_l（Dⁱ_lt表示拉損傷；Dⁱ_lc表示壓損傷，其中i表示當(dāng)前迭代步，l表示當(dāng)前荷載步）判斷該單元是否發(fā)生過(guò)損傷，即Dⁱ_l是否大于D⁰_l，其中，0是上一荷載步的最后迭代步，若“否”，則說(shuō)明該單元此前未發(fā)生損傷，進(jìn)入③；若“是”，則說(shuō)明該單元此前發(fā)生過(guò)損傷，進(jìn)入④。

③由②所求應(yīng)變s檢驗(yàn)是否發(fā)生拉損傷，如發(fā)生了則進(jìn)入⑤，若未發(fā)生拉損傷，則檢驗(yàn)是否為壓損傷，如發(fā)生了則進(jìn)入⑥，未發(fā)生ANSYS數(shù)據(jù)庫(kù)中，計(jì)算整體結(jié)構(gòu)響應(yīng)。

④判斷發(fā)生損傷的類(lèi)型：若為拉損傷，計(jì)算壓損傷變量Dⁱ_lc。

⑤計(jì)算相應(yīng)的拉損傷變量Dⁱ_lt以及弱化后彈性模量E。

⑥計(jì)算相應(yīng)的壓損傷變量Dⁱ_lc以及弱化后彈性模量E。

⑦計(jì)算單元矩陣，包括單元?jiǎng)偠汝?、阻尼陣、質(zhì)量陣以及單元等效荷載向量，然后將上面得到數(shù)據(jù)返回到ANSYS數(shù)據(jù)庫(kù)中，計(jì)算整體結(jié)構(gòu)響應(yīng)。如果計(jì)算結(jié)果不滿(mǎn)足收斂條件，則進(jìn)入下一平衡迭代步，若滿(mǎn)足則進(jìn)入下一荷載步，直至地震分析結(jié)束。

3 巖質(zhì)邊坡?lián)p傷演化模擬

3.1 邊坡計(jì)算參數(shù)及有限元模型

計(jì)算參數(shù)選?。簭椥阅Ａ縀=3000GPa，泊松比v=0.23，質(zhì)量密度ρ=2300kg/m³，動(dòng)抗拉強(qiáng)度f(wàn)_t=1.5×10⁵。選取材料的隨機(jī)力學(xué)參數(shù)為：彈性模量和抗拉強(qiáng)度的均質(zhì)度m=2，泊松比的均質(zhì)度m=100。選取瑞利阻尼進(jìn)行計(jì)算，計(jì)算時(shí)取第1階與第10階的固有頻率，陣型阻尼比為0.05。本文所建立的邊坡模型坡度為45°，坡高為20m，坡前計(jì)算長(zhǎng)度取20m，坡后計(jì)算長(zhǎng)度取44m，坡腳以下取20m。為研究隧道洞口段巖質(zhì)邊坡在地震作用下的損傷破壞過(guò)程，在平面應(yīng)變狀態(tài)下，采用等效損傷單元（幾何形狀如四節(jié)點(diǎn)平面應(yīng)變單元）對(duì)邊坡模型進(jìn)行有限元網(wǎng)格離散，離散后的單元尺寸約為0.4m×0.4m，模型共分18589個(gè)單元、18941個(gè)節(jié)點(diǎn)。模型左右兩邊施加水平約束，底部施加水平和豎向約束，邊坡有限元模型見(jiàn)圖4。分水平與豎直2個(gè)方向輸入地震波，其中水平峰值加速度為0.474g，豎向加速度峰值為0.312g，地震動(dòng)總的持續(xù)時(shí)間為10.56s，加速度時(shí)程曲線(xiàn)見(jiàn)圖5。

3.2 邊坡地震損傷演化分析

本文從地震荷載作用整個(gè)過(guò)程中，選取6個(gè)典型的時(shí)間點(diǎn)（2.64s、2.92s、5.78s、6.96s、8.98s、9.42s）來(lái)分析隧道洞口段巖質(zhì)邊坡從開(kāi)始出現(xiàn)微裂紋到形成清晰可見(jiàn)的裂縫，裂縫不斷擴(kuò)展至形成最終損傷區(qū)域總過(guò)程，其損傷演化過(guò)程及其主拉應(yīng)力的分布見(jiàn)圖6。隧道洞口段巖質(zhì)邊坡?lián)p傷演化分析過(guò)程分析中：在靜力荷載（本文僅施加了重力荷載）作用的基礎(chǔ)上施加地震荷載。地震荷載激勵(lì)的初期階段，邊坡絕大部分區(qū)域處于受壓狀態(tài)下，雖然結(jié)構(gòu)沒(méi)有出現(xiàn)嚴(yán)重的裂縫，但由于單元材料隨機(jī)分布的影響，少數(shù)材料強(qiáng)度較弱的單元依舊出現(xiàn)了不同程度的損傷，坡面開(kāi)始出現(xiàn)微裂紋，但邊坡結(jié)構(gòu)整體上處于完好的狀態(tài)。隨著地震激勵(lì)時(shí)間的增加，已出現(xiàn)的微裂紋不斷擴(kuò)展，新的微裂紋也不斷出現(xiàn)，這些微裂紋逐漸擴(kuò)展、延伸，最后在坡后形成數(shù)條清晰可見(jiàn)的裂縫，如圖6a、b所示；隨著地震動(dòng)的繼續(xù)震蕩，坡后裂縫范圍繼續(xù)擴(kuò)展、延伸，坡面發(fā)生破壞區(qū)域逐漸擴(kuò)大；洞口拱頂上方和坡面中間開(kāi)始出現(xiàn)裂紋并擴(kuò)展成較大裂縫，形成多處貫通的損傷區(qū)域，同時(shí)造成了主拉應(yīng)力的重新分布，如圖6c、d所示；最后，隨著地震動(dòng)的繼續(xù)震蕩，拱頂上方裂縫繼續(xù)擴(kuò)展，最終形成了貫穿整個(gè)拱頂上方的損傷區(qū)域，此時(shí)可以認(rèn)為拱頂上方已完全破壞，會(huì)給隧道的安全運(yùn)行帶來(lái)影響；坡面中間形成的自上而下裂縫會(huì)繼續(xù)擴(kuò)展、延伸，最終形成了一條清晰可見(jiàn)的裂縫，可以近似將該裂縫看成邊坡在強(qiáng)震作用下產(chǎn)生的滑坡趨勢(shì)，如圖6e、f所示。

4 結(jié)論

本文通過(guò)引入服從Weibull函數(shù)的細(xì)觀(guān)損傷微元體，運(yùn)用連續(xù)損傷介質(zhì)理論，理論上推導(dǎo)出了可反映巖石非均質(zhì)性和宏觀(guān)破裂與失穩(wěn)過(guò)程的巖石細(xì)觀(guān)損傷本構(gòu)關(guān)系式。

地震作用下巖質(zhì)邊坡的破壞是材料細(xì)觀(guān)損傷不斷產(chǎn)生、累積的結(jié)果。在地震荷載激勵(lì)的初期階段，邊坡絕大部分區(qū)域處于受壓狀態(tài)下，雖然結(jié)構(gòu)沒(méi)有出現(xiàn)嚴(yán)重的裂縫，但由于單元材料隨機(jī)分布的影響，少數(shù)材料強(qiáng)度較弱的單元依舊出現(xiàn)了不同程度的損傷，坡后開(kāi)始出現(xiàn)微裂紋，這些微裂紋逐漸貫通并形成數(shù)條清晰可見(jiàn)的裂縫。隨著地震動(dòng)的繼續(xù)震蕩，在洞口拱頂上方和坡面中間開(kāi)始出現(xiàn)新的裂縫。拱頂上方裂縫最終形成一個(gè)連續(xù)貫通的損傷區(qū)域，此時(shí)可以認(rèn)為拱頂上方已完全破壞，這會(huì)給隧道的安全運(yùn)行帶來(lái)影響，因此在洞口設(shè)計(jì)、施工時(shí)要對(duì)拱頂上方進(jìn)行必要的加固處理；坡面中間形成的自上而下的裂縫，可以看成是邊坡結(jié)構(gòu)在強(qiáng)震作用下產(chǎn)生的滑坡趨勢(shì)，設(shè)計(jì)、施工也要給與重視。