基于變形體離散元方法的裂縫模型及變形分析

基于變形體離散元方法的裂縫模型及變形分析

混凝土結(jié)構(gòu)的負(fù)荷影響下的損害和破壞，巖質(zhì)邊坡的變形和穩(wěn)定，以及不連續(xù)的地下孔室圍巖的穩(wěn)定性是當(dāng)前防災(zāi)和民防工程中必須解決的問題。這些問題存在一個(gè)顯著的特點(diǎn),系統(tǒng)在外載作用下從連續(xù)或斷續(xù)狀態(tài)逐漸轉(zhuǎn)變?yōu)榉沁B續(xù)狀態(tài)而形成獨(dú)立塊體,系統(tǒng)的力學(xué)行為最終由獨(dú)立塊體間的相互作用決定。采用數(shù)值方法描述系統(tǒng)的力學(xué)行為,實(shí)現(xiàn)連續(xù)介質(zhì)的損傷斷裂并向非連續(xù)介質(zhì)的轉(zhuǎn)化以致失穩(wěn)是當(dāng)今計(jì)算固體力學(xué)迫待研究的課題。Rashid于1968年提出彌散裂縫模型模擬材料開裂行為,Bazant發(fā)展的斷裂帶理論實(shí)現(xiàn)了斷裂過程的數(shù)值模擬,Rots采用固定裂縫模型和Jirasek采用旋轉(zhuǎn)裂縫模型研究了混凝土的Ⅰ型開裂,但以有限元為基礎(chǔ)的連續(xù)介質(zhì)力學(xué)方法限制于小位移問題,不能處理大變形失穩(wěn)問題,基于非連續(xù)介質(zhì)的離散元法在模擬連續(xù)介質(zhì)損傷斷裂力學(xué)方面研究較少。本文將準(zhǔn)脆性材料的開裂模型嵌入到離散單元法中,提出一種將損傷斷裂模型與變形離散元結(jié)合的方法以研究連續(xù)-非連續(xù)系統(tǒng)的過渡過程,文中首先討論材料開裂模型,不連續(xù)界面的引入以及裂縫表面的接觸模擬,離散塊體的運(yùn)動(dòng)求解等問題,文末以三點(diǎn)彎梁斷裂問題作數(shù)值解與試驗(yàn)比較。1三維非線性衰減理論模型依據(jù)裂縫方向的不同處理方法,基于斷裂帶理論的本構(gòu)關(guān)系,將斷裂力學(xué)模型分為固定裂縫和旋轉(zhuǎn)裂縫兩種。1.1土體裂裂變形模型固定裂縫模型模擬材料的一應(yīng)力點(diǎn)的開裂方向由該點(diǎn)在加載過程中達(dá)到抗拉強(qiáng)度時(shí)主應(yīng)力的方向所決定,在后續(xù)計(jì)算中不再發(fā)生改變?？紤]到裂縫面上由于骨料互瑣而產(chǎn)生的殘余剪切剛度,通常引入一隨著裂縫的張開而逐漸減小或常剪切抗力系數(shù),采用應(yīng)變分解方法,認(rèn)為混凝土總應(yīng)變?chǔ)庞苫炷敛牧献陨響?yīng)變?chǔ)與o和裂縫張開引起的應(yīng)變?chǔ)與r組成,即ε=εco+εcr.(1)以裂縫面法向?yàn)閤軸建立局部坐標(biāo)系,根據(jù)坐標(biāo)轉(zhuǎn)換關(guān)系有:εcr=Νecr?(2)s=ΝΤσ.(3)式中:N為坐標(biāo)轉(zhuǎn)換矩陣,開裂應(yīng)變ecr包含垂直裂縫面方向的應(yīng)變分量及與裂縫面平行的剪應(yīng)變分量,s為同一局部坐標(biāo)系下的應(yīng)力。假設(shè)混凝土材料自身滿足線彈性本構(gòu)關(guān)系,在局部坐標(biāo)系下裂縫應(yīng)變與應(yīng)力的關(guān)系為s=Dcrecr?(4)式中Dcr為綜合考慮裂縫拉伸和剪切軟化后構(gòu)成的3×3矩陣。綜合式(1)—(4),可得開裂混凝土材料的應(yīng)力-應(yīng)變?nèi)筷P(guān)系為σ=[Dco-DcoΝ(Dcr+ΝΤDcoΝ)-1ΝΤDco]ε?(5)式中Dco為混凝土彈性剛度矩陣。固定裂縫模型裂縫方向由達(dá)到抗拉強(qiáng)度時(shí)的主應(yīng)力方向決定,在后續(xù)的計(jì)算中不發(fā)生改變,隨著裂縫的擴(kuò)展,一點(diǎn)的材料主軸可能發(fā)生偏轉(zhuǎn),因此采用固定裂縫模型計(jì)算一點(diǎn)的力學(xué)行為通常與材料發(fā)生開裂后的實(shí)際力學(xué)特征不符。其后,學(xué)者們發(fā)展了固定裂縫模型,定義了一個(gè)開裂后連續(xù)旋轉(zhuǎn)的應(yīng)力主軸,在旋轉(zhuǎn)的主軸空間中建立應(yīng)力應(yīng)變關(guān)系,提出了旋轉(zhuǎn)裂縫模型,并應(yīng)用該模型在模擬混凝土開裂方面做了大量工作。1.2局部應(yīng)力應(yīng)變動(dòng)力系統(tǒng)的求解旋轉(zhuǎn)裂縫模型考慮當(dāng)一點(diǎn)主應(yīng)力分量達(dá)到抗拉強(qiáng)度時(shí),材料發(fā)生開裂,在開裂方向的應(yīng)力應(yīng)變分量滿足材料的軟化力學(xué)行為,發(fā)生開裂后的不同軟化階段,裂縫方向始終與當(dāng)前主應(yīng)力方向垂直,只需定義主軸空間中的應(yīng)力應(yīng)變關(guān)系,即:εcr=Τecr?(6)ecr=gs.(7)式中:ecr為局部坐標(biāo)系下3個(gè)主方向的開裂應(yīng)變,s為局部坐標(biāo)系下的應(yīng)力,g為裂縫特性矩陣,T為坐標(biāo)轉(zhuǎn)換矩陣。根據(jù)旋轉(zhuǎn)裂縫模型的假設(shè),在求解過程中需要不斷更新T矩陣,使得局部坐標(biāo)系隨主應(yīng)變方向變化而旋轉(zhuǎn)。以各點(diǎn)有一個(gè)主方向發(fā)生開裂為例:ecr1=g1s1?(8)g1=(1-μ)μE?(0＜μ≤1).(9)式中μ為材料割線模量相對(duì)于線彈性模量E的折減系數(shù)。綜合式(1)、(6)—(9),假定混凝土材料自身滿足線彈性本構(gòu)關(guān)系,在由主應(yīng)力方向建立的局部坐標(biāo)系中應(yīng)力應(yīng)變?nèi)筷P(guān)系為:{s1s2s3}=E1-ν-2μν2[μ(1-ν)μνμνμν1-μν21+νν(1+μν)1+νμνν(1+μν)1+ν1-μν21+ν]{e1e2e3}?(10)式中:si、ei(i=1,2,3)分別為3個(gè)主應(yīng)力和主應(yīng)變,E為彈性模量,ν為Poisson比。根據(jù)當(dāng)前的應(yīng)變狀態(tài)確定折減系數(shù)μ,以表征不同的開裂程度。1.3斷裂帶及材料參數(shù)對(duì)力學(xué)性能的影響在復(fù)雜應(yīng)力路徑作用下,準(zhǔn)脆性材料的開裂行為趨向于以Ⅰ型裂紋為主的方向。因此,混凝土的張拉軟化曲線成為上述裂縫模型的基礎(chǔ),依據(jù)斷裂帶理論,遵循斷裂能Gf唯一準(zhǔn)則,根據(jù)材料抗拉強(qiáng)度ft和基于網(wǎng)格尺寸確定的斷裂帶寬度h,調(diào)整材料應(yīng)變軟化曲線參數(shù)以準(zhǔn)確描述材料的斷裂能釋放,通常采用線性和雙線性軟化曲線描述準(zhǔn)脆性材料的力學(xué)響應(yīng)。兩種軟化曲線的臨界開裂應(yīng)變?yōu)?εf1=2Gf/fth?(11)εf2=3.6Gf/fth.(12)對(duì)于雙線性軟化曲線,軟化段中折點(diǎn)位置由式(13)和式(14)確定:εm=2εf2/9?(13)fm=ft/3.(14)2局部阻尼分析方法文提出的離散單元法將為不連續(xù)面所切割的系統(tǒng)差分為四面體常應(yīng)變單元以模擬塊體的變形,不連續(xù)界面通過接縫本構(gòu)模型規(guī)定界面上的力和位移關(guān)系,作用力被集中在差分節(jié)點(diǎn)上。該模型基于Newton第二定律,采用動(dòng)態(tài)松弛法求解節(jié)點(diǎn)運(yùn)動(dòng)方程,使用顯式中心差分方法求解節(jié)點(diǎn)的位移和單元應(yīng)變,各節(jié)點(diǎn)的動(dòng)力平衡方程為m??ui=∫sσijnjds+Fexti?(15)式中:m為分配在節(jié)點(diǎn)上的質(zhì)量,Fexti為節(jié)點(diǎn)所受的外力,s為包圍積分質(zhì)點(diǎn)域的外邊界,nj為邊界的單位外法向分量。由于采用動(dòng)力松弛解法,對(duì)于穩(wěn)態(tài)問題,通常采用粘性阻尼以消散系統(tǒng)的能量。對(duì)于有荷載突變或系統(tǒng)逐漸破壞的情況,采用局部阻尼,阻尼力與節(jié)點(diǎn)的不平衡力成正比,即Fdi=α|∑F(t)i|sgn(˙u(t-Δt/2)i)?(16)式中:Fdi為節(jié)點(diǎn)阻尼力,∑F(t)i為除阻尼力外節(jié)點(diǎn)所受的合力,α為小于1的常數(shù)。對(duì)運(yùn)動(dòng)方程采用差分格式得:˙u(t+Δt/2)i=˙u(t-Δt/2)i+(∑F(t)i-Fdi)Δtm?(17)u(t+Δt)i=uti+˙u(t+Δt/2)iΔt.(18)對(duì)于已有的不連續(xù)面以及材料開裂后形成的新的不連續(xù)面,通過法向和切向的彈簧傳遞接觸面之間的相互作用力,接觸模型法向采用拉斷準(zhǔn)則,切向遵循Mohr-Coulomb準(zhǔn)則:ΔFn=-ΚnΔunAc?(19a)ΔFs=-ΚsΔusAc?(19b)Fn=0?Fs=0?當(dāng)Fn＞ftAc?(20a)Fs=Fs(f|Fn|+CAc)/|Fs|?當(dāng)|Fs|＞f|Fn|+CAc?(20b)式中:Kn、Ks是接觸法向和切向彈簧的剛度,Δun、Δus為相對(duì)位移增量的法向和切向分量,ΔFn、ΔFs為法向和切向彈簧力變化量,Fn、Fs為切向和法向彈簧力的大小,Ac為子接觸面積,f為接觸面的摩擦因數(shù),C為粘聚力?？勺冃坞x散元能夠很好模擬不連續(xù)介質(zhì)的靜動(dòng)力力學(xué)行為,由于其采用增量解法,可以較容易將彌散裂縫模型引入到材料本構(gòu)模型中以模擬材料開裂行為。3材料的拉拔試驗(yàn)對(duì)于準(zhǔn)脆性材料,定義開裂系數(shù)Ff即沿裂縫法向的非彈性開裂應(yīng)變?chǔ)與rn和最大開裂應(yīng)變?chǔ)舊的比來表征開裂程度,即Ff=εcrn/εf.(21)當(dāng)Ff為0時(shí),材料處于線彈性階段,當(dāng)0<Ff<1時(shí),材料處于軟化階段,當(dāng)Ff=1時(shí),沿裂縫面的抗拉強(qiáng)度降為零,此時(shí)通過常應(yīng)變四面體單元形心插入一離散裂縫,在該離散裂縫面上形成新的接觸,接觸面的抗拉強(qiáng)度和粘聚力均為0,f采用本體材料的摩擦因數(shù),接觸的法向和切向仍遵循Mohr-Coulomb準(zhǔn)則。對(duì)于固定裂縫模型,離散裂縫方向垂直于裂紋起裂時(shí)的最大主應(yīng)力方向,對(duì)于旋轉(zhuǎn)裂縫模型,離散裂縫方向垂直于Ff=1時(shí)的最大主應(yīng)力方向,而忽略開裂過程中裂縫方向的旋轉(zhuǎn)。為此,根據(jù)開裂系數(shù)Ff的大小可以確定彌散裂縫張開的程度以及離散裂縫插入的順序,最終形成非連續(xù)或斷續(xù)系統(tǒng)。4網(wǎng)格劃分的確定選取Peterson三點(diǎn)彎梁試驗(yàn),梁的尺寸為2000mm×200mm×50mm,切口深度為半梁高,力學(xué)計(jì)算參數(shù)為:混凝土彈性模量E=30GPa,Poisson比ν=0.2,抗拉強(qiáng)度ft=3.33MPa,斷裂能Gf=124N/m。計(jì)算采用兩套網(wǎng)格,其主要差別在于切口寬度和切口處的網(wǎng)格劃分,細(xì)網(wǎng)格的單元尺寸為粗網(wǎng)格(切口寬度20mm)的單元尺寸的一半,從正面看粗網(wǎng)格劃分情況如圖1所示。本文采用三維固定和旋轉(zhuǎn)裂縫兩種模型,分別應(yīng)用線性、雙線性軟化曲線計(jì)算施加的總荷載F與切口頂部節(jié)點(diǎn)的位移D關(guān)系曲線如圖2—圖6所示。5結(jié)構(gòu)承載能力分析依據(jù)圖2—圖5,按照Bazant斷裂帶理論,每個(gè)模型中不同粗細(xì)的網(wǎng)格剖分計(jì)算得到的荷載位移關(guān)系曲線基本相同,說明在非線性斷裂力學(xué)模型中,當(dāng)計(jì)算遵循斷裂能守恒準(zhǔn)則,根據(jù)數(shù)值模擬的斷裂帶寬度的大小來調(diào)整應(yīng)變軟化參數(shù),可以克服網(wǎng)格剖分尺寸對(duì)計(jì)算結(jié)果的影響,同時(shí)也說明了細(xì)網(wǎng)格并不能克服由于網(wǎng)格線與開裂方向不一致而帶來的剪切鎖死現(xiàn)象,即目前的局部化模型存在一定程度的網(wǎng)格方向依賴性。對(duì)比試驗(yàn)結(jié)果的上下包絡(luò)線與數(shù)值模擬結(jié)果,可知采用線性軟化曲線較雙線性軟化曲線計(jì)算的峰值荷載偏大,荷載位移曲線下降較陡,雙線性軟化模型后期承載力偏大,這是由于對(duì)于同一斷裂帶寬度,在保證斷裂能守恒的情況下,雙線性軟化曲線在拐點(diǎn)以前的下降段較線性軟化曲線陡,在拐點(diǎn)后下降段相比較緩所致。從圖6中可以看出,在裂縫發(fā)生的初始階段,固定裂縫模型和旋轉(zhuǎn)裂縫模型的計(jì)算結(jié)果基本相同,而后旋轉(zhuǎn)裂縫模型計(jì)算的殘余承載能力偏小,更接近于試驗(yàn)值。分析原因知隨著裂縫擴(kuò)展,裂縫應(yīng)變?cè)诳倯?yīng)變分布中占主導(dǎo)地位,由于開裂方向與網(wǎng)格線的不平行造成彌散表征裂縫的單元內(nèi)產(chǎn)生了虛假的剪切應(yīng)力,裂縫單元主應(yīng)力方向發(fā)生一定程度的偏轉(zhuǎn)。采用固定裂縫模型使得該虛假剪應(yīng)力的作用隨著裂縫的進(jìn)一步擴(kuò)展而增強(qiáng),系統(tǒng)剛度較實(shí)際情況偏大;采用旋轉(zhuǎn)裂縫模型由于裂縫法向始終與主應(yīng)力方向共軸,一定程度上避免了主應(yīng)力方向偏離初始開裂方向而帶來的剪切鎖死現(xiàn)象,計(jì)算的結(jié)構(gòu)承載能力更符合實(shí)際情況。圖7給出了粗網(wǎng)格劃分情況下不同斷裂模型的裂縫開裂方向,旋轉(zhuǎn)裂縫模型與固定裂縫模型相比開裂方向與梁所承受的對(duì)稱荷載更為相符。6支護(hù)結(jié)構(gòu)的變形1)本文采用非線性斷裂力學(xué)模型,考慮混凝土、巖石等準(zhǔn)脆性材料的受拉開裂,基于變形體離散元方法分析了材料受拉破壞行為。模型能夠定量預(yù)測(cè)初始裂紋時(shí)的應(yīng)力狀態(tài)及開裂過程中的變形,獲得了合理的結(jié)構(gòu)荷載位移響應(yīng)曲線,為連續(xù)-非