基于X-SBFEM的非線性斷裂數值模型研究

傅 興 安, 李 建 波*, 林 皋

( 1.大連理工大學 海岸和近海工程國家重點實驗室， 遼寧 大連 116024；2.大連理工大學 建設工程學部 工程抗震研究所， 遼寧 大連 116024 )

基于X-SBFEM的非線性斷裂數值模型研究

傅 興 安1,2, 李 建 波*1,2, 林 皋1,2

( 1.大連理工大學 海岸和近海工程國家重點實驗室， 遼寧 大連 116024；2.大連理工大學 建設工程學部 工程抗震研究所， 遼寧 大連 116024 )

從非線性斷裂力學模型的角度，開展準脆性材料(混凝土)裂紋過程區的有效模擬，是當前的研究熱點之一．擴展比例邊界有限元法(X-SBFEM)兼有擴展有限元法(XFEM)和比例邊界有限元法(SBFEM)兩種方法的優勢，利用SBFEM求解裂尖段應力奇異性問題，利用XFEM模擬非裂尖段位移場不連續．為在X-SBFEM中增加非線性斷裂模型，提出采用side-face力的形式，基于黏聚力模型，通過線性疊加迭代法來模擬準脆性材料(混凝土)裂紋過程區．最后，以數值算例——單邊缺口的三點彎曲梁和四點剪切梁——模擬裂紋過程區能量耗散影響，驗證了所提方法的精度與應用效果．

擴展比例邊界有限元法；應力強度因子；黏聚力；side-face力；裂紋過程區

0 引 言

材料的斷裂行為以及斷裂過程的數值模擬是評估結構承載能力的重要依據．近幾年從斷裂力學角度提出的擴展比例邊界有限元法(X-SBFEM)是在擴展有限元法[1-2](XFEM)和比例邊界有限元法[3-4](SBFEM)的基礎上發展而來的．X-SBFEM[5-6]兼有兩種方法的優勢，在裂紋主體段用XFEM模擬，在裂尖段用SBFEM模擬，最終建立裂紋體的總體平衡方程并進行求解．其中，XFEM的虛自由度與SBFEM的全實自由度的銜接是重點環節．X-SBFEM克服了擴展有限元在域內裂紋尖端近場需要增強函數來描述，且某些情況下，該增強函數將是非多項式或不連續的，造成在形成剛度矩陣時需要特殊處理．而進一步考慮模型裂紋過程區(FPZ)的非線性影響，在X-SBFEM的基礎上，引入非線性斷裂力學，對于推進模型的應用具有重要工程意義．一般的，準脆性材料裂紋的模擬主要有線彈性的方法(LEFM)和非線性的方法(NFM)，用LEFM還是NFM取決于FPZ相對于結構的尺寸[7]．現在線彈性力學只能應用在大尺度結構上(如混凝土大壩等)的觀點已被廣泛接受．為得到更高的計算精度，對于正常尺寸的結構(如混凝土梁)應用基于FPZ的非線性斷裂力學模型是有必要的[7]．

Yang等[8]、Shi等[9]近年采用FEM-SBFEM耦合的方法或P-SBFEM模擬FPZ的非線性行為時，都通過插入黏聚界面單元(CIEs)來進行分析．這種方法在模擬裂紋擴展過程中，新生成裂紋邊界變得扭曲復雜,對于復雜結構模型應用起來有一定的難度[7]．

本文提出在X-SBFEM基礎上采用線性漸進疊加假設的迭代法求解裂紋面黏聚力來考慮裂紋過程區的非線性影響，并通過對裂紋擴展問題采用基于LEFM的最大環向拉應力準則來預測裂紋路徑，模擬裂紋尖端FPZ的能量耗散．

1 擴展比例邊界有限元法(X-SBFEM)基本理論

X-SBFEM的主要工作是通過建立裂尖超單元利用SBFEM半解析的性質來模擬裂尖近場的非連續行為，而在裂紋非裂尖段則采用XFEM來描述，關鍵是單元接觸界面處的數值處理．圖1是模型域內拓撲關系的示意圖[6]．

圖1 X-SBFEM含裂紋模型域內拓撲關系

1.1 擴展有限元法(XFEM)原理

基于單元分解思想，裂紋的主體段用擴展有限元法進行模擬，其位移場的一般形式為[10]

(1)

式中：Nfem、Nc分別為普通有限元、因被內部裂紋貫穿而增強的節點．N表示形函數．qI為普通自由度，而aJ為與Q所關聯的廣義自由度(如圖1中方形節點所示)，Q(x)為Heaviside階躍函數．擴展有限元平衡方程為

