三維斷裂彈性動力學的改進無單元Galerkin法

景永強 彭妙娟

摘要：為提高斷裂彈性動力學問題數(shù)值計算的精度，避免出現(xiàn)病態(tài)或奇異方程組，基于改進的移動最小二乘法建立三維彈性動力學問題的積分弱形式，采用罰函數(shù)法施加位移邊界條件，引入隱式時間積分并且結(jié)合三維斷裂力學的形函數(shù)考慮裂紋尖端的奇異性，探究將改進的無單元Galerkin（improved element-free Galerkin， IEFG）法用于斷裂彈性動力學問題的數(shù)值計算。通過懸臂梁、柱和矩形板等3個算例，討論節(jié)點分布、影響域比例參數(shù)、罰因子和時間步長等參數(shù)對計算精度的影響，證明IEFG法用于求解三維斷裂彈性動力學問題的正確性和有效性。

關鍵詞：

彈性動力學; 斷裂力學; 裂紋; 改進的移動最小二乘法; 形函數(shù); 隱式時間積分

中圖分類號：O241.82; TB115.1

文獻標志碼：B

Improved element-free Galerkin method for

three-dimensional fracture elasto-dynamics

JING Yongqiang， PENG Miaojuan

（Department of Civil Engineering， School of Mechanics and Engineering Science， Shanghai University， Shanghai 200444， China）

Abstract：

To improve the accuracy of the numerical calculation of fracture elasto-dynamics， and to avoid the ill conditioned or singular equations， the improved element-free Galerkin（IEFG） method is applied to the numerical calculation of fracture elasto-dynamics. Based on the improved moving least square method， the integral weak form of the three-dimensional elasto-dynamics problem is established. The penalty function method is used to set the displacement boundary conditions. The singularity of crack tip is considered by introducing the implicit time integration and the shape function of three-dimensional fracture mechanics. Three examples（include the cantilever beam， the column and the square plate） are selected to analyze the influence of computing parameters（include the node distribution， the influence domain scale parameter， the penalty factor and the time step） on the calculation accuracy. It is proved that the IEFG method is correct and effective for solving three-dimensional fracture elasto-dynamics problems.

Key words：

elasto-dynamics; fracture mechanics; crack; improved moving least square method; shape function; implicit time integration

0 引 言

斷裂彈性動力學問題是斷裂力學與彈性動力學耦合的問題，傳統(tǒng)的解析方法難以求解，目前可采用的數(shù)值方法主要有有限元法和無單元法等。采用有限元法求解斷裂力學問題時，在裂紋的發(fā)生擴展過程中需要重新劃分網(wǎng)格，計算精度較低。無單元法無須劃分網(wǎng)格，只需要節(jié)點信息。對于含有裂紋的構(gòu)件，為提高計算精度，無單元法可以在裂紋處布置加密節(jié)點，避免對網(wǎng)格和單元依賴的問題，在求解斷裂彈性動力學問題時具有一定的優(yōu)勢。

BELYTSCHKO等[1]采用無單元Galerkin法研究三維幾何和材料非線性的動力學問題。CHEN等[2]采用無單元Galerkin法研究材料高速撞擊發(fā)生彈塑性大變形的問題。LEE等[3]采用無單元Galerkin法對裂紋的奇點進行不連續(xù)性建模，提出一種改進的裂紋分析技術。RAO等[4]提出一種改進的無單元法，用于均勻、等向、非線性彈性二維固體裂縫的斷裂分析。

