改進(jìn)鍵型近場動力學(xué)方法下的多裂紋板破壞分析

孫 杰, 徐業(yè)鵬

(河海大學(xué)力學(xué)與材料學(xué)院，南京 211100)

裂縫的存在通常會削弱工程結(jié)構(gòu)的承載能力，甚至導(dǎo)致結(jié)構(gòu)失效，因此含缺陷材料和結(jié)構(gòu)的斷裂問題及其數(shù)值模擬一直都是計(jì)算力學(xué)的研究重點(diǎn)。近年來，近場動力學(xué)方法以其在模擬裂紋擴(kuò)展以及分析結(jié)構(gòu)斷裂方面的獨(dú)特優(yōu)勢日益受到國內(nèi)外學(xué)者的青睞。

近場動力學(xué)(peridynamics，PD)是基于非局部作用思想建立模型并通過求解空間積分方程來描述物質(zhì)力學(xué)行為的一種方法。該方法由Silling[1]于2000年提出，隨后，為證明求解不存在裂尖應(yīng)力奇異性，Silling等[2]又對PD基于鍵的理論進(jìn)行更深一步的推導(dǎo)，成功應(yīng)用于研究微彈性材料[3-4]和混凝土材料[5]中的復(fù)雜裂紋問題。常規(guī)微彈性脆性模型(prototype microelastic brittle,PMB)在應(yīng)用時存在限制，為此Gerstle等[6]在開展多晶斷裂和PD斷裂模擬的同時，引入了彎矩密度概念來解除泊松比的限制，引入了雙參數(shù)微極模型。Huang等[7]、顧鑫等[8]為提高計(jì)算精度，又引入核函數(shù)修正項(xiàng)。黃丹等[9-10]進(jìn)一步分析改進(jìn)了核函數(shù)修正項(xiàng)對PD的作用，并開發(fā)出相應(yīng)的靜力動力的算法。錢劍等[11]在此基礎(chǔ)上進(jìn)行了動載作用下的復(fù)合型裂紋擴(kuò)展模擬。

與此同時，多裂紋擴(kuò)展問題中裂紋間的相互影響和連接一直都是較為復(fù)雜的難題，眾多學(xué)者針對多裂紋問題進(jìn)行過大量研究。Silling[12]用沖擊物撞擊了帶有平行裂紋的圓盤，利用近場動力學(xué)成功模擬此算例，在試驗(yàn)中觀察到了裂紋擴(kuò)展的角度。Zhou等[13]研究了裂紋陣列對多裂紋擴(kuò)展和聚結(jié)過程的影響。Zeng等[14]用擴(kuò)展有限元方法，觀察脆性巖石類材料的裂紋對其彈性性能和強(qiáng)度的影響情況。Vazic等[15]研究了小裂紋對宏觀裂紋動力學(xué)擴(kuò)展的影響，考慮了不同數(shù)量、不同位置和密度的小裂紋的各種組合，其大小取決于小裂紋的位置、密度和數(shù)量。張振南等[16]基于VMIB方法，通過對兩條不同排列的平行裂紋開裂過程的模擬，探究多裂紋之間的相互影響?；诟倪M(jìn)后的近場動力學(xué)方法，對含雙裂紋的脆性板進(jìn)行單軸拉伸破壞模擬，并將所得結(jié)果與已有文獻(xiàn)結(jié)果進(jìn)行對比，驗(yàn)證本文提出的本構(gòu)模型和算法的可靠性。隨后通過模擬不同初始形態(tài)的多裂紋板的裂紋擴(kuò)展路徑和分析承載能力，探究裂紋初始狀態(tài)對結(jié)構(gòu)破壞型式、擴(kuò)展路徑和承載能力的影響規(guī)律。

1 近場動力學(xué)基本理論

1.1 鍵型PD基本思想

近場動力學(xué)的基本思想如下：在假設(shè)物質(zhì)點(diǎn)系統(tǒng)做剛體運(yùn)動，保持物質(zhì)系統(tǒng)構(gòu)型不變并且物質(zhì)點(diǎn)對之間的作用與時間無關(guān)的前提下，考慮某一個空間R內(nèi)的所有物質(zhì)點(diǎn)，領(lǐng)域內(nèi)任意物質(zhì)點(diǎn)x和一定范圍大小δ(近場范圍)內(nèi)的其他物質(zhì)點(diǎn)x′∈H:{|x′-x|≤δ,H∈R}存在關(guān)于位移u=u(x,t)的相互作用力f(圖1)，那么對于該物質(zhì)點(diǎn)就有如下關(guān)系：

f=f[x,x′,u(x,t),u(x′,t),t]

(1)

(2)

式中：H為近場范圍；ρ為材料密度；b為外力密度。

圖1 物質(zhì)點(diǎn)對之間的相互作用

對于均質(zhì)材料，點(diǎn)對力函數(shù)可以繼續(xù)簡化為與參考構(gòu)型中兩物質(zhì)點(diǎn)的相對位置ξ=x′-x和相對位移η=u′-u相關(guān)的函數(shù)：

f=f(x′-x,u′-u)=f(ξ,η)