(2)

1.2 比例邊界有限元法(SBFEM)原理

如圖2所示，考慮裂尖SBFEM超單元內裂紋面有side-face力．首先在不考慮體力的情況下，SBFEM的位移場和應力場為[3]

(3)

(4)

式中：N(η)為邊界的形函數，φi和λi分別為特征值及位移模態，ci表示積分常數，D為彈性模量矩陣，L為線性算子矩陣，B1(η)和B2(η)決定于SBFEM超單元的幾何特性．根據虛功原理可推得：

(5)

(6)

式中：P為超單元邊界等效節點力，Φ為由位移模態為列所組成的矩陣，uh為邊界節點位移，E0、E1和E2分別為位移控制方程的系數陣．于是應力強度因子的公式可表達為[6,12]

(7)

式中：KⅠ、KⅡ為兩型應力強度因子，L0為裂尖到沿裂紋面延長線與邊界相交點的距離．ψyy、ψxy為應力模態．

圖2 SBFEM有限域的邊界離散及相似坐標變換

Fig.2 Discretization on boundary with element and scaled transformation of coordinates of bounded media of SBFEM

1.3 X-SBFEM耦合模型邊界模式

圖3描述了一個典型的用XFEM和SBFEM模擬裂紋的方法．SBFEM超單元直接通過單元兩條邊界來描述裂尖裂紋面，XFEM則由階躍函數引入附加的自由度來描述裂紋非裂尖段非連續位移場[6]．為了保證兩種不同單元間的位移連續，則需引入特殊的矩陣T使SBFEM和XFEM銜接節點的位移協調．在把SBFEM得到的未知節點位移(uE、uF、uA和uB)轉換成XFEM未知節點位移(q2、q3、a2和a3)的過程中可以根據之前的公式推導出一個轉換矩陣，如下式：

(8)

其中I為單位矩陣．為保證位移協調并且把單剛矩陣集入總體剛度陣，需把SBFEM超單元的位移向量和剛度矩陣重新排列，如下式：

(9)

其中usb為SBFEM單元在非共用邊界上的節點位移，uxf為XFEM單元在共用邊界上的節點位移．轉換矩陣T只與SBFEM域和XFEM域共用邊界上的裂紋開口處的形函數有關．K表示剛度矩陣，其中下標a表示SBFEM域和XFEM域單元之間共用節點自由度，下標b表示非共用節點自由度．

圖3 擴展有限元域與比例邊界有限元域的耦合

2 非線性斷裂力學模型

2.1 線性漸進疊加假設

Hillerborg等[13]提出的P-a曲線(external load-crack mouth opening displacement curve)是由混凝土結構的裂紋過程區(FPZ)的黏性行為引起的．為簡化Xu等[14]求解黏性裂紋擴展的方法，基于線性漸進疊加假設[7]，把FPZ當作線彈性裂紋進行求解．

這種線性簡化方法中假設塑性變形為0，忽略圖4(a)中abp曲線非線彈性部分的影響，即把圖4(a)簡化成圖4(b)．由此，圖4(b)中B點是線彈性點．這樣線彈性斷裂力學將可以適用于虛擬的OB段加載的分析．由此，相同材料及幾何特性條件下，不同預設裂紋長度的一系列構件將會得到一系列類似的B點．此時，P-a曲線可以當成這一系列B點的包絡線．

2.2 開裂準則

混凝土開裂往往呈現復合斷裂形態，從數值角度，結構所承受的總荷載如圖5(c)所示，除承受外荷載外，在裂紋過程區內，如裂紋面之間相對的位移，包括裂紋面的張裂位移(COD)和滑裂位移(CSD)，沒有超過圖6所示的限值wc和sc時，則還承受虛擬裂紋面內的黏聚力．但裂紋面相對位移超過限值的部分令黏聚力為0．對于有黏聚力的情況，以Ⅰ型應力強度因子為例，應力強度因子由兩部分組成[14]：

(a) 真實卸載/重加載軌跡

(b) 虛擬卸載/重加載軌跡

圖4 線性漸進疊加假設

(a) 只有外力作用

(b) 只有黏聚力作用

(c) 外力與裂紋面黏聚力疊加作用

圖5 疊加法計算應力強度因子KⅠ,Ⅱ

Fig.5 Superposition method for calculating stress intensity factorKⅠ,Ⅱ