無單元Galerkin法利用移動最小二乘法構(gòu)造試函數(shù)，因此具有較高的計算精度，但計算過程中形成的方程組有時是病態(tài)的或奇異的。陳美娟等[5]提出改進的移動最小二乘法，基函數(shù)采用正交函數(shù)，可有效避免病態(tài)或奇異方程組的產(chǎn)生，不需要求解逆矩陣，可提高計算效率。ZHANG等[6]提出改進的無單元Galerkin（improved element-free Galerkin， IEFG）法，并用于解決三維勢問題、波動方程、斷裂力學和瞬態(tài)熱傳導問題。CHENG等[7]提出針對三維彈性力學、彈性動力學等問題的IEFG法。SUN等[8]研究共平面方形裂紋，對等向彈性體中二維和三維共平面方形裂紋的相互作用進行數(shù)值模擬。PENG等[9]在改進的移動最小二乘法的基礎上，提出針對三維黏彈性問題的IEFG法。ZHANG等[10]基于改進的移動最小二乘法，在裂縫尖端使用基函數(shù)，提出一種可以解決二維斷裂問題的IEFG法。魏慶節(jié)等[11]基于改進的移動最小二乘法建立三維彈性動力學問題的形函數(shù)，結(jié)合三維彈性動力學的Galerkin積分弱形式，采用罰函數(shù)法施加位移邊界條件，引入隱式時間積分，并提出三維彈性動力學的IEFG法。蔡小杰等[12]基于改進的移動最小二乘法，提出針對彈塑性大變形問題的IEFG法，討論權(quán)函數(shù)、影響域比例參數(shù)、罰因子、節(jié)點數(shù)和迭代步數(shù)等參數(shù)對計算精度的影響。鄒詩瑩等[13]在IEFG法的基礎上，將能反映裂紋尖端附近應力奇異性的特征項引入改進的移動最小二乘法的基函數(shù)中，將斷裂力學與IEFG法相結(jié)合，用于研究線彈性斷裂力學問題，并對含反射裂縫的機場復合路面層狀體系結(jié)構(gòu)進行數(shù)值分析。

劉學聰?shù)萚14]基于擴展有限元法提出一種新的裂紋尖端加強函數(shù)，組合傳統(tǒng)的函數(shù)基，繼承傳統(tǒng)附加函數(shù)的特性，使節(jié)點的奇異附加自由度減少為2個，可減小總矩陣的規(guī)模，提高計算效率。石路楊等[15]引入裂紋交叉匯合加強函數(shù)，分析多裂紋的交叉匯合過程，在裂紋附近區(qū)域使用廣義形函數(shù)，并引入線性遞增函數(shù)消除混合單元，建立求解多裂紋擴展的擴展有限元法，可有效提高裂紋附近的分析精度。楊永濤等[16]基于包含裂紋尖端增強函數(shù)的數(shù)值流形方法，采用Newmark隱式動力學算法進行時間積分，重點研究動力載荷作用下裂紋動態(tài)應力強度因子的求解方法。張洪武等[17]基于精細網(wǎng)格的傳統(tǒng)有限元分析方法，針對非均質(zhì)飽和多孔介質(zhì)彈塑性動力問題，提出一種廣義耦合擴展多尺度有限元法。

分析以上研究發(fā)現(xiàn)，在求解斷裂力學問題時IEFG法優(yōu)于有限元法。目前，IEFG法在三維彈性力學、彈性動力學、彈塑性大變形、黏彈性和線彈性斷裂力學等方面的研究較多，針對斷裂彈性動力學問題，IEFG法的研究很少。本文基于改進的移動最小二乘法，建立三維彈性動力學問題的積分弱形式，采用罰函數(shù)法施加位移邊界條件，引入隱式時間積分，并且結(jié)合三維斷裂力學的形函數(shù)考慮裂紋尖端的奇異性。采用改進的移動最小二乘法可避免病態(tài)或奇異方程組，保證IEFG法的計算精度。三維斷裂彈性動力學IEFG法不僅要考慮三維彈性動力學中加速度和阻尼的影響，還要考慮斷裂力學裂紋尖端的奇異性。本文結(jié)合改進的移動最小二乘法、彈性動力學和斷裂力學的無單元Galerkin法，建立三維斷裂彈性動力學的IEFG法的離散方程，并推導相應的計算公式，采用懸臂梁、柱和方形板等3個數(shù)值算例證明三維斷裂彈性動力學IEFG法的正確性和有效性。

1 三維斷裂彈性動力學基本方程

三維線彈性動力學的控制方程為

LTσ+b=ρu··+μu·， x∈Ω

