基于ABAQUS 的兩種應(yīng)力強(qiáng)度因子計(jì)算方法對(duì)比

姜 旭，寇園園，石朝龍

（1.西安石油大學(xué)石油工程學(xué)院，陜西西安 710065；2.西安石油大學(xué)陜西省油氣井及儲(chǔ)層滲流與巖石力學(xué)重點(diǎn)實(shí)驗(yàn)室，陜西西安 710065）

斷裂力學(xué)[1]在鉆井工程領(lǐng)域中運(yùn)用甚廣，其研究日趨復(fù)雜化、多樣化，由于工程實(shí)踐的復(fù)雜性，大量的解析解和經(jīng)驗(yàn)公式充斥在實(shí)踐的各個(gè)環(huán)節(jié)[2]。表征裂紋尖端應(yīng)力場(chǎng)強(qiáng)度的“應(yīng)力強(qiáng)度因子（SIF）”是井壁裂紋分析的核心內(nèi)容[3]。裂紋尖端應(yīng)力場(chǎng)的奇異性明顯，很大程度上造成了其求解誤差，解析解和經(jīng)驗(yàn)公式都給工程時(shí)間帶來(lái)了巨大的偏差與困擾[4]。因而，建立高效、高精度的計(jì)算方法至關(guān)重要。在ABAQUS 中內(nèi)置有SIF的求解方法，本文將通過(guò)算例對(duì)ABAQUS 內(nèi)置求解方法與傳統(tǒng)“外推法”進(jìn)行對(duì)比，得出若干結(jié)論。

1 基于單元應(yīng)力的外推法計(jì)算SIF 的理論基礎(chǔ)

基于應(yīng)力的外推法是計(jì)算應(yīng)力強(qiáng)度因子最直接的方法，由應(yīng)力強(qiáng)度因子的極限定義式（1）便知，KI是對(duì)應(yīng)于裂尖（r=0）時(shí)的極限值，但是直接計(jì)算無(wú)法使得r=0。為了得到KI，便可使用合理外推的方式來(lái)計(jì)算KI，方法如下：

在運(yùn)用ABAQUS 中易直接查詢(xún)（Inquery）到裂尖前端單元積分點(diǎn)上的應(yīng)力值σy和對(duì)應(yīng)的積分點(diǎn)的坐標(biāo)值r，且其可以被輸出（Output）。同樣可以在繪圖軟件中畫(huà)出σy和r 的關(guān)系，不難得到應(yīng)力分布曲線（見(jiàn)圖1）。這時(shí)，當(dāng)單元被細(xì)化，其應(yīng)力值將趨向于無(wú)窮，文獻(xiàn)中大都稱(chēng)其為應(yīng)力奇異。

圖1 裂紋尖端的應(yīng)力分布示意圖

為了得到KI值，用ri來(lái)表示應(yīng)力積分點(diǎn)與裂尖之間的距離，與其相對(duì)應(yīng)的KIi共同組成數(shù)據(jù)對(duì)（ri，KIi）。在得到數(shù)據(jù)對(duì)之后，運(yùn)用最小二乘法來(lái)擬合。

此時(shí)，令：

設(shè)KI=Ar+B，當(dāng)r=0 時(shí)，KI≈KI（r=0）=B，顯然，易得每個(gè)數(shù)據(jù)點(diǎn)的偏差為KIi-KIi。

根據(jù)最小二乘法，擬合出的最優(yōu)結(jié)果滿足：

對(duì)A 和B 求偏導(dǎo)得到線性方程組，再進(jìn)一步用ri和KIi相關(guān)的因子來(lái)表示A（斜率）和B（截距）的表達(dá)式如下：

因此，只要做出反映σy和r 的關(guān)系曲線，則其截距為應(yīng)力強(qiáng)度因子。

2 ABAQUS 中兩種SIF 計(jì)算方法的實(shí)現(xiàn)[5-8]