(10)

2.3 FPZ軟化曲線

裂紋過程區內的正黏聚力σ(x)決定于圖6(a)σ-w雙線性的軟化曲線或圖6(b)σ-w單線性的軟化曲線．切分量τ(x)由圖6(c)τ-s曲線決定．圖6(a)(b)曲線以下的面積為Ⅰ型斷裂能GfI．圖6(c)曲線間的面積為Ⅱ型斷裂能GfⅡ的2倍[8]．假設圖6中曲線關系都為不可逆的卸載路徑，即當裂紋面閉合時，一種彈性的卸載將隨著w或s的減小而沿著圖中的割線發生．實心箭頭表示加載，空心箭頭表示卸載．

(a) σ-w雙線性關系

(b)σ-w單線性關系

(c)τ-s關系

圖6 裂紋相對位移與黏聚力的關系

Fig.6 Relations between relative displacements of crack and cohesive tractions

2.4 以迭代形式模擬黏聚力模型

本文工作的關鍵是根據裂紋面的相對位移采用線性疊加的迭代法求解裂紋面黏聚力．

步驟1假設結構只受外力F作用，這樣可以根據X-SBFEM基于線彈性假定求出超單元裂紋面相對位移Δui，再根據圖6即可得出相應黏聚力ti．

步驟2如圖5，把外荷載F和步驟1中得到的黏聚力ti一起作用在結構上．其中黏聚力以side-face力形式施加在裂紋面上．

side-face力分布荷載[5]對應的虛功項為

(11)

易得SBFEM非齊次控制方程：

(12)

假定荷載及位移模式能由冪級數表達，即

(13)

(14)

將式(14)代入式(12)和式(13)可得

(15)

(16)

因此，邊界節點完整位移和等效節點力分別為

(17)

(18)

其中Φ和Q分別為式(6)解出的位移模態矩陣和力的模態矩陣．綜合式(17)和式(18)可得

(19)

故由side-face力形式加載的黏聚力產生的SBFEM等效邊界節點力為

(20)

(21)

則可得到SBFEM超單元的位移場為

(22)

式中：φt為SBFEM解出的應力模態．由式(22)求解出相對位移Δui+1．

步驟3重復步驟1、2直至ti與Δui+1的關系能很好地符合圖6中的關系．

3 算例與分析

3.1 單邊缺口三點彎曲梁

Hillerborg等曾通過試驗研究單邊缺口三點彎曲梁模型[13]，梁的材料和幾何尺寸數據等如圖7所示．抗拉強度ft=3.33 MPa，Ⅰ型斷裂能GfⅠ=137 N/m．本算例基于LEFM最大環向拉應力準則來預測裂紋擴展路徑．單線性軟化曲線(圖6(b)) 被用來做分析．由Ⅰ型斷裂能GfⅠ所得出的單線性軟化曲線限值為wc=0.082 3 mm．計算了在20×200網格密度條件下3種裂紋擴展步長10、20和30 mm的結果，并與基于線彈性方法[6]的結果進行對比．

圖7 單邊缺口三點彎曲梁

圖8是考慮FPZ非線性的3種不同裂紋擴展步長時的荷載-加載點位移曲線．由圖中可以看出，陳白斌等[6](基于X-SBFEM) 和Yang等[8](基于SBFEM)基于LEFM結果數據與Hillerborg等的試驗數據[13]有很大的差別，尤其是他們的荷載峰值都比試驗的峰值高出較多．這是因為基于線彈性斷裂力學方法無法模擬FPZ的能量耗散．而從圖中可以看出本文方法計算得到的結果與試驗結果[13]是較為符合的，說明本文中所使用的基于X-SBFEM 通過線性疊加的迭代法求解黏聚力的方法能很好地模擬FPZ的能量耗散．而且從圖中可以看出3種裂紋擴展步長的結果都很符合試驗曲線，這說明不同裂縫擴展步長對計算結果影響不大．

圖8 基于非線性斷裂力學不同裂紋擴展步長時三點彎曲梁的荷載-加載點位移曲線

Fig.8 NFM-based load-LPD curves for three-point bending beam with different crack increment lengths

3.2 單邊缺口四點剪切梁