（1）

式中：L為微分算子矩陣;σ為應力向量;b為外力向量;ρ為質(zhì)量密度;u為位移向量;μ為阻尼系數(shù);x為點的集合。L、σ、b和u分別定義為

L=x1000x2000x3x2x10x30x10x3x2

（2）

σ=[σ11 σ22 σ33 σ12 σ13 σ23]T

（3）

b=[b1（t） b2（t） b3（t）]T

（4）

u=[u1（t） u2（t） u3（t）]T

（5）

式中：xi（i=1，2，3）為點x緊支域內(nèi)的點;σij（i，j=1，2，3）為不同方向的應力;t為時間。

三維線彈性動力學控制方程的邊界條件為

ui-i=0， i=1，2，3， x∈Γu

（6）

njσij-ti=0， i，j=1，2，3， x∈Γt

（7）

式中：i為邊界Γu上的已知位移;nj為邊界Γt的外法線方向的余弦;ti為邊界Γt上的已知面力。Γu∪Γt=Γ，Γ為求解域Ω的邊界。

三維線彈性動力學控制方程的初始條件為

ui（x，t0）=ui0（x）， x∈Ω（8）

u·i（x，t0）=vi0（x）， x∈Ω （9）

式中：vi0和ui0分別為點x的初始速度和初始位移;t0為初始時間。

三維線彈性動力學控制方程的本構(gòu)方程為

σ=Dε

（10）

其中：

ε=Lu

（11）

D=E（1-v）（1+v）（1-2v）1v1-vv1-v000v1-v0v1-v000v1-vv1-v00000001-2v2（1-v）0000001-2v2（1-v）0000001-2v2（1-v）

（12）

式（11）即為三維線彈性動力學控制方程的幾何方程。

為得到彈性動力學的Galerkin弱形式，采用罰函數(shù)法施加位移邊界條件，即

∫ΩδuT·ρu··dΩ+ ∫ΩδuT·μu·dΩ+∫ΩδεT·σdΩ+

∫ΓtδuT·bdΩ-∫ΓtδuT·tdΩ+

α∫ΓtδuT·S（u-）dΓ=0

（13）

式中：α為罰因子;

=[1 2 3]T （14）

t=[t1 t2 t3]T（15）

將式（10）和（11）代入式（13），計算得到罰函數(shù)法的可行域為

S=S1000S2000S3

（16）

∫ΩδuT·ρu··dΩ+∫ΩδuT·μu·dΩ+∫Ωδ（Lu）T·（Lu）dΩ+∫ΩδuT·bdΩ-∫ΓtδuT·tdΩ+α∫ΓtδuT·S（u-）dΓ=0（17）

當x1、x2、x3有位移約束時，對應的S1、S2、S3等于1，否則等于0。

2 三維斷裂彈性動力學的IEFG法

利用改進的移動最小二乘法建立逼近函數(shù)。由文獻[13]中裂紋尖端的應力場表達式可知，裂紋尖端附近的位移與特征項r（r為點x到裂紋尖端的距離）成正比，應力具有1/r階奇異性。為有效計算裂紋尖端場的應力數(shù)值解，在基函數(shù)中加入擴展函數(shù)項，由改進的移動最小二乘法得到形函數(shù)。基于三維彈性動力學IEFG法，在基函數(shù)中加入文獻[13]中裂紋尖端位移場表達式中的重要項和重要的梯度項，使IEFG法的計算結(jié)果更精確。

結(jié)合文獻[13]得到三維斷裂力學IEFG法的試函數(shù)為

pT（x）=1 x1 x2 x3 rcos θ2 rsin θ2

rsin θ2sin θ rcos θ2sin θ

（18）

PI（x）=[w（x-x1）p（x1） w（x-x2）p（x2） …

w（x-xn）p（xn）]

（19）

式中：θ為裂紋的擴展角;n為影響域覆蓋點x的節(jié)點數(shù);xn為點x的緊支域內(nèi)的節(jié)點;w（x-xn）為具有緊支集特性的權(quán)函數(shù)。