(3)

均勻材料的本構(gòu)力函數(shù)的一般形式為

f(η,ξ)=F(ξ)(η+ξ)

(4)

根據(jù)Silling[1]以及黃丹等[17],本構(gòu)力(點(diǎn)對力)可表示為

(5)

式(5)中：μ判斷物質(zhì)點(diǎn)對破壞情況。

(6)

式(6)中：s0為物質(zhì)點(diǎn)對的臨界伸長率，可取s0=[4πG0/(9Eδ)]1/2。

物質(zhì)點(diǎn)對勢能與傳統(tǒng)應(yīng)變能等效[1]求解微模量系數(shù)c，結(jié)果如下：

平面應(yīng)力問題：

(7)

平面應(yīng)變問題：

(8)

式中：E為彈性模量；ν為泊松比。

常規(guī)PMB模型在平面應(yīng)力問題求解時受到泊松比限制，為突破該限制，Gerstle等[6]引入物質(zhì)點(diǎn)對的轉(zhuǎn)動自由度，提出雙參數(shù)微極模型：

f(η,ξ)=D(ξ)η

(9)

令ξ=|ξ|，其中，

(10)

通過物質(zhì)點(diǎn)對勢能與傳統(tǒng)應(yīng)變能等效[10]，求出雙參數(shù)微極模型的參數(shù)：

(11)

為了提高計(jì)算精度，加入考慮了長程力的尺寸效應(yīng)，加入核函數(shù)修正項(xiàng)g=[1-(ξ/δ)2]2。將原來的常規(guī)微彈性脆性模型系數(shù)修正為c(ξ)=c(0,δ)g(ξ)。

核函數(shù)修正項(xiàng)必須滿足長程力的變化規(guī)律：

(12)

引入核函數(shù)修正項(xiàng)g=[1-(ξ/δ)2]2后得到修正的微模量系數(shù):

(13)

1.2 求解體系

在任意時刻，x由近場范圍內(nèi)物質(zhì)點(diǎn)對作用力產(chǎn)生的內(nèi)部合力可以表示為

?x∈R,t≥0

(14)

PD是將固體離散為一系列具有信息的物質(zhì)點(diǎn)，設(shè)物質(zhì)點(diǎn)間距為|Δx|，則將式(2)改寫為

(15)

通常在求解時加入人工阻尼來保證獲得靜力解，基本方程如下：

(16)

采用中心差分格式

(17)

按時間順序進(jìn)行離散,代入(15)化簡，得到：

(18)

由Fortran程序完成求解。

2 算例分析

2.1 模型驗(yàn)證

為驗(yàn)證本文本構(gòu)模型和算法的正確性，首先對含雙裂紋混凝土板進(jìn)行了單軸拉伸破壞模擬。具體模型如圖2所示，彈性模量E=30 GPa，泊松比ν=0.33，材料密度為ρ=2 400 kg/m3，離散為10 200個物質(zhì)點(diǎn)，離散間距為0.000 5 m,近場尺寸δ=3Δx=0.001 5 m。在模擬過程中采用位移加載控制方式對其進(jìn)行加載，每一步荷載增量為0.04 mm，迭代步長為Δt=1×10-7s。

為模擬不同裂紋初始布置下板的斷裂破壞模式，采用如下三種方案(均保持兩條裂紋平行)。

(1)方案I：2b/2a=1，即b=5 mm，α=45°。

(2)方案II：2b/2a=2，即b=10 mm，α=45°。

(3)方案III：2b/2a=1，即b=5 mm，α=60°。

模擬結(jié)果如圖3(方案I～I(xiàn)II的裂紋擴(kuò)展過程)所示，裂紋沿著最易匯合的方向擴(kuò)展，一起沿水平方向開裂，在兩條裂紋之間形成破裂區(qū)，進(jìn)而貫通，同時兩條裂紋也向板的兩側(cè)擴(kuò)展，最后形成一條貫通性裂紋。三個方案中PD所得裂紋擴(kuò)展路徑與VMIB[16]方法所得裂紋擴(kuò)展路徑基本一致。

圖2 雙裂紋板模型

圖3 裂紋擴(kuò)展路徑比對

2.2 含三裂紋脆性板單軸拉伸破壞分析

主要研究在含三裂紋脆性板在受到單軸拉伸荷載作用時，裂紋的初始布置是否會改變裂紋擴(kuò)展的路徑以及對臨界破壞荷載是否有影響。本模型具體材料參數(shù)和2.1節(jié)一致，q=1 000 kPa，迭代時間步長Δt=1×10-7s。模型如圖4所示。

圖4 三裂紋板模型

觀察裂紋的擴(kuò)展路徑主要有兩種形式，即以中心裂紋為主的開裂與以兩邊裂紋為主的開裂，最后都是沿著豎直方向擴(kuò)展直到板完全破壞。以α=30°為例，當(dāng)β=0°時，脆性板以兩邊裂紋為主開裂，直至破壞，當(dāng)β逐漸增大后，裂紋擴(kuò)展路徑開始發(fā)生變化，由以周邊裂紋為主的開裂變?yōu)橐灾行牧鸭y為主的開裂，如圖5(a)所示。