Arrea和Ingraffea曾試驗和分析了單邊缺口四點剪切梁模型[15]．其幾何尺寸及材料參數如圖9所示．假設為平面應力狀態．抗拉強度ft=3.00 MPa，Ⅰ型斷裂能GfⅠ=100 N/m，Ⅱ型斷裂能GfⅡ=10 N/m．裂紋路徑預測基于LEFM最大環向拉應力準則．單線性的軟化曲線(圖6(b))和τ-s曲線(圖6(c))被用來做分析．由Ⅰ型斷裂能GfⅠ=100 N/m計算得到的COD限值wc=0.067 mm，由Ⅱ型斷裂能GfⅡ求得的CSD限值sc=0.02 mm．計算了在20×200網格密度條件下3種裂紋擴展步長20、30和40 mm的結果，并與基于線彈性方法[6]的結果進行對比．

圖9 單邊缺口四點剪切梁

圖10和圖11分別展示了基于本文方法計算得到的不同裂紋擴展步長下荷載與裂紋口滑移位移的關系曲線、荷載與其加載點位移關系的曲線．其中，由圖10可以看出，本文計算結果與Arrea和Ingraffea試驗得到的數據[15]及Yang等[8](基于NFM)的數值解擬合較好；從圖11可以看出，本文計算結果與Yang等[8](基于NFM)的數值解較為接近，峰值都明顯低于陳白斌等[6](基于X-SBFEM、LEFM)所得的數值解結果，說明本文方法較好地基于X-SBFEM算法模擬了裂紋尖端存在FPZ能量耗散的影響．且由圖看出3種裂紋擴展步長的結果接近，即裂紋擴展步長對計算結果的影響不是很大，從而驗證了X-SBFEM算法對于復合型裂紋擴展問題的適用性．此外，圖11可看成本文方法較好地表現出荷載-加載點位移曲線snap-back的現象．圖12所示為當a=20 mm時F=36.7 kN作用下的裂紋軌跡及正黏聚力沿裂紋面分布圖．

圖10 基于非線性斷裂力學不同裂紋擴展步長下的荷載-滑移位移曲線

Fig.10 NFM-based load-CMSD curves with the change of crack increment length

圖11 基于非線性斷裂力學不同裂紋擴展步長時四點剪切梁的荷載-加載點位移曲線

(a) 變形幾何圖

4 結 語

本文提出基于擴展比例邊界有限元法(X-SBFEM)采用線性疊加迭代法模擬準脆性材料(混凝土)裂紋過程區的非線性行為，并將其結果與基于X-SBFEM采用LEFM方法得到的結果進行對比，可以看出其模擬裂紋尖端FPZ的能量耗散的有效性．本文發展了X-SBFEM在非線性斷裂過程模擬中的應用．此外，采用該方法不需要插入CIEs和引入影子域(shadow domain)就可以得到較高精度的結果．

[1] RABCZUK T, BORDAS S, ZI G. On three-dimensional modelling of crack growth using partition of unity methods [J].ComputersandStructures, 2010,88(23/24):1391-1411.

[2] BELYTSCHKO T, GRACIE R, VENTURA G. A review of extended/generalized finite element methods for material modeling [J].ModellingandSimulationinMaterialsScienceandEngineering, 2009,17(4):043001.

[3] SONG Chongmin, WOLF J P. Semi-analytical representation of stress singularities as occurring in cracks in anisotropic multi-materials with the scaled boundary finite element method [J].ComputersandStructures, 2002,80(2):183-197.

[4] SONG Chongmin, TIN-LOI F, GAO Wei. A definition and evaluation procedure of generalized stress intensity factors at cracks and multi-material wedges [J].EngineeringFractureMechanics, 2010,77(12):2316-2336.

[5] 陳白斌,李建波,林 皋. 無需裂尖增強函數的擴展比例邊界有限元法[J]. 水利學報, 2015,46(4):489-496, 504.

CHEN Baibin, LI Jianbo, LIN Gao. An extended scaled boundary finite element method without asymptotic enrichment of the crack tip [J].JournalofHydraulicEngineering, 2015,46(4):489-496,504. (in Chinese)

[6] 陳白斌,李建波,林 皋. 基于X-SBFEM的裂紋體非網格重剖分耦合模型研究[J]. 工程力學, 2015,32(3):15-21.

CHEN Baibin, LI Jianbo, LIN Gao. Study on the coupling model of crack without remeshing based on X-SBFEM [J].EngineeringMechanics, 2015,32(3):15-21. (in Chinese)

[7] 朱朝磊. 基于比例邊界有限元方法的混凝土結構靜動態斷裂模擬[D]. 大連:大連理工大學, 2014.