類似文獻[11]中移動最小二乘法的形函數(shù)推導，可得到三維斷裂彈性動力學的逼近函數(shù)為

uα，h（x）=nI=1pT（x）A-1（x）PI（x）uα（xI）

（20）

式中：A-1（x）為由pi構(gòu)成的對角矩陣。

三維斷裂彈性動力學的形函數(shù)為

ΦI（x）=pT（x）A-1（x）PI（x）

（21）

將式（21）代入式（20），可得

uα，h（x）=nI=1ΦI（x）uα（xI）

（22）

采用局部基函數(shù)擴展法，在基函數(shù)中加入特征項r，即

pT（x）=[1 x1 x2 x3 r]

（23）

將求解域離散為有限個節(jié)點，節(jié)點總數(shù)為M。利用改進的移動最小二乘法建立逼近函數(shù)，可以得到位移的逼近函數(shù)為

ui（x，t）=nI=1ΦI（x）uI（xI，t）， i=1，2，3

（24）

即 u（x，t）=nI=1ΦI（x）uI

（25）

將式（24）和（25）代入式（17）可得

∫ΩδnI=1ΦIuITρnJ=1ΦJu··JdΩ+∫ΩδnI=1ΦIuITμnJ=1ΦJu·JdΩ+∫ΩδnI=1LuITDnJ=1ΦJLuJdΩ-

∫ΩδnI=1ΦIuITbdΩ-∫ΓtδnI=1ΦIuITtdΓ+∫ΓtδnI=1ΦIuITSnJ=1ΦJu··JdΓ+∫ΓtδnI=1ΦIuITSu··dΓ=0

（26）

對式（26）進行積分運算，可以得到三維斷裂彈性動力學的離散方程為

MU··（t）+CU·（t）+（K+Kα）U（t）=F（t）+Fα

（27）

其中：

M=（Mij），Mij=∫ΩΦiρΦjdΩ， i，j=1，2，…，M

（28）

U=[uT1（t） uT2（t） … uTM（t）]T

（29）

C=（Cij），Cij=∫ΩΦiμΦjdΩ， i，j=1，2，…，M

（30）

K=（Kij），Kij=∫ΩBTiDBjdΩ， i，j=1，2，…，M

（31）

Bi=Φi，2000Φi，2000Φi，3Φi，2Φi，10Φi，30Φi，10Φi，3Φi，2， i=1，2，…，M

（32）

Kα=（Kij，α），Kij，α=α∫ΓtΦiSΦjdΓ，

i，j=1，2，…，M

（33）

F=（Fij），F(xiàn)ij=∫ΩΦibdΩ+∫ΓtΦjtdΓ，

i，j=1，2，…，M

（34）

Fα=（Fi，α），F(xiàn)i，α=α∫ΓtΦiSdΓ，i=1，2，…，M

（35）

3 隱式時間積分

隱式時間積分方法可實現(xiàn)求解目標在大時間步長內(nèi)的無條件穩(wěn)定，被廣泛應用于各種動力學問題中。本文采用Newmark-β法對運動方程進行時間離散。假設時刻t的時間域[0，T]被離散為n個時間步（即Δt=T/n），位移Ut及其對應導數(shù)U·t、U··t的值均已知，對U·t+Δt和U··t+Δt進行二階Taylor展開可得

Ut+Δt=Ut+ΔtU·t+12（1-β2）Δt2U··t+12β2Δt2U··t+Δt

（36）

U·t+Δt=U·t+（1-β1）ΔtU··t+β1ΔtU··t+Δt

（37）

U··t+Δt=2β2Δt2（Ut+Δt-Ut）-2β2ΔtU·t-1β2-1U··t

（38）

式中：β1和β2均為Newmark參數(shù)，分別取β1=3/2，β2=8/5。

令α1=2β2Δt2、α2=2β2Δt、α3=1β2-1，可得

U··t+Δt=α1（Ut+Δt-Ut）-α2U·t-α3U··t

（39）

將式（39）代入式（37）可得

