基于二次插值重構(gòu)有限元法的動態(tài)裂紋擴展模擬

王玉光 +李忠文 吳圣川 臧曉蕾 王燕

摘要：

基于二次插值重構(gòu)有限元法（Twiceinterpolation Finite Element Method， TFEM）分析動態(tài)斷裂力學(xué)問題并進行數(shù)值實驗，考察TFEM在裂紋動態(tài)擴展模擬中的準(zhǔn)確性和可靠性.由于TFEM保證節(jié)點處梯度場的連續(xù)性，因此裂尖附近的應(yīng)力場可以得到較好的逼近.把該算法成功移植到自主開發(fā)的三維裂紋擴展仿真軟件（ZonCrack）中.利用ZonCrack進行的裂紋擴展，分析結(jié)果表明：TFEM得到裂尖應(yīng)力強度因子（Stress Intensity Factor， SIF）與解析解基本一致；裂紋擴展的模擬結(jié)果與實驗值吻合良好.

關(guān)鍵詞：

插值重構(gòu)有限元； 斷裂力學(xué)； 疲勞裂紋擴展； 裂紋尖端場； 損傷容限； ZonCrack軟件

中圖分類號： TB115；O242.21

文獻標(biāo)志碼： B

0 引 言

動態(tài)疲勞裂紋的模擬分析對工程結(jié)構(gòu)服役安全性和剩余壽命評價有重要的工程指導(dǎo)意義.[1]為模擬裂紋擴展，一般采用節(jié)點松弛方法，即在計算中把節(jié)點沿著裂縫路徑一分為二.[2]由于這種技術(shù)適用性差，學(xué)者們便又提出一些新型計算方法，如無網(wǎng)格[3]、邊界單元法[4]，以及有效的單元分割策略[5].這些方法雖然可在一定程度上提高計算精度，但計算效率卻降低，顯然不利于工程應(yīng)用.

為此，ZHENG等[6]提出一種二次插值重構(gòu)有限元法（Twiceinterpolation Finite Elements Method， TFEM）求解線彈性斷裂力學(xué)問題.為提高缺陷體的分析精度，通過2次連續(xù)插值得到形函數(shù)，其主要特點是把標(biāo)準(zhǔn)有限元法（Finite Elements Method， FEM）插值函數(shù)作為第二階段的插值參數(shù)，或者說把傳統(tǒng)FEM不連續(xù)的節(jié)點梯度（節(jié)點導(dǎo)數(shù)）作為插值參數(shù)，這樣，新的形函數(shù)便具有連續(xù)的節(jié)點梯度場.結(jié)果發(fā)現(xiàn)，TFEM比經(jīng)典FEM有更好的計算精度和收斂性.

該算法植入新一代基于自適應(yīng)加密的擴展有限元（eXtend FEM，XFEM）的斷裂仿真軟件ZonCrack.ZonCrack基于概率疲勞斷裂力學(xué)理論，對含裂紋工程的結(jié)構(gòu)進行裂紋擴展仿真，獲取無損探傷數(shù)據(jù)，為工程師提供科學(xué)的基礎(chǔ)數(shù)據(jù)和判據(jù).該算法能夠預(yù)測靜載荷或疲勞載荷作用下的開裂行為，確定工程結(jié)構(gòu)損傷容限，不僅可以進行傳統(tǒng)有限元強度分析，而且提供疲勞斷裂及擴展壽命的分析功能.

1 基于平面三角形的TFEM

選取對復(fù)雜幾何形狀適應(yīng)性好的三角形單元構(gòu)建插值函數(shù)以求解斷裂力學(xué)問題.

1.1 一次插值函數(shù)

給定一個典型的常應(yīng)變的單純形三角形，單元由節(jié)點i，j，m逆時針排列構(gòu)成.插值重構(gòu)法在二維單純形單元中的示意見圖1.

如圖1所示，Si，Sj和Sm為有公共節(jié)點i，j和k的關(guān)系單元集合，組成支撐域.對u（x）有貢獻的節(jié)點集（無單元法中一般稱為支撐點）正是Si，Sj和Sm所涵蓋的全部節(jié)點.設(shè)qs為支撐節(jié)點的位移向量，

1.2 二次插值函數(shù)

根據(jù)經(jīng)典FEM得到節(jié)點導(dǎo)數(shù)后，再進入二次插值獲取一個新的形函數(shù)

1.3 關(guān)于C0節(jié)點

從式（16）看出，二次插值得到的試函數(shù)在節(jié)點處具有C1連續(xù)性.這一特點具有提高計算精度的優(yōu)勢，但位于材料界面和位移邊界上的節(jié)點給求解帶來困難.為此，設(shè)節(jié)點i具有C0連續(xù)性，則對于插值點x位于單元e內(nèi)的節(jié)點，可令

