廣義微極近場動力學(xué)模型及三維斷裂行為模擬

陳希卓，禹海濤，朱建波，劉建鋒

（1.同濟大學(xué)土木工程學(xué)院，上海 200092；2.四川大學(xué)深地科學(xué)與工程教育部重點實驗室，四川成都 610065；3.同濟大學(xué)土木工程防災(zāi)國家重點實驗室，上海 200092；4.深圳大學(xué)土木與交通工程學(xué)院，廣東深圳 518061）

材料的裂紋擴展和破壞問題一直都是固體力學(xué)研究和工程領(lǐng)域面臨的關(guān)鍵難題[1］。目前學(xué)者們已經(jīng)提出了眾多的力學(xué)模型和數(shù)值方法來模擬固體材料斷裂過程，如擴展有限元[2-3］、內(nèi)聚力模型[4-5］和邊界元方法[6］等。然而，這些方法都需要引入特定的附加函數(shù)和斷裂準(zhǔn)則來描述不連續(xù)行為，無法有效模擬裂紋的自發(fā)萌生和擴展過程[7-8］。此外，由于傳統(tǒng)的數(shù)值方法大都是基于連續(xù)性假設(shè)，因此在處理斷裂等非連續(xù)行為時必須面對導(dǎo)數(shù)奇異性問題，從本質(zhì)上難以準(zhǔn)確模擬破壞問題。

為了解決這一難題，Silling[7］提出了近場動力學(xué)（PD）理論，采用考慮非局部作用的積分模型代替?zhèn)鹘y(tǒng)理論的微分模型，有效避免了傳統(tǒng)連續(xù)介質(zhì)力學(xué)方法在處理斷裂行為時的導(dǎo)數(shù)奇異性問題。該方法實現(xiàn)了由連續(xù)到不連續(xù)、微觀到宏觀力學(xué)作用的統(tǒng)一描述，尤其在分析裂紋擴展等非連續(xù)力學(xué)問題時具有顯著優(yōu)勢。近些年來，PD方法已逐漸成為計算力學(xué)及相關(guān)領(lǐng)域的關(guān)注熱點，其主要可以分為鍵基和態(tài)基2種模型[8］。對比態(tài)基近場動力學(xué)模型，鍵基模型的計算量更小，同時不存在零能模態(tài)的問題，因此廣泛應(yīng)用于巖石和混凝土等準(zhǔn)脆性材料的斷裂行為模擬。朱其志等[9］運用近場動力學(xué)方法對含有預(yù)制裂隙的巖石類材料試件的單軸壓縮試驗進(jìn)行數(shù)值模擬，分析了不同裂隙傾角對裂紋擴展模式的影響。Rabczuk和Ren[10］基于帶對偶-近場作用的近場動力學(xué)（DHPD）模型實現(xiàn)了巖石內(nèi)部裂紋擴展和準(zhǔn)脆性斷裂過程模擬。Wang等[11］發(fā)展了一種新型共軛鍵基模型，模擬了單軸壓縮條件下預(yù)先含有裂隙的巖石試樣中裂隙的分叉和合并行為，揭示了單軸壓縮條件下巖石的裂紋擴展機制。Gerstle等[12］提出了微極近場動力學(xué)（MPPD）模型，采用Euler-Bernoulli梁模型描述物質(zhì)點之間的相互作用，有效模擬了混凝土構(gòu)件失穩(wěn)和破壞過程。Diana等[13］運用微極近場動力學(xué)分析了砂巖在三點彎曲試驗中的力學(xué)特性，有效捕捉到了試件的I型和混合型斷裂過程。

上述近場動力學(xué)模型在準(zhǔn)脆性材料斷裂過程的模擬中均取得了不錯的效果，但其主要還是集中于二維斷裂問題。實際上，準(zhǔn)脆性材料在三維條件下會出現(xiàn)更復(fù)雜的體破壞現(xiàn)象，裂紋的空間擴展模式也會呈現(xiàn)較大差別[14］。而傳統(tǒng)的微極近場動力學(xué)模型對復(fù)雜條件下材料三維斷裂行為的模擬精度不高，主要原因在于其無法保證不均勻應(yīng)變場下近場動力學(xué)和傳統(tǒng)連續(xù)介質(zhì)力學(xué)的能量一致性以及缺乏相應(yīng)的三維斷裂準(zhǔn)則。

本文提出一種廣義微極近場動力學(xué)（GMPD）模型，在傳統(tǒng)微極近場動力學(xué)（MPPD）模型基礎(chǔ)上進(jìn)一步考慮三維狀態(tài)下鍵的軸向變形、切向變形、相對轉(zhuǎn)角之間的耦聯(lián)作用，可有效提高三維斷裂問題的模擬精度。通過引入對應(yīng)于鍵拉、剪、彎力學(xué)行為的3種近場動力學(xué)參數(shù)，以保證復(fù)雜荷載條件下近場動力學(xué)與傳統(tǒng)連續(xù)介質(zhì)力學(xué)應(yīng)變能的一致性。另外建立基于能量的新型斷裂準(zhǔn)則，旨在實現(xiàn)準(zhǔn)脆性材料的三維斷裂過程模擬與準(zhǔn)確預(yù)測。

1 廣義微極近場動力學(xué)模型

基于非局部作用的思想，近場動力學(xué)理論假定物體內(nèi)的任一物質(zhì)點位置矢量xi與其周圍一定區(qū)域Hx內(nèi)的其他任意物質(zhì)點位置矢量xj之間存在相互作用，這種作用可以理解為鍵，而力通過鍵在物質(zhì)點間進(jìn)行傳遞。為了消除泊松比的限制，微極近場動力學(xué)模型引入了Euler-Bernoulli梁模型來描述物質(zhì)點之間的相互作用，但是該模型忽略了鍵的剪切變形和轉(zhuǎn)動慣量，無法對復(fù)雜加載條件下的力學(xué)行為尤其是三維斷裂問題進(jìn)行準(zhǔn)確模擬。鑒于此，提出了一種廣義微極近場動力學(xué)模型，引入了Timoshenko梁來模擬物質(zhì)點間的相互作用，進(jìn)一步提高模型的模擬精度與適用性。

1.1 控制方程

基于Timoshenko梁理論，xi在t時刻的運動控制方程可表示為

式中：ρ為物質(zhì)密度；u（xi，t）為xi的位移向量；Hx為xi的作用域；f為xi和xj間的相互作用力；θi和θj分別為xi和xj的相對轉(zhuǎn)角；b（xi，t）為作用于xi上的外部體力；I為轉(zhuǎn)動慣量；A為鍵的截面面積；m為xi和xj間的彎矩相互作用；n（xi，t）為作用于xi上的外部彎矩；η和ξ分別為物質(zhì)點間的相對位置向量和相對位移向量，可表示為

基于局部坐標(biāo)系，鍵的對偶力函數(shù)f與對偶力矩m可分解為

式中：e1、e2和e3分別為局部坐標(biāo)系下沿x、y和z軸的單位向量。

引入分別對應(yīng)拉伸、剪切和彎曲剛度的鍵參數(shù)CN、Cθ和CM，并結(jié)合Timoshenko梁單元受力方程，可得到單根鍵上的力和彎矩在局部三維坐標(biāo)系下的表達(dá)為

1.2 應(yīng)變能與鍵參數(shù)

通過對鍵的微勢能w進(jìn)行域內(nèi)的積分，可以得到該模型的應(yīng)變能密度WPD。

式中：d為局部坐標(biāo)系下鍵的變形分量的統(tǒng)一表達(dá)；K為鍵的剛度矩陣。

建立如圖1所示的球面坐標(biāo)系OX′Y′Z′，三維條件下GMPD模型的應(yīng)變能密度可以進(jìn)一步表示為

圖1 球坐標(biāo)系下鍵的變形構(gòu)型Fig.1 A bond with the deformed configuration in global spherical coordinate system

式中：δ為作用域的半徑大??；α為鍵與球坐標(biāo)系中z′軸的夾角；β為鍵在Ox′y′平面的投影與x′軸的夾角。

基于近場動力學(xué)應(yīng)變能WPD和傳統(tǒng)連續(xù)介質(zhì)力學(xué)應(yīng)變能WCM的等價關(guān)系可以得到三維條件下鍵參數(shù)的表達(dá)為

將式（5）和式（6）代入式（8），并通過局部到整體的坐標(biāo)變換可以得到GMPD模型的應(yīng)變能密度表達(dá)，如式（9）：

而傳統(tǒng)連續(xù)介質(zhì)力學(xué)的應(yīng)變能密度WCM為

1.3 斷裂準(zhǔn)則

為了有效捕捉準(zhǔn)脆性材料的破壞過程，提出了一種基于能量的新型鍵基斷裂準(zhǔn)則。在該模型中，假定當(dāng)穿越某一平面的鍵AB全部斷開時，物質(zhì)點A的域內(nèi)形成相應(yīng)的宏觀裂紋，如圖2所示，因此材料的平面斷裂能應(yīng)等于所有斷開鍵的應(yīng)變能總和，如式（13）：

圖2 斷裂演化過程示意Fig.2 Evaluation of fracture energy

展開求解式（13），可以得到鍵的軸向變形、剪切變形以及相對轉(zhuǎn)角的限值s0、γ0和θ0，分別為

當(dāng)鍵的變形達(dá)到相應(yīng)的峰值后會自動斷開，即鍵兩端的物質(zhì)點不再產(chǎn)生相互作用，從而在宏觀上形成裂紋。這里引入標(biāo)量函數(shù)μ（xi，xj，t）描述物質(zhì)點間鍵是否發(fā)生破壞，如式（15）：

基于上述標(biāo)量函數(shù)，xi的損傷度d（xi，t）可以表示為

損傷度為0，代表物質(zhì)點完全無損傷；損傷度為1，表明物質(zhì)點周圍的鍵完全斷裂。

2 數(shù)值實現(xiàn)

為了求解GMPD模型中的控制方程，首先在空間上對求解域進(jìn)行離散化處理，在這里將求解域離散成物質(zhì)點，每個物質(zhì)點都占據(jù)一定空間體積且擁有一定物理性質(zhì)，離散后的運動控制方程可以表示為

式中：ΔVj為空間離散化處理后xi域內(nèi)的xj所占據(jù)的空間體積。

采用顯示中心差分方法求解上述控制方程，用有限差分近似代替位移對時間的導(dǎo)數(shù)，可以得到物質(zhì)點的加速度和速度的表達(dá)為

將式（19）代入式（17），可以得到位移的遞推求解格式為

將式（20）代入式（18），可以得到轉(zhuǎn)角的遞推求解格式為

該模型的邊界條件主要包括位移邊界和力邊界。其中外力P（x）是通過轉(zhuǎn)化為體力b（x）作用于最外層的物質(zhì)點Lr上，如圖3所示，其中Δ為物質(zhì)點間距。具體的轉(zhuǎn)換關(guān)系為

位移邊界主要是作用于物質(zhì)點外的虛擬邊界上，通過虛擬邊界層帶動內(nèi)部物質(zhì)點一起運動。如圖3所示，虛擬邊界的尺寸和作用域半徑δ一致。

圖3 力和位移邊界條件Fig.3 Force and displacement boundary conditions

GMPD模型的完整數(shù)值求解流程如圖4所示。

圖4 數(shù)值求解流程Fig.4 Computational flowchart of proposed model

3 三維斷裂模擬

為了說明GMPD模型對于三維斷裂問題的適用性，基于該模型模擬單軸壓縮條件下準(zhǔn)脆性材料試件的準(zhǔn)靜態(tài)破壞過程，并通過與傳統(tǒng)微極模型以及文獻(xiàn)中已有試驗結(jié)果進(jìn)行對比來驗證模型的有效性。

選取含單條初始裂隙的立方體試件為對象，如圖5所示?；谖墨I(xiàn)[14］的試驗參數(shù)，試件的彈性模量取E=3.81×109Pa，密度ρ=1.21g·cm-3，泊松比v=0.4，表面斷裂能G0=452N·m-1，初始裂紋半徑a=5mm，裂紋傾角α=30°，在兩側(cè)施加速率為2mm·min-1的軸向壓縮荷載。在近場動力學(xué)模擬中，假定試件為均勻和各向同性的，采用間距Δ=0.001m的物質(zhì)點對試件進(jìn)行離散，作用域的大小取δ=3Δ。

圖5 壓縮荷載下含初始裂隙試件[14］（單位：mm）Fig.5 Single flawed samples under uniaxial compression[14](unit:mm)

圖6為加載過程中本模型預(yù)測的試件裂紋擴展?fàn)顟B(tài)。由圖6a可知，預(yù)制裂紋邊緣最先開始發(fā)生破壞并逐漸形成翼型裂紋。在翼型裂紋不斷延伸的同時，反翼型裂紋（方向和翼型裂紋相反的拉伸裂紋）也開始在預(yù)制裂紋尖端萌生并擴展（見圖6b、圖6c）。隨著荷載的逐漸增加，翼型裂紋邊緣產(chǎn)生了次生裂紋，裂紋沿水平方向延伸至試件邊緣造成材料破壞（見圖6d）。本文模擬所得裂紋類型及其擴展路徑與試驗結(jié)果[14］基本一致（如圖7），有效捕捉到了三維空 間內(nèi)產(chǎn)生的翼型裂紋、反翼型裂紋以及次生裂紋。

圖6 加載過程中GMPD模型預(yù)測的裂紋擴展過程Fig.6 Predicted crack propagation in specimen under axial compression

圖7 壓縮荷載下試驗觀測的試件裂紋擴展過程[14］Fig.7 Experimental observation of crack development of specimen under axial compression[14]

為了進(jìn)一步說明二維斷裂與三維斷裂模擬的差異性，圖8給出了相同荷載作用下傳統(tǒng)MPPD模型預(yù)測的二維裂紋的擴展模式。

由圖8可見，二維裂紋發(fā)展主要以翼型裂紋為主，且擴展方向與荷載施加方向基本一致。進(jìn)而，通過與圖6對比可以發(fā)現(xiàn)：三維條件下裂紋的擴展過程更為復(fù)雜，會演化生成更多不同類型的裂紋，從而進(jìn)一步證明本文所建立的GMPD模型在模擬準(zhǔn)脆性材料三維裂紋擴展問題方面的準(zhǔn)確性和優(yōu)越性。

圖8 傳統(tǒng)MPPD模型預(yù)測的最終裂紋擴展模式Fig.8 Final growth path obtained from the original MPPD model

4 結(jié)語

提出了一種廣義微極近場動力學(xué)（GMPD）模型，可以有效模擬準(zhǔn)脆性材料的三維斷裂行為。模型采用Timoshenko梁來模擬物質(zhì)點間的相互作用，從而充分考慮鍵在三維受力條件下的軸向變形、切向變形、相對轉(zhuǎn)角以及三者之間的耦合作用。引入了可分別表征鍵的拉伸、剪切以及彎曲剛度的3個鍵參數(shù)，實現(xiàn)了任意變形場下近場動力學(xué)和傳統(tǒng)連續(xù)介質(zhì)力學(xué)的能量一致。同時基于能量提出了一種新型斷裂準(zhǔn)則，給出了鍵軸向變形、切向變形以及相對轉(zhuǎn)角的臨界值，實現(xiàn)準(zhǔn)脆性材料的三維破壞模擬。

基于所提出的GMPD模型模擬了單軸壓縮荷載作用下試件的三維破壞過程，結(jié)果表明該模型可以有效捕捉復(fù)雜荷載條件下不同類型裂紋的萌生和擴展，包括翼型裂紋、反翼型裂紋以及次生裂紋的發(fā)展過程。通過模型預(yù)測結(jié)果與試驗結(jié)果對比，進(jìn)一步驗證了模型的準(zhǔn)確性，表明本模型可以有效應(yīng)用于準(zhǔn)脆性材料的三維斷裂行為模擬，可為工程破壞問題提供分析依據(jù)。

作者貢獻(xiàn)聲明：

陳希卓：模型構(gòu)建、數(shù)據(jù)分析呈現(xiàn)及論文撰寫。

禹海濤：項目負(fù)責(zé)人、論文構(gòu)思、指導(dǎo)模型構(gòu)建及數(shù)據(jù)分析，論文修改。

朱建波：論文修改，提供試驗數(shù)據(jù)。

劉建鋒：論文修改。