U·t+Δt=β1α2（Ut+Δt-Ut）+1-2β1β2U·t+

1-β1β2ΔtU··t

（40）

根據(jù)隱式時間積分方法，式（27）可寫成

MU··t+Δt+CU·t+Δt+（K+Kα）Ut+Δt=Ft+Δt+Fα

（41）

將式（39）和（40）代入式（41）可得

（α1M+β1α2C+K+Kα）Ut+Δt=

Ft+Δt+Fα+M（α1Ut+α2U·t+α3U··t）+

Cβ1α1Ut+2β1β2-1U·t+β1β2-1ΔtU··t

（42）

4 數(shù)值算例

選擇懸臂梁、柱和矩形板等3個算例，驗證IEFG法的正確性和有效性。3個算例的權(quán)函數(shù)均為三次樣條權(quán)函數(shù)，將計算結(jié)果與Abaqus軟件有限元法計算結(jié)果進行比較。

為分析IEFG法的精確度，定義均方差為

eI=1NNi=1（uIEFG，i-uAbaqus，i）

（43）

式中：uIEFG，i為IEFG法的計算結(jié)果;uAbaqus，i為Abaqus的計算結(jié)果;N為節(jié)點數(shù)。

為討論Abaqus解的準確性和收斂性，定義方差

eA=1NNk=1（uJ，k-IJ，k）2+（uI，k-IJ，k）2

（44）

式中：uI，k和uJ，k分別為節(jié)點分布I和J下第k個節(jié)點的數(shù)值解;IJ，k為uI，k和uJ，k的平均值。

4.1 受突加線性分布載荷作用的懸臂梁

選取某典型懸臂梁，上、下表面自由，一端固定，另一端受突加線性分布載荷作用，見圖1，其中a=120 mm、b=120 mm、c=360 mm、r=30 mm。材料參數(shù)分別取ρ=2 405 kg/m3、E=20.69 GPa、ν=0.15，g取9.8 m/s2。在x3=-60 mm和x3=60 mm處施加突加載荷P=1×103 N。

利用有限元軟件Abaqus對該算例進行計算，有限元網(wǎng)格采用20節(jié)點二次六面體單元C3D20R，分別取9×9×15、15×15×31、15×15×51、17×17×65和33×33×81等5種節(jié)點分布，懸臂梁自由端Q點（x1=360 mm、x2=-60 mm）位移數(shù)值解的方差見圖2。由此可知，隨著網(wǎng)格的細化，方差逐漸趨于0，說明Abaqus的數(shù)值解收斂。與IEFG方法的數(shù)值解進行對比時選擇節(jié)點分布為33×33×81，Abaqus中懸臂梁的網(wǎng)格劃分結(jié)果見圖3。

根據(jù)三維斷裂彈性動力學IEFG法的離散方程，采用與Abaqus計算時同樣的節(jié)點分布，利用MATLAB對懸臂梁Q點位移進行數(shù)值計算，不同節(jié)點分布時懸臂梁位移的均方差結(jié)果見圖4。由此可知，隨著節(jié)點數(shù)增加，IEFG法的均方差逐漸減小，計算精度逐漸提高。

影響域比例參數(shù)dmax在1.5～3.5范圍內(nèi)取不同值時，懸臂梁位移均方差的變化見圖5。由此可知，當dmax在2.5～3.0范圍內(nèi)時，懸臂梁位移的均方差較小，說明此時數(shù)值解的計算精度較高。

給定影響域比例參數(shù)dmax，取不同時間步長Δt，計算懸臂梁位移的均方差，結(jié)果見圖6。由此可知，隨著時間步長的增大，均方差逐漸增大。當時間步長為1.0×10-6 s時均方差較小，說明此時數(shù)值解的計算精度較高。

給定影響域比例參數(shù)dmax和時間步長，取不同罰因子α，計算懸臂梁位移的均方差，結(jié)果見圖7。由此可知，隨著罰因子的增大，均方差呈類似拋物線分布。當罰因子α=2.5×1015時均方差較小，說明此時數(shù)值解的計算精度較高。