分析發(fā)現(xiàn)，當(dāng)問題域內(nèi)所有節(jié)點進行C0連續(xù)性處理后，

TFEM型函數(shù)將退化為經(jīng)典FEM型函數(shù)，即在某些情況下2種算法可以耦合.也就是說，對于含缺陷體需要更高插值精度的局部區(qū)域可采用TFEM，而在遠離裂紋前緣區(qū)則可采用一般精度FEM，從而實現(xiàn)計算精度和效率的平衡考慮.

2 裂紋擴展計算關(guān)鍵技術(shù)

發(fā)展TFEM方法用于求解斷裂力學(xué)問題，并植入ZonCrack軟件中.專業(yè)斷裂仿真軟件ZonCrack是一個面向動態(tài)斷裂過程的仿真實驗平臺，其以成熟的疲勞斷裂力學(xué)為理論基礎(chǔ)，對含缺陷工程部件進行仿真模擬.[7]為減少對經(jīng)典FEM源代碼的改動，引入節(jié)點松弛技術(shù)模擬材料的劈開過程.

2.1 應(yīng)力強度因子

采用區(qū)域互作用積分計算混合型裂紋的應(yīng)力強度因子（Stress Intensity Factor， SIF）解[8]，以最大周向應(yīng)力判斷裂紋方向[9]，則極坐標(biāo)形式的裂尖應(yīng)力場分量可有I型和II型裂紋應(yīng)力強度因子表示

2.2 節(jié)點投影線技術(shù)

研究發(fā)現(xiàn)，網(wǎng)格質(zhì)量對插值重構(gòu)有限元法的計算精度影響較小（網(wǎng)格敏感性低），也就是說調(diào)整裂紋面附近的網(wǎng)格對SIF的影響也較小.

當(dāng)前求解斷裂力學(xué)問題的通用處理方法是確保裂紋面/線與單元邊/節(jié)點的重合.本文采用節(jié)點投影方法準(zhǔn)確模擬裂紋的擴展路徑，具體操作過程為：網(wǎng)格整體上不調(diào)整，僅對裂紋附近節(jié)點投影到裂紋面/線上（見圖2）.一般來說，裂紋與三角元有3種關(guān)系：第一，裂紋與2個邊相交，把相交邊的公共節(jié)點移至裂紋；第二，裂紋與1個邊相交，把離裂紋近的節(jié)點移至裂紋；第三，裂尖位于單元內(nèi)，則移動單

元上距裂紋尖端最近的節(jié)點.

采用上述處理后，問題域內(nèi)所有單元邊/節(jié)點都會位于裂紋面/線上.自此，采用經(jīng)典節(jié)點分離法人為把裂紋面/線上的節(jié)點一分為二，并賦予不同節(jié)點號.通過這樣處理，裂紋擴展路徑上就出現(xiàn)一個強不連續(xù)性位移，即實現(xiàn)裂紋張開模擬.

3 算例驗證

為考察新型二次插值有限元的求解精度和穩(wěn)定性，定義2范數(shù)位移誤差和能量誤差為

式中：下標(biāo)e表示精確解或者理論解；下標(biāo)n表示數(shù)值方法參考解或者TFEM解.在TFEM實現(xiàn)中，單元剛度矩陣積分采用4點Hammer積分方法.

在ZonCrack中調(diào)用TFEM求解器分析一些典型的斷裂力學(xué)問題，以測試TFEM法準(zhǔn)確性與合理性.文中所用變量均采用國際單位制.

3.1 單邊直裂紋

一個單邊直裂紋的無限大平板，兩端受均勻拉力，假設(shè)平面應(yīng)變狀態(tài).模型的幾何及加載見圖3，且有H=L=P=1.力學(xué)性能參數(shù)E=1010，ν=0.25.解析解為

應(yīng)力強度因子解對比曲線見圖4.由此可知，對于不同擴展長度的裂紋，TFEM數(shù)值解與理論解和經(jīng)典FEM解吻合良好，且TFEM法的計算精度高于傳統(tǒng)FEM法.TFEM與傳統(tǒng)FEM的裂紋尖端應(yīng)力場分布下見圖5.由此可見，TFEM結(jié)果與理論解給出的裂紋簡單應(yīng)力場吻合更好，而傳統(tǒng)有限元法卻存在著明顯的應(yīng)力跳躍現(xiàn)象（即不連續(xù)）.

3.2 中心斜裂紋

為進一步考察TFEM求解SIF的準(zhǔn)確性與可靠性，以中心斜裂紋板為例研究，見圖6.在無限大平板x方向和y方向分別搜到σ1和σ2大小的拉應(yīng)力，β為斜裂紋與x軸的夾角，顯然，KI和KII為關(guān)于裂紋傾角β的函數(shù)，理論解為