觀察有無周邊裂紋時臨界破壞荷載，如“β=0°且無周邊裂紋”與“β=α=0°”對比可得，周邊裂紋的出現(xiàn)明顯削弱了板的承載能力，但同時對比“β=15°且無周邊裂紋”與“α=0°，β=15°”，板承載能力又有稍許提升，可見裂紋的增加對板的承載能力存在復(fù)雜的影響。進(jìn)一步觀察圖5(b)可知，臨界破壞荷載隨著中心裂紋角度與周邊裂紋角度的增加而減小。當(dāng)三條裂紋角度都在0°時，板的承載能力最強(qiáng)，即此時板處在最穩(wěn)定的狀態(tài)，而當(dāng)有一條裂紋達(dá)到90°時，臨界破壞荷載達(dá)到最小值。在周邊裂紋角度值一定的情況下，板的承載能力隨著中心裂紋角度的增加而減小。同樣，在中心裂紋角度一定的情況下。臨界破壞荷載值也隨著α的增大而減小。即板整體承載能力隨著中心裂紋角度與周邊裂紋角度的增加而減小。

圖5 PD模擬結(jié)果

2.3 含多裂紋脆性板單軸拉伸破壞分析

含初始裂紋板的破壞問題[18-20]在實(shí)際工程與航空航天領(lǐng)域都有重要的作用，在文獻(xiàn)[21]中，由于航空事故對航空問題中的鋁板進(jìn)行了裂紋破壞問題分析，而該類問題在脆性板中也具有深遠(yuǎn)研究意義。采用的模型如圖6所示，由于初始裂紋布置對板的承載能力有著不可忽略的影響，通過改變角度α與中心裂紋角度β，觀察初始裂紋布置對臨界破壞荷載的影響。材料數(shù)據(jù)如下：彈性模量E=0.69 GPa，泊松比ν=0.3，材料密度為ρ=2 400 kg/m3，取物質(zhì)點(diǎn)間距為0.001 m，物質(zhì)點(diǎn)總數(shù)為30 200，近場尺寸δ=3Δx=0.003 m。對其進(jìn)行均布拉伸，q=1 000 kPa，迭代步長為Δt=1×10-7s。具體模型如圖6所示。

圖6 多裂紋板模型

由圖7易知，在α與β同時達(dá)到90°時，臨界破壞荷載達(dá)到最大值，幾乎接近板的抗拉強(qiáng)度。而在α=β=0°時，臨界破壞荷載最小，即板處在最脆弱的狀態(tài)。

在角度β不變的情況下，臨界破壞荷載隨著周邊裂紋角度α的增加而增加，同樣，在α不變的情況下，臨界破壞荷載也幾乎隨著β的增加而增加，僅在β=45°時，α由15°～30°時臨界破壞荷載出現(xiàn)了較為明顯的下降，整體規(guī)律為隨著α與β的增大，臨界破壞荷載越大，即板的承載能力越強(qiáng)，越穩(wěn)定。反之當(dāng)α與β的越小，臨界破壞荷載越小，即板越易破壞。

此外，在不同階段，臨界破壞荷載隨裂紋角度變化幅度也有差異。在中心裂紋角度β<45°的時候，板的承載能力受周邊裂紋角度α影響較小，例如在β=15°時，圖7中折線較為平穩(wěn)，即臨界破壞荷載隨角度α變化波動不大，當(dāng)β>45°的情況，整體折線相對較陡，即臨界破壞荷載隨角度α變化波動較為明顯。

圖7 臨界破壞荷載示意圖

3 結(jié)論

(1)傳統(tǒng)鍵型PD理論在處理裂紋問題時具有獨(dú)特的優(yōu)勢，在常規(guī)微彈脆性近場動力學(xué)本構(gòu)模型中引入表征非局部長程作用力強(qiáng)度尺寸效應(yīng)的核函數(shù)修正項(xiàng)，構(gòu)建雙參數(shù)微彈脆性近場動力學(xué)本構(gòu)模型，對含雙裂紋脆性板進(jìn)行拉伸破壞模擬，將所得結(jié)果與VMIB方法所得結(jié)果進(jìn)行對比，證明該模型與方法的可行性。

(2)對于含中心裂紋的三裂紋板，裂紋的初始角度對裂紋擴(kuò)展路徑與臨界破壞荷載有一定影響。臨界破壞荷載隨著中心裂紋角度與周邊裂紋角度的增加而減小，當(dāng)三條裂紋角度都在0°時，板的承載能力最強(qiáng)，即此時板處在最穩(wěn)定的狀態(tài)，而當(dāng)有一條裂紋達(dá)到90°時，臨界破壞荷載達(dá)到最小。

(3)多裂紋板算例中，在α與β同時到達(dá)90°時，臨界破壞荷載達(dá)到最大，幾乎接近板的抗拉強(qiáng)度。而在α=β=0°時，臨界破壞荷載最小，即板處在最脆弱的狀態(tài)。在中心裂紋角度β<45°的時候，板的承載能力受周邊裂紋角度α影響較小，反之影響較大。