取dmax=2.7、Δt=1.0×10-6 s、α=2.5×1015，Abaqus中節(jié)點分布取33×33×81，計算圖1中自由端Q點的位移u1、u3和應力σ隨時間的變化，見圖8～10。由此可知，IEFG法與有限元法的計算結(jié)果吻合，證明IEFG法解決三維斷裂彈性動力問題的正確性和有效性。

4.2 受突加載荷作用的柱

選取單邊帶有半圓型裂紋的柱形試件，裂紋位于試件中部，裂紋半徑r=5 cm，試件下端固定，上端施加單向突加載荷作用，見圖11。突加載荷P=20 kPa，材料和幾何參數(shù)取E=210 GPa、ν=0.35、a=30 cm、b=30 cm、c=120 cm。

利用Abaqus對柱位移進行計算，有限元網(wǎng)格采用20節(jié)點二次六面體單元C3D20R，取9×9×17、9×9×31、17×17×64、21×21×71和41×41×81等5種節(jié)點分布，計算得到Q點（x1=-15 cm、x3=120 cm）位移的方差，見圖12。隨著節(jié)點數(shù)增加，Abaqus數(shù)值解的方差趨于0，說明Abaqus的數(shù)值解收斂。與IEFG法的數(shù)值解進行對比時選擇節(jié)點分布為41×41×81，Abaqus中柱的網(wǎng)格劃分結(jié)果見圖13。

根據(jù)三維斷裂彈性動力學IEFG法的離散方程，取dmax=2.8、Δt=1.0×10-6 s、α=2.5×1015，利用MATLAB對柱位移進行數(shù)值計算，圖11中Q點的位移u1、u3和應力σ隨時間的變化見圖14～16。由此可知，IEFG法與有限元法計算結(jié)果吻合，證明IEFG法的正確性和有效性。

4.3 受突加載荷作用的矩形板

選取帶有半圓型裂紋的矩形板，下表面固定，上表面受突加均布載荷作用，見圖17，其中r=0.03 m、a=3 m、b=3 m、c=1 m。矩形板的彈性模量E=210 GPa、泊松比ν=0.15，突加載荷P=100 MPa。

利用Abaqus計算矩形板位移，有限元網(wǎng)格采用20節(jié)點二次六面體單元C3D20R，取21×21×11，31×31×21，41×41×31，61×61×41和81×81×61等5種節(jié)點分布，圖17中矩形板Q點（x1=1.5 m、x3=1.0 m）位移的方差見圖18。隨著節(jié)點數(shù)增加，方差趨于0，說明Abaqus的數(shù)值解收斂。與IEFG法數(shù)值解對比時選擇節(jié)點分布為81×81×61，Abaqus中矩形板的網(wǎng)格劃分結(jié)果見圖19。

根據(jù)三維斷裂彈性動力學IEFG法的離散方程，利用MATLAB對矩形板位移進行數(shù)值計算，取dmax=2.8、Δt=1.0×10-6 s、α=2.5×1015，圖17中Q點的位移u1、u3和應力σ隨時間的變化見圖20～22。IEFG法與有限元法計算結(jié)果吻合，證明IEFG法的正確性和有效性。

5 結(jié)束語

基于改進的移動最小二乘法，提出三維斷裂彈性動力學的IEFG法。通過懸臂梁、柱和矩形板3個算例，對比分析節(jié)點分布、影響域比例參數(shù)、罰因子和時間步長等參數(shù)對計算精度的影響。本文提出的三維斷裂力學的IEFG法計算結(jié)果與有限元法計算結(jié)果吻合，證明該方法的正確性和有效性。

參考文獻：

[1] BELYTSCHKO T， KRYSL P， KRONGAUZ Y. A three-dimensional explicit element-free Galerkin method[J]. International Journal for Numerical Methods in Fluids， 1997， 24（12）： 1253-1270. DOI： 10.1002/（SICI）1097-0363（199706）24：12〈1253：：AID-FLD558〉3.0.CO;2-Z.