ZHU Chaolei. Static and dynamic fracture simulation of concrete structures based on scaled boundary finite element method [D]. Dalian:Dalian University of Technology, 2014. (in Chinese)

[8] YANG Z J, DEEKS A J. Fully-automatic modelling of cohesive crack growth using a finite element-scaled boundary finite element coupled method [J].EngineeringFractureMechanics, 2007,74(16):2547-2573.

[9] SHI Mingguang, ZHONG Hong, OOI E T,etal. Modelling of crack propagation of gravity dams by scaled boundary polygons and cohesive crack model [J].InternationalJournalofFracture, 2013,183(1):29-48.

[10] 茹忠亮,朱傳銳,趙洪波. 裂紋擴展問題的改進XFEM算法[J]. 工程力學, 2012,29(7):12-16.

RU Zhongliang, ZHU Chuanrui, ZHAO Hongbo. An improved algorithm of XFEM for the crack propagation problems [J].EngineeringMechanics, 2012,29(7):12-16. (in Chinese)

[11] 張曉東,丁 勇,任旭春. 混凝土裂紋擴展過程模擬的擴展有限元法研究[J]. 工程力學, 2013,30(7):14-21.

ZHANG Xiaodong, DING Yong, REN Xuchun. Simulation of the concrete crack propagation process with the extended finite element method [J].EngineeringMechanics, 2013,30(7):14-21. (in Chinese)

[12] 劉鈞玉,林 皋,胡志強. 裂紋面荷載作用下多裂紋應力強度因子計算[J]. 工程力學, 2011,28(4):7-12.

LIU Junyu, LIN Gao, HU Zhiqiang. The calculation of stress intensity factors of multiple cracks under surface tractions [J].EngineeringMechanics, 2011,28(4):7-12. (in Chinese)

[13] HILLERBORG A, MODEER M, PETERSSON P E. Analysis of crack formation and crack growth in concrete by means of fracture mechanics and finite elements [J].CementandConcreteResearch, 1967,6(6):773-781.

[14] XU Shilang, REINHARDT H W. Determination of double-Kcriterion for crack propagation in quasi-brittle fracture. Part Ⅱ: Analytical evaluating and practical measuring methods for three-point bending notched beams [J].InternationalJournalofFracture, 1999,98(2):151-177.

[15] ARREA M, INGRAFFEA A. Mixed-mode crack propagation in mortar and concrete: 81-13 [R]. Ithaca: Department of Structural Engineering, Cornell University, 1982.

StudyofnumericalmodelofnonlinearfracturebasedonX-SBFEM

FU Xing′an1,2, LI Jianbo*1,2, LIN Gao1,2

( 1.State Key Laboratory of Coastal and Offshore Engineering, Dalian University of Technology, Dalian 116024, China； 2.Institute of Earthquake Engineering, Faculty of Infrastructure Engineering, Dalian University of Technology, Dalian 116024, China )

The effective simulation of the fracture process zone (FPZ) of quasi-brittle material (concrete) by using nonlinear fracture mechanics model is one of the hot spots in the current research. Extended scaled boundary finite element method (X-SBFEM) has advantages of both extended finite element method (XFEM) and scaled boundary finite element method (SBFEM), making full use of XFEM to describe discontinuous displacement field in non-crack-tip area and SBFEM to precisely solve stress singular problems in the crack-tip area. For the application of nonlinear fracture models in X-SBFEM, based on the model of the cohesive traction, the FPZ of the quasi-brittle material (concrete) is simulated using the iterative method of linear superposition by the form of side-face traction. Finally, numerical examples of the three-point single notched bending beam and the four-point single notched shear beam are used to simulate the energy dissipation effect in the FPZ, which verify the accuracy and application effect of the proposed method.

extended scaled boundary finite element method (X-SBFEM); stress intensity factor; cohesive traction; side-face traction; fracture process zone (FPZ)

2016-09-18；

2017-05-26．

國家重點研發計劃資助項目(2016YFB0201000)；國家自然科學基金資助項目(51779222)；國家自然科學基金委創新研究群體項目(51421064)；中央高校基本科研業務費專項資金資助項目(DUT17LK16).

傅興安(1991-)，男，碩士生，E-mail:fuxingan@mail.dlut.edu.cn；李建波*(1977-)，男，博士，副教授，E-mail:jianboli@dlut.edu.cn;林 皋(1929-)，男，教授，中國科學院院士，E-mail:gaolin@dlut.edu.cn.

1000-8608(2017)05-0494-07

TV313

A

10.7511/dllgxb201705009