2.1 ABAQUS 內(nèi)置方法計(jì)算SIF

ABAQUS 中建立半寬W=100 mm，半高H=200 mm，厚度1 mm 的平板模型（見(jiàn)圖2），該平板中心預(yù)制有一條裂紋，其半裂紋長(zhǎng)度為20 mm；建立彈性模量為2×105MPa，泊松比為0.25 的“彈性各向異性”固體材料；經(jīng)過(guò)裝配、網(wǎng)格劃分、指定裂紋和裂紋擴(kuò)展方向以及裂紋尖端場(chǎng)的奇異性設(shè)置等操作。進(jìn)行邊界條件設(shè)置和應(yīng)力加載。（施加均布拉力30 MPa，下端兩角固定）提交作業(yè)后，裂紋尖端應(yīng)力（見(jiàn)圖3）。

圖2 裂紋平板模型

圖3 Job 提交后的裂紋尖端沿Y 方向的正應(yīng)力分布

在后處理中輸出的Data 文件中，ABAQUS 內(nèi)置方法計(jì)算出的SIF 值（見(jiàn)圖4）。

圖4 Data 文件中輸出的SIF 計(jì)算結(jié)果

2.2 ABAQUS 中外推法的實(shí)現(xiàn)過(guò)程

在后處理中，運(yùn)用“單元應(yīng)力的外推法”，通過(guò)定義“path”輸出（20，0）處至（25，0）處對(duì)應(yīng)的應(yīng)力值，輸出rpt 文件。將該文件中的數(shù)據(jù)導(dǎo)入表格中進(jìn)行處理（見(jiàn)表1）。

表1 基于rpt 文件的數(shù)據(jù)處理

表1 基于rpt 文件的數(shù)據(jù)處理（續(xù)表）

基于表1 中數(shù)據(jù)繪制KI的變化趨勢(shì)圖（見(jiàn)圖5）。其截距238.26 即為KI值，該值與ABAQUS 求解器所求平均值244.13 的相對(duì)誤差僅為2.46%。

圖5 KI 隨ri 變化的趨勢(shì)圖

略去表1 中前21 組數(shù)據(jù)，排除裂尖的奇異性干擾?？汕蟮酶咏咏腒I值：243.36（見(jiàn)圖6），其與ABAQUS 求解器所求平均值244.13 的相對(duì)誤差縮至0.32%，相關(guān)性系數(shù)從0.158 2 上升至0.989 2。

圖6 KI 隨ri 變化的趨勢(shì)圖

3 ABAQUS 求解器求解與外推法的求解值的對(duì)比

基于以上計(jì)算過(guò)程，進(jìn)行對(duì)比（見(jiàn)表2）。

表2 SIF 求解結(jié)果對(duì)比

4 結(jié)論

（1）ABAQUS 求出的值與應(yīng)力外推法的求解值相近，說(shuō)明了ABAQUS 內(nèi)置的算法已經(jīng)相當(dāng)精確。

（2）在削弱裂尖部位的奇異性影響后，可使得其與ABAQUS 內(nèi)置算法逼近，說(shuō)明ABAQUS 的計(jì)算精度是可靠的。數(shù)據(jù)的剔除過(guò)程隨意性較大，客觀上造成計(jì)算結(jié)果的不確定性。

（3）基于應(yīng)力的兩次外推過(guò)程，其差異化結(jié)果印證了裂紋尖端場(chǎng)所具有的奇異性特點(diǎn)，去掉裂紋尖端附近的點(diǎn)可有力減弱奇異性對(duì)應(yīng)力強(qiáng)度因子的影響。

（4）裂紋尖端應(yīng)力強(qiáng)度因子的研究對(duì)象是“裂紋尖端”，而實(shí)際計(jì)算中卻因尖端數(shù)據(jù)不收斂而舍棄了該部分?jǐn)?shù)據(jù)，這造成了邏輯上的悖論。