[2] CHEN Y P， ESKANDARIAN A， OSKARD M， et al. Meshless analysis of high-speed impact[J]. Theoretical and Applied Fracture Mechanics， 2005， 44（3）： 201-207. DOI： 10.1016/j.tafmec.2005.09.007.

[3] LEE S H， YOON Y C. An improved crack analysis technique by element-free Galerkin method with auxiliary supports[J]. International Journal for Numerical Methods in Engineering， 2003， 56（9）： 1291-1314. DOI： 10.1002/nme.611.

[4] RAO B N， RAHMAN S. An enriched meshless method for non-linear fracture mechanics[J]. International Journal for Numerical Methods in Engineering， 2004， 59（2）： 197-223. DOI： 10.1002/nme.868.

[5] 陳美娟， 程玉民. 改進的移動最小二乘法[J]. 力學季刊， 2003， 24（2）： 266-272. DOI： 10.3969/j.issn.0254-0053.2003.02.019.

[6] ZHANG Z， WANG J F， CHENG Y M， et al. Improved element-free Galerkin method for three-dimensional transient heat conduction problems[J]. Science China Physics， Mechanics and Astronomy， 2013， 56： 1568-1580. DOI： 10.1007/s11433-013-5135-0.

[7] CHENG Y M， PENG M J. Boundary element-free method for elastodynamics[J]. Science in China Series G： Physics and Astronomy， 2005， 48： 641-657. DOI： 10.1360/142004-25.

[8] SUN Y Z， LIEW K M. Analyzing interaction between coplanar square cracks using an efficient boundary element-free method[J]. International Journal for Numerical Methods in Engineering， 2012， 91（11）： 1184-1198. DOI： 10.1002/nme.4309.

[9] PENG M J， LI R X， CHENG Y M. Analyzing three-dimensional viscoelasticity problems via improved element-free Galerkin（IEFG） method[J]. Engineering Analysis with Boundary Elements， 2014， 40： 104-113. DOI： 10.1016/j.enganabound.2013.11.018.

[10] ZHANG Z， LIEW K M， CHENG Y M， et al. Analyzing 2D fracture problems with improved element-free Galerkin method[J]. Engineering Analysis with Boundary Elements， 2008， 32（3）： 241-250. DOI： 10.1016/j.enganabound.2007.08.012.

[11] 魏慶節(jié)， 彭妙娟. 三維彈性動力學的改進的無單元Galerkin方法[J]. 應用力學學報， 2019， 36（1）： 22-27.

[12] 蔡小杰， 彭妙娟， 程玉民. 彈塑性大變形問題的改進的無單元Galerkin方法[J]. 中國科學， 2018， 48（2）： 53-62. DOI： 10.1360/SSPMA2017-00231.

[13] 鄒詩瑩， 席偉成， 彭妙娟， 等. 運用改進的無單元Galerkin方法分析機場道面斷裂力學問題[J]. 物理學報， 2017， 66（12）： 36-44. DOI： 10.7498/aps.66.120204.

[14] 劉學聰， 章青， 夏曉舟. 一種新型裂尖加強函數(shù)的顯式動態(tài)擴展有限元法[J]. 工程力學， 2017， 34（10）： 10-18. DOI： 10.6052/j.issn.1000-4750.2016.05.0399.

[15] 石路楊， 余天堂. 多裂紋擴展的擴展有限元法分析[J]. 巖土力學， 2014， 35（1）： 263-272. DOI： 10.16285/j.rsm.2014.01.001.

[16] 楊永濤， 徐棟棟， 鄭宏. 動載下裂紋應力強度因子計算的數(shù)值流形元法[J]. 力學學報， 2014， 46（5）： 730-738. DOI： 10.6052/0459-1879-14-024.

[17] 張洪武， 盧夢凱， 鄭勇剛. 非均質(zhì)飽和多孔介質(zhì)彈塑性動力分析的廣義耦合擴展多尺度有限元法[J].計算力學學報， 2016， 33（4）： 454-461. DOI： 10.7511/jslx201604005.

（編輯 武曉英）