有限元計算選取的計算域為邊長L=10的方板，所含裂紋的長度a=0.25，假設(shè)軸向拉力P=1.考慮到裂紋長度遠小于矩形板的長度，所以可以視為無限大平板的中心斜裂紋問題.

3.3 4點剪切含裂紋矩形板

4點剪切受載問題的裂紋擴展研究較多，幾何模型及受載情況見圖7，L=4，b=1，P=1.材料參數(shù)為E=200 MPa，ν=0.25，假設(shè)平面應(yīng)變條件.第12擴展步后的裂紋路徑見圖8.由此可見，插值重構(gòu)法的計算結(jié)果與實驗基本吻合.

4 結(jié) 論

將新型TFEM植入高級三維裂紋擴展仿真與分析系統(tǒng)ZonCrack中用于求解動態(tài)斷裂問題.通過標(biāo)準(zhǔn)算例分析，得出以下結(jié)論.

（1）在ZonCrack中調(diào)用TFEM求解器，能夠獲得與理論解可比的I型、II型及混合型SIF解，且在同一單元網(wǎng)格下，計算精度高于傳統(tǒng)FEM.

（2）與傳統(tǒng)FEM解相比較，TFEM可在裂尖區(qū)域得到更加準(zhǔn)確和光滑連續(xù)的應(yīng)力分布.

（3）TFEM法預(yù)測的裂紋路徑與實測解基本一致，表明TFEM算法的可靠性.

必須指出，與傳統(tǒng)有限元比較，雖然插值重構(gòu)有限元法使用經(jīng)典FEM的導(dǎo)數(shù)進行二次插值，但并未因此顯著增加問題的總求解自由度.

參考文獻：

[1] 楊小彬， 莊茁， 莊傳晶. 富氣輸送管道裂紋動態(tài)擴展的數(shù)值模擬[J]. 清華大學(xué)學(xué)報（自然科學(xué)版）， 2008， 48（8）： 13551358. DOI： 10.3321/j.issn：10000054.2008.08.032.

YANG X B， ZHUANG Z， ZHUANG C J. Numerical study of crack propagation in rich gas transporting pipelines[J]. Journal of Tsinghua University （Science and Technology）， 2008， 48（8）： 13551358. DOI： 10.3321/j.issn：10000054.2008.08.032.

[2] 解德， 錢勤， 李長安. 斷裂力學(xué)中的數(shù)值計算方法及工程應(yīng)用[M]. 北京： 科學(xué)出版社， 2009.

[3] BELYTSCHKO T， LU Y Y， GU L. Crack propagation by elementfree Galerkin methods[J]. Engineering Fracture Mechanics， 1995， 51（2）： 295315. DOI： 10.1016/00137944（94）001539.

[4] 李俊， 馮偉哲， 高效偉. 一種基于直接計算高階奇異積分的斷裂力學(xué)雙邊界積分方程分析法[J]. 力學(xué)學(xué)報， 2016， 48（2）： 387398.DOI： 10.6052/0459187915342.

LI J， FENG W Z， GAO X W. A dual boundary integral equation method based on direct evaluation of higher order singular integral for crack problems[J]. Chinese Journal of Theoretical and Applied Mechanics， 2006， 48（2）： 387398.DOI： 10.6052/0459187915342.

[5] MELENK J M， BABUS K. The partition of unity finite element method： basic theory and applications[J]. Computer Methods in Applied Mechanics and Engineering， 1996， 139（14）： 289314. DOI： 10.1016/S00457825（96）010870.

[6] ZHENG C， WU S C， TANG X H， et al. A novel twice interpolation finite element method for solid mechanics problems[J]. Acta Mechanica Sinica， 2010， 26（2）： 265278. DOI： 10.1007/s1040900902653.

[7] 吳圣川， 吳玉程. 斷裂分析及CAE軟件的現(xiàn)狀與發(fā)展[J]. 計算機輔助工程， 2011， 20（1）： 12.

WU S C， WU Y C. Current status and development of fracture mechanics based CAE software[J]. Computer Aided Engineering， 2011， 20（1）： 12.

[8] BUDIANSKY B， RICE J R. Conservation laws and energyrelease rates[J]. Journal of Applied Mechanics， 1973， 40（2）： 201203. DOI： 10.1115/1.3422926.

[9] FLEMING M， CHU Y A， MORAN B， et al. Enriched elementfree Galerkin methods for cracktip fields[J]. International Journal for Numerical Methods in Engineering， 1997， 40（8）： 14831504. DOI： 10.1002/（SICI）10970207（19970430）40：8<1483：：AIDNME123>3.0.CO；26.

（編輯 武曉英）