改進(jìn)的鍵基正交各向異性近場動(dòng)力學(xué)模型

鄭國君 陳瑞 申國哲 夏陽

摘要：針對經(jīng)典的鍵基近場動(dòng)力學(xué)（bond?based peridynamic， BPD）模型受固定泊松比限制的問題，提出一種改進(jìn)BPD模型。該模型可解除泊松比限制，并可用于分析正交各向異性單向板的變形和裂紋擴(kuò)展問題。在改進(jìn)BPD模型中，每根鍵受到軸向和橫向成對力的作用，額外增加的節(jié)點(diǎn)轉(zhuǎn)動(dòng)可消除由橫向力引起的附加彎矩，從而確保該模型滿足角動(dòng)量守恒條件。仿真結(jié)果驗(yàn)證所提出的改進(jìn)BPD模型的精度，并展示其預(yù)測碳纖維復(fù)合材料變形和裂紋擴(kuò)展的能力。

關(guān)鍵詞：

近場動(dòng)力學(xué); 鍵基; 剪切影響系數(shù); 泊松比; 顆粒旋轉(zhuǎn); 角動(dòng)量守恒; 正交各向異性; 碳纖維復(fù)合材料

中圖分類號(hào)：TB334;TP391.99

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

An improved bond?based peridynamic orthotropic model

ZHENG Guojun， CHEN Rui， SHEN Guozhe， XIA Yang

（School of Vehicle Engineering， Dalian University of Technology， Dalian 116024， Liaoning， China）

Abstract：

As to the problem that the classical bond?based peridynamic（BPD） model is limited by the fixed Poisson′s ratio， an improved BPD model is proposed， which removes the limitation of Poisson′s ratio and can be used to analyze the deformation and crack propagation of orthotropic anisotropic plates. In the improved BPD model， each bond is subjected to axial and transverse forces in pairs， and the additional rotation of particles can eliminate the additional bending momentum caused by transverse forces， so that the balance of angular momentum can be ensured. The simulation results can verify the accuracy of the improved BPD model， and show its ability to predict the deformation and crack propagation in carbon fiber composite material.

Key words：

peridynamic; bond?based; shear influence coefficient; Poisson′s ratio; particle rotation; angular momentum balance; orthotropic; carbon fiber composite material

引?言

復(fù)合材料是指由2種或者2種以上不同性質(zhì)的材料通過物理或者化學(xué)的方法結(jié)合而成的具有新性能的材料。復(fù)合材料具有質(zhì)量輕、強(qiáng)度高、抗疲勞性好等優(yōu)點(diǎn)[1]，在航空航天、汽車、風(fēng)力發(fā)電、醫(yī)療器械等領(lǐng)域的應(yīng)用越來越廣泛，如空中客車A380、波音787等客機(jī)的主要結(jié)構(gòu)（整體機(jī)身、機(jī)翼和艙門等）都采用復(fù)合材料[2]。隨著復(fù)合材料應(yīng)用領(lǐng)域越來越廣泛，其載荷環(huán)境日益復(fù)雜，損傷失效問題日益突出，對其失效模式的分析也成為工程中的重點(diǎn)和難點(diǎn)。

傳統(tǒng)的連續(xù)介質(zhì)理論具有微分形式的運(yùn)動(dòng)方程，在面對損傷破壞等不連續(xù)問題時(shí)，由于位移場的不連續(xù)性，在不連續(xù)的區(qū)域難以得到位移場的偏微分方程，因此傳統(tǒng)連續(xù)介質(zhì)理論在求解不連續(xù)問題時(shí)會(huì)遇到困難[3?4]。目前，針對復(fù)合材料結(jié)構(gòu)損傷等問題主要采用有限元法，但是有限元法也是以連續(xù)介質(zhì)理論為基礎(chǔ)的，因此使用有限元法求解損傷失效問題也極具挑戰(zhàn)性，需要借助一些附加的失效準(zhǔn)則重新劃分網(wǎng)格并進(jìn)行求解[5?7]，而且損傷只能沿某些特定的方向傳播。為解決傳統(tǒng)有限元在求解損傷裂紋時(shí)的缺陷，有學(xué)者提出擴(kuò)展有限元和非連續(xù)有限元等思想[8?14]，該思想主要使用獨(dú)立于網(wǎng)格劃分的思想解決裂紋擴(kuò)展問題，不需要對結(jié)構(gòu)內(nèi)部存在的裂紋等缺陷進(jìn)行網(wǎng)格劃分。擴(kuò)展有限元法可提高模型描述復(fù)雜位移場的能力，避免網(wǎng)格的重新劃分，解決大量的斷裂問題。但是，在解決不連續(xù)問題時(shí)，擴(kuò)展有限元法仍然需要額外的斷裂準(zhǔn)則。COX等[15]提出分子動(dòng)力學(xué)理論解決傳統(tǒng)有限元存在的缺陷，但卻增大計(jì)算的消耗。KRNER等[16]根據(jù)原子間的作用力，考慮長程力的影響，提出非局部連續(xù)理論，該理論可不區(qū)分不連續(xù)性，但由于裂紋的存在，其位移導(dǎo)數(shù)仍然是不存在的。

為解決傳統(tǒng)數(shù)值方法在求解不連續(xù)問題時(shí)存在的困難，SILLING[17]和SILLING等[18]提出近場動(dòng)力學(xué)理論（peridynamic，PD），基于非局部思想，通過求解空間積分方程描述物質(zhì)的力學(xué)行為。

鍵基近場動(dòng)力學(xué)（bond?based PD，BPD）模型使用成對力函數(shù)描述節(jié)點(diǎn)間的相互作用，因此對于各向同性材料來說，鍵力僅與兩個(gè)節(jié)點(diǎn)間的相對伸長率有關(guān)，在平面應(yīng)力狀態(tài)下泊松比被限制為1/3，在平面應(yīng)變狀態(tài)下泊松比被限制為1/4[19]。GERSTLE等[20]通過引入桿單元提出微極模型，通過增加成對力矩模擬具有不同泊松比的線彈性材料，這樣可解除泊松比的限制，但是未考慮鍵的轉(zhuǎn)動(dòng)，因此不滿足角動(dòng)量守恒。REN等[21]提出考慮剪切變形的線彈性固體PD模型，從總變形中減去剛體的轉(zhuǎn)動(dòng)部分，但是未具體處理泊松比的限制。GHAJARI等[22]引入勒讓德多項(xiàng)式，提出正交各向異性材料的連續(xù)微模量函數(shù)模型，能夠預(yù)測復(fù)雜的斷裂現(xiàn)象，但是其泊松比仍然是受限制的。ZHU等[23]考慮鍵轉(zhuǎn)動(dòng)的影響，重構(gòu)近場動(dòng)力學(xué)鍵基模型，消除泊松比的限制。ZHOU等[24]提出共軛鍵線彈性模型，鍵能不僅與其法向的伸長有關(guān)，也與一對共軛鍵的旋轉(zhuǎn)角度有關(guān)，從而克服泊松比的限制。SILLING等[25]提出基于狀態(tài)的PD模型，通過引入節(jié)點(diǎn)的變形狀態(tài)，徹底解除泊松比的限制，但其計(jì)算過程比較復(fù)雜。

對于正交各向異性材料，雖然有些PD模型可解除泊松比情況限制，但未考慮鍵轉(zhuǎn)動(dòng)的影響，不滿足角動(dòng)量守恒，在剪切模擬時(shí)誤差較大。基于狀態(tài)的PD模型雖然可有效解除泊松比的限制，但是基于狀態(tài)的PD模型計(jì)算過程比基于鍵的PD模型更復(fù)雜，應(yīng)用更不便。

本文基于Timoshenko梁理論提出一種改進(jìn)的鍵基梁模型，考慮節(jié)點(diǎn)的旋轉(zhuǎn)效應(yīng)和剪切變形，不僅可以解除泊松比的限制，也可以使橫向剛度計(jì)算更加精確。

1?BPD基礎(chǔ)理論

在經(jīng)典的BPD模型中，PD模型的參考構(gòu)型和當(dāng)前構(gòu)型[17]見圖1。參考構(gòu)型中某節(jié)點(diǎn)在位置xi的動(dòng)力學(xué)方程可以寫為

式中：ρ為質(zhì)量密度;u為位移矢量場;f為成對力函數(shù)，表示節(jié)點(diǎn)xj施加到節(jié)點(diǎn)xi上的每體積的力;V為節(jié)點(diǎn)xi所占據(jù)的體積;b為體力矢量場;Ηxi為節(jié)點(diǎn)xi的作用域。

為滿足線動(dòng)量和角動(dòng)量守恒，對參考構(gòu)型中的任意相對位置矢量ξ和相對位移矢量η有

f（η，ξ）=-f（-η，-ξ），η，ξ

（η+ξ）×f（η，ξ）=0，η，ξ

（2）

其中：

ξ=xj-xi

η=u（xj，t）-u（xi，t）

（3）

根據(jù)式（2）和（3）可知，兩個(gè)節(jié)點(diǎn)彼此施加的力大小相等、方向相反，且與當(dāng)前構(gòu)型中的相對位置矢量平行。

在微彈性材料模型中，成對力函數(shù)可以由標(biāo)量微勢能函數(shù)w（η，ξ）推導(dǎo)[17]，即

f（η，ξ，t）=w（η，ξ）η

（4）

單個(gè)鍵力與相對伸長率之間的線性關(guān)系可以根據(jù)線微彈性理論假設(shè)和微勢能函數(shù)推導(dǎo)，即

w（η，ξ）=c（η，ξ）s2ξ/2（5）

式中：c（η，ξ）為鍵剛度常數(shù);s為鍵的伸長;ξ為ξ的模。s的定義為

s=（ξ+η-ξ）/ξ

（6）

對于微彈性材料來說，給定節(jié)點(diǎn)xi處的應(yīng)變能密度可通過對作用域內(nèi)的微勢能函數(shù)w（η，ξ）積分得到，即

WPD（xi）=12∫Ηxiw（η，ξ）dVxj

（7）

式中的1/2是因?yàn)槊總€(gè)節(jié)點(diǎn)的鍵能只有總鍵能的1/2。相同載荷作用下PD得到的應(yīng)變能密度應(yīng)等于經(jīng)典彈性力學(xué)理論得到的應(yīng)變能密度，即

WPD（xi）=WCL（xi）

（8）

從而可以得到鍵剛度系數(shù)與彈性模量之間的關(guān)系。

2?改進(jìn)的BPD模型

2.1?求解鍵常數(shù)

對于改進(jìn)的BPD模型，由節(jié)點(diǎn)xi和在其近場區(qū)域Ηxi內(nèi)的相鄰點(diǎn)xj組成PD鍵，見圖2。

與經(jīng)典的BPD模型相同，假設(shè)鍵的軸向力只與鍵的拉伸變形和軸向力密度有關(guān)，即

f^xij=-f^xji=cNs

（9）

式中：f^xij為在近場區(qū)域Ηxi內(nèi)的節(jié)點(diǎn)xj施加在節(jié)點(diǎn)xi上的軸向力密度;f^xji為在近場區(qū)域Ηxj內(nèi)的節(jié)點(diǎn)xi施加在節(jié)點(diǎn)xj上的軸向力密度;cN為要求解的PD參數(shù)。

考慮用鍵的橫向力消除泊松比的限制。鍵的橫向力附著在兩個(gè)節(jié)點(diǎn)上，會(huì)出現(xiàn)附加的彎矩，必須對其進(jìn)行消除以滿足角動(dòng)量守恒。類似于有限元模型中的梁單元，可通過增加兩個(gè)節(jié)點(diǎn)的旋轉(zhuǎn)自由度抵消由橫向力引起的彎矩。鍵的橫向位移、力和節(jié)點(diǎn)的旋轉(zhuǎn)角度及彎矩關(guān)系見圖2。與軸向力相同，在變形構(gòu)型上建立橫向力和彎矩。假設(shè)PD鍵是具有長度ζ和高度Δ的梁模型，高度相對于長度的比值Δ/ζ沒有小到足以忽略剪切變形的影響。Timoshenko梁理論可描述剪切變形和轉(zhuǎn)動(dòng)慣量對各種細(xì)長梁的影響。[26]剪切影響系數(shù)可解釋剪切應(yīng)力的變化，適于描述短梁的剪切變形性能。[27]在局部坐標(biāo)系中，基于Timoshenko梁理論的有限元方程為

p=Kd

（10）

式中：p為附著在材料節(jié)點(diǎn)xi和xj之間的鍵的橫向力和彎矩分量;d為鍵的位移分量;K為Timoshenko梁單元的局部剛度矩陣[28?29]。p、d和K分別定義為

p=[fixfiymizfjxfjymjz]T（11）

d=[uiviθiujvjθj]T

（12）

K=EAξ00-EAξ00012EI（1+b）ξ36EI（1+b）ξ20-12EI（1+b）ξ36EI（1+b）ξ206EI（1+b）ξ2（4+b）EI（1+b）ξ0-6EI（1+b）ξ2（2-b）EI（1+b）ξ-EAξ00EAξ000-12EI（1+b）ξ3-6EI（1+b）ξ2012EI（1+b）ξ3-6EI（1+b）ξ206EI（1+b）ξ2（2-b）EI（1+b）ξ0-6EI（1+b）ξ2（4+b）EI（1+b）ξ

（13）

式中：E為彈性模量;G為剪切模量;I為梁的轉(zhuǎn)動(dòng)慣量;b為剪切影響系數(shù)。b可使梁的橫向剛度更加準(zhǔn)確，其定義為

b=6EΔ2/5Gξ2（14）

引入節(jié)點(diǎn)的附加旋轉(zhuǎn)角度和彎矩以計(jì)算鍵的應(yīng)變能，基于Timoshenko梁理論的PD模型見圖3。

由鍵連接的兩個(gè)節(jié)點(diǎn)之間的位移是相互獨(dú)立的，初始的相對位置矢量為ξ，相對位移矢量為η，x′軸與全局坐標(biāo)系的x軸夾角為。

在單個(gè)鍵中，PD的軸向力密度、橫向力密度和彎矩密度可以表示為

f^y=12d（1+b）ξ2ξ-θ^1+θ^22

M^z=M^1z+M^2z2=6c（1+b）ξξ-θ^1+θ^22

（15）

式中：c和d為要求解的PD參數(shù);和分別為軸向和橫向位移;θ^1和θ^2為局部坐標(biāo)系中節(jié)點(diǎn)的旋轉(zhuǎn)角度;b為剪切影響系數(shù)。因此，

=2-1

=2-1

（16）

那么，存儲(chǔ)在由相鄰兩個(gè)節(jié)點(diǎn)xi和xj形成的PD鍵內(nèi)的應(yīng)變能密度為

w（η，ξ）=dTKd/2=

wf^x（η，ξ）+

wf^y（η，ξ）+wM^z（η，ξ）

（17）

式中：wf^x（η，ξ）為由軸向變形引起的應(yīng)變能密度;wf^y（η，ξ）為由剪切變形引起的應(yīng)變能密度;wM^z（η，ξ）為由節(jié)點(diǎn)旋轉(zhuǎn)力矩引起的應(yīng)變能密度。

wf^x（η，ξ）=c2ξ2

wf^y（η，ξ）=12d2（1+b）ξ22ξ-（θ^1+θ^2）2

wM^z（η，ξ）=d2（1+b）ξ-6ξ（θ^1+θ^2）+

（4+b）（θ^21+θ^22）+（4-2b）θ^1θ^2

（18）

節(jié)點(diǎn)的旋轉(zhuǎn)角度與坐標(biāo)系無關(guān)，因此定義

θ=（θ^1+θ^2）/2 （19）

那么應(yīng)變能密度方程可以化簡為

w（η，ξ）=

12cξ2+12d（1+b）ξξ-θ2+dξ（θ^1-θ^2）2

（20）

在局部小變形過程中，可以假設(shè)兩個(gè)節(jié)點(diǎn)的旋轉(zhuǎn)角度是相等的，即

θ^1=θ^2（21）

因此，式（20）可化簡為

w（η，ξ）=12cξ2+12d（1+b）ξξ-θ2

（22）

在x′和y′方向上的應(yīng)變分別為ε1和ε2，那么在局部坐標(biāo)系下的應(yīng)變定義為

ε1ε2γ12=cos2 sin2 sin cos sin2 cos2 -sin cos -2sin cos 2sin cos cos2 -sin2 εxεyγxy

（23）

由局部坐標(biāo)系與全局坐標(biāo)系關(guān)系可知

12=ξε1=ξ（εxcos2 θ+εysin2 θ+γxysin θcos θ）

12=ξ12γ12+θ=ξ-εxsin θcos θ+εysin θcos θ+γxy（cos2 θ-sin2 θ）2+θ

（24）

能量密度

WPD=12∫H12cξε21+12d（1+b）ξ12γ122tdH

（25）

正交各向異性單向板材料可以作為二維問題處理，其PD模型見圖4，厚度方向僅考慮單層物質(zhì)點(diǎn)，面內(nèi)分為纖維方向和基體方向，纖維鍵與基體鍵在拉伸方向的鍵剛度分別設(shè)為cf和cm，纖維鍵與基體鍵在剪切方向的鍵剛度分別設(shè)為df和dm，任意角度的鍵剛度滿足

c=cf+cm，θ=

c=cm，θ≠

d=df+dm，θ=

d=dm，θ≠

（26）

式中：θ為纖維方向與x軸的夾角。

因此，式（25）中的應(yīng)變能密度可以近似表示為

WPD=12∑Qq=112cfξε21+12df（1+b2）ξ12γ122tVq+

12∫H12cmξε21+12dm（1+b1）ξ12γ122tdH

（27）

式中：b1和b2分別為基體方向和纖維方向的剪切影響因子。

在經(jīng)典連續(xù)介質(zhì)力學(xué)中，正交各向異性材料的平面應(yīng)力狀態(tài)主方向的應(yīng)力應(yīng)變關(guān)系為

σ1σ2τ12=C11C120C21C22000C66ε1ε2γ12

（28）

式中：Cij為剛度矩陣。因此，經(jīng)典連續(xù)介質(zhì)理論的應(yīng)變能密度為

WCL=12Cijεiεj=12C11ε21+C12ε1ε2+

12C22ε22+12C66γ212

（29）

式中：i和j均取1、2、6，ε6=γ12。

PD得到的應(yīng)變能密度應(yīng)等于經(jīng)典連續(xù)介質(zhì)理論的應(yīng)變能密度，可求得cm、cf、dm和df分別為

cf=2（E1-E2）（1-ν12ν21）Qq=1ξVq

df=2G12（1-ν12ν21）-（1-ν12）E2Qq=15G12ξ25G12ξ2+6E1Δ2×3（1-ν12ν21）Qq=1Vqξ

cm=6（1+ν12）E2（1-ν12ν21）tπδ3

dm=λ（1-3ν12）E26πδt（1-ν12ν21）

（30）

式中：ν為泊松比;λ為調(diào)整系數(shù)，

λ=ΦΦ-arctan Φ，?Φ=δΔ5G126E2

（31）

2.2?失效準(zhǔn)則

對于二維平面問題，破壞每單位斷裂面積的所有鍵所需要的功等于材料的臨界應(yīng)變能釋放率，由此可以計(jì)算臨界拉伸。[30]對于正交各向異性單向板PD模型，鍵失效定義為纖維鍵和基體鍵失效。在纖維方向和基體方向的鍵的臨界伸長示意見圖5。對于纖維方向的鍵，如果鍵伸長超過臨界伸長sft，則纖維鍵發(fā)生斷裂;對于其他方向的所有鍵，如果鍵伸長超過臨界伸長smt，則基體鍵發(fā)生斷裂。令均勻介質(zhì)中單位斷裂區(qū)域的所有鍵斷裂所需要的功等于應(yīng)變能釋放率GI，Cr，1（即垂直于纖維方向的Ι型斷裂應(yīng)變能釋放率[31]），可計(jì)算臨界伸長sft;等效相應(yīng)的功和應(yīng)變能釋放率GI，Cr，2（即沿著纖維方向的Ι型斷裂應(yīng)變能釋放率[31]），可計(jì)算臨界伸長smt。臨界伸長sft和smt可以表示為

sft=4πGI，Cr，19E1δ，?smt=4πGI，Cr，29E2δ

（32）

節(jié)點(diǎn)xi與其近場區(qū)域內(nèi)的任意點(diǎn)xj之間有鍵的相互作用，當(dāng)鍵的伸長超過臨界伸長s0時(shí)，鍵會(huì)發(fā)生斷裂。定義標(biāo)量函數(shù)μ描述鍵是否發(fā)生斷裂，見式（33）。當(dāng)鍵伸長小于臨界伸長時(shí)μ=1，意味著鍵未發(fā)生斷裂，否則鍵斷裂并且相應(yīng)的鍵力變?yōu)?。

μ（ξ，t）=1，s（ξ，t′）

（33）

將節(jié)點(diǎn)xi的損傷定義為函數(shù)φ，見式（34）。若在其近場區(qū)域內(nèi)沒有鍵的斷裂，則損傷值為0;當(dāng)與節(jié)點(diǎn)xi連接的所有鍵都斷開時(shí)，節(jié)點(diǎn)的損傷值為1;當(dāng)損傷函數(shù)達(dá)到0.5時(shí)產(chǎn)生裂縫。

φ（x，t）=1-∫H μ（ξ，t）dV′∫HdV′

（34）

3?數(shù)值結(jié)果

3.1?碳纖維單向板拉伸

為驗(yàn)證模型的有效性，對不同纖維方向的單向板進(jìn)行數(shù)值模擬，分別考慮纖維方向θ為0°、15°、30°、45°和90°的單向板，板長為100 mm、寬為20 mm、厚為1 mm。初始坐標(biāo)設(shè)置在薄板左邊的中心點(diǎn)，在左側(cè)建立寬度為3Δ的區(qū)域作為約束區(qū)域，右端施加靜力拉伸載荷，力的總和為F=42 000 N。纖維方向的彈性模量E1=41.0 GPa，基體方向的彈性模量E2=10.4 GPa，泊松比ν12=0.28，剪切模量G12=4.3 GPa，密度ρ=1 970 kg/m3。近場區(qū)域均勻離散為1 mm×1 mm正方形物質(zhì)點(diǎn)的集合，近場區(qū)域設(shè)置為δ=3Δ。使用OpitiStruct進(jìn)行數(shù)值仿真模擬，有限元計(jì)算同樣也采用1 mm×1 mm正方形網(wǎng)格，將均勻分布在矩形板上的7個(gè)作用點(diǎn)作為測量點(diǎn)，見圖6。

為驗(yàn)證模型的有效性，繪制具有不同纖維方向的正交各向異性單向板沿x軸的位移云圖，與有限元分析進(jìn)行對比，結(jié)果見圖7～11。

PD模型位移云圖與有限元模型位移云圖結(jié)果吻合較好，說明該模型準(zhǔn)確。

3.2?碳纖維單向板裂紋擴(kuò)展

基于本文提出的BPD模型，引入含有中心裂紋缺陷的碳纖維單向板，模擬含缺陷碳纖維單向板的漸進(jìn)損傷，通過研究單向板的損傷擴(kuò)展路徑和最終破壞模式，驗(yàn)證該方法的可行性和有效性。

考慮長為70.0 mm、寬為40.0 mm、厚為0.6 mm的單向板，初始坐標(biāo)設(shè)置在薄板的中心點(diǎn)，纖維方向的彈性模量E1=105 GPa，基體方向的彈性模量E2=8.4 GPa，泊松比ν12=0.32，剪切模量G12=4 GPa，密度ρ=1 800 kg/m3。近場區(qū)域設(shè)置為δ=3.015Δ，均勻離散為0.5 mm×0.5 mm正方形物質(zhì)點(diǎn)的集合。缺陷在試件的幾何中心，裂紋長度2a=8.0 mm，在兩側(cè)施加2.0 mm的位移載荷。含線裂紋碳纖維單向板受拉伸載荷作用示意見圖12。纖維鍵的臨界伸長sft=0.03，基體鍵的臨界伸長smt=0.02。纖維角度為45°時(shí)，其漸進(jìn)損傷結(jié)果對比見圖13。

當(dāng)纖維方向?yàn)?5°時(shí)，基體首先損傷且損傷開始于初始裂紋的兩端;隨著位移的增大，裂紋沿著纖維方向向兩端擴(kuò)展，當(dāng)其達(dá)到纖維的臨界伸長時(shí)，纖維發(fā)生斷裂。PD仿真結(jié)果與試驗(yàn)結(jié)果[32]較吻合，可驗(yàn)證該模型的有效性和可行性，同時(shí)也表明PD可以模擬碳纖維單向板裂紋損傷，凸顯PD理論在

模擬裂紋損傷方面的優(yōu)越性。

4?結(jié)?論

提出一種改進(jìn)的BPD模型用于克服泊松比受到限制的問題，假設(shè)節(jié)點(diǎn)受到軸向和橫向力的作用，考慮節(jié)點(diǎn)的旋轉(zhuǎn)以消除由橫向力引起的額外彎矩，從而確保該模型滿足角動(dòng)量守恒條件。由于鍵的高度與長度的比值沒有小到足以忽略的程度，因此在該模型中考慮鍵的剪切影響系數(shù)，提高其精度。

通過模擬不同纖維方向的碳纖維單向板的單向拉伸，驗(yàn)證改進(jìn)的碳纖維單向板BPD模型的準(zhǔn)確性和解除泊松比限制的特點(diǎn)。通過模擬不同纖維角度矩形板的裂紋傳播路經(jīng)，驗(yàn)證模型在裂紋損傷模擬中的可行性和準(zhǔn)確性。對比結(jié)果表明，數(shù)值模擬結(jié)果和有限元仿真結(jié)果均與試驗(yàn)結(jié)果一致，可得出以下結(jié)論：

（1）改進(jìn)的鍵基PD模型通過在鍵的兩個(gè)節(jié)點(diǎn)上加入成對的橫向力，可消除泊松比的限制。這是因?yàn)橛?個(gè)獨(dú)立的PD參數(shù)對應(yīng)于復(fù)合材料4個(gè)宏觀材料常數(shù)，分別為彈性模量、泊松比、剪切模量和密度。

（2）通過考慮節(jié)點(diǎn)的旋轉(zhuǎn)，確保該模型滿足角動(dòng)量守恒，是該模型能夠精確模擬變形的關(guān)鍵。

（3）鍵的剪切影響系數(shù)可以提高單軸拉伸下板件的橫向變形精度，對模擬裂紋損傷和擴(kuò)展至關(guān)重要。

參考文獻(xiàn)：

[1]?杜善義. 先進(jìn)復(fù)合材料與航空航天[J]. 復(fù)合材料學(xué)報(bào)， 2007， 24（1）： 1?12. DOI： 10.3321/j.issn：1000?3851.2007.01.001.

[2]?任鵬. 復(fù)合材料層合板沖擊損傷及沖擊后疲勞壽命研究[D]. 哈爾濱： 哈爾濱工業(yè)大學(xué)， 2014.

[3]?LIEW K M， ZHAO X， FERREIRA A J M. A review of meshless methods for laminated and functionally graded plates and shells[J]. Composite Structures， 2011， 93（8）： 2031?2041. DOI： 10.1016/j.compstruct.2011.02.018.

[4]?NGUYEN V P， RABCZUK T， BORDAS S， et al. Meshless methods： A review and computer implementation aspects[J]. Mathematics and Computers in Simulation， 2008， 79（3）： 763?813. DOI： 10.1016/j.matcom.2008.01.003.

[5]?MUROTANI K， YAGAWA G， CHOI J B. Adaptive finite elements using hierarchical mesh and its application to crack propagation analysis[J]. Computer Methods in Applied Mechanics and Engineering， 2013， 253： 1?14. DOI： 10.1016/j.cma.2012.07.024.

[6]?SCHLLMANN M， FULLAND M， RICHARD H A. Development of a new software for adaptive crack growth simulations in 3D structures[J]. Engineering Fracture Mechanics， 2003， 70（2）： 249?268. DOI： 10.1016/S0013?7944（02）00028?0.

[7]?EL KHAOULANI R， BOUCHARD P O. An anisotropic mesh adaptation strategy for damage and failure in ductile materials[J]. Finite Elements in Analysis and Design， 2012， 59： 1?10. DOI： 10.1016/j.finel.2012.04.006.

[8]?BELYTSCHKO T， BLACK T. Elastic crack growth in finite elements with minimal remeshing[J]. International Journal for Numerical Methods in Engineering， 1999， 45（5）： 601?620.

[9]?COMBESCURE A， GRAVOUIL A， GRéGOIRE D， et al. X?FEM a good candidate for energy conservation in simulation of brittle dynamic crack propagation[J]. Computer Methods in Applied Mechanics and Engineering， 2008， 197（5）： 309?318. DOI： 10.1016/j.cma.2007.04.011.

[10]?SHI J X， CHOPP D， LUA J， et al. Abaqus implementation of extended finite element method using a level set representation for three?dimensional fatigue crack growth and life predictions[J]. Engineering Fracture Mechanics， 2010， 77（14）： 2840?2863. DOI： 10.1016/j.engfracmech.2010.06.009.

[11]?GENIAUT S， GALENNE E. A simple method for crack growth in mixed mode with X?FEM[J]. International Journal of Solids and Structures， 2012， 49（15/16）： 2094?2106. DOI： 10.1016/j.ijsolstr.2012.04.015.

[12]?GRGOIRE D， MAIGRE H， RTHOR J， et al. Dynamic crack propagation under mixed?mode loading： Comparison between experiments and X?FEM simulations[J]. International Journal of Solids and Structures， 2007， 44（20）： 6517?6534. DOI： 10.1016/j.ijsolstr.2007.02.044.

[13]?PRABEL B， MARIE S， COMBESCURE A. Using X?FEM method to model dynamic propagation and arrest of cleavage cracks in ferritic steel[J]. Engineering Fracture Mechanics， 2008， 75 （10）： 2984?3009. DOI： 10.1016/j.engfracmech.2008.01.008.

[14]?MENK A， BORDAS S P A. Crack growth calculations in solder joints based on microstructural phenomena with X?FEM[J]. Computational Materials Science， 2011， 50（3）： 1145?1156. DOI： 10.1016/j.commatsci.2010.11.014.

[15]?COX B N， GAO H J， GROSS D， et al. Modern topics and challenges in dynamic fracture[J]. Journal of Mechanics and Physics of Solids， 2005， 53（3）： 565?596. DOI： 10.1016/j.jmps.2004.09.002.

[16]?KRNER E. Elasticity theory of materials with long range cohesive forces[J]. International Journal of Solids and Structures， 1967， 3（5）： 731?742. DOI： 10.1016/0020?7683（67）90049?2.

[17]?SILLING S A. Reformulation of elasticity theory for discontinuities and long?range forces[J]. Journal of Mechanics and Physics of Solids， 2000， 48（1）： 175?209. DOI： 10.1016/S0022?5096（99）00029?0.

[18]?SILLING S A， ASKARI E. A meshfree method based on peridynamic model of solid mechanics[J]. Computers & Structures， 2005， 83（17/18）： 1526?1535. DOI： 10.1016/j.compstruc.2004.11.026.

[19]?WALTER G， NICOLAS S， STEWART S. Peridynamic modeling of plain and reinforced concrete structures[C]// Proceedings of 18th International Conference on Structural Mechanics in Reactor Technology（SMiRT 18）. Bejing： IASMiRT， 2005.

[20]?GERSTLE W， SAU N， SILLING S. Peridynamic modeling of concrete structures[J]. Nuclear Engineering and Design， 2007， 237（12/13）： 1250?1258. DOI： 10.1016/j.nucengdes.2006.10.002.

[21]?REN H L， ZHUANG X Y， RABCZUK T. A new peridynamic formulation with shear deformation for elastic solid[J]. Journal of Micromechanics and Molecular Physics， 2016， 1（2）： 1?28. DOI： 10.1142/S2424913016500090.

[22]?GHAJARI M， IANNUCCI L， CURTIS P. A peridynamic material model for analysis of dynamic crack propagation in orthotropic media[J]. Computer Methods in Applied Mechanics and Engineering， 2014， 276： 431?452. DOI： 10.1016/j.cma.2014.04.002.

[23]?ZHU Q Z， NI T. Peridynamic formulations enriched with bond rotation effects[J]. International Journal of Engineering Science， 2017， 121： 118?129. DOI： 10.1016/j.ijengsci.2017.09.004.

[24]?ZHOU X P， WANG Y T， SHOU Y D， et al. A novel conjugated bond linear elastic model in bond?based peridynamics for fracture problems under dynamic loads[J]. Engineering Fracture Mechanics， 2018， 188： 151?183. DOI： 10.1016/j.engfracmech.2017.07.031.

[25]?SILLING S A， EPTON M， WECKNER O， et al. Peridynamic states and constitutive modeling[J]. Journal of Elasticity， 2007， 88（2）： 151?184. DOI： 10.1007/s10659?007?9125?1.

[26]?TIMOSHENKO S P. LXVI： On correction for shear of differential equation for transverse vibrations of prismatic bars[J]. Philosophical Magazine， 1921， 41（245）： 744?746. DOI： 10.1080/14786442108636264.

[27]?KIM J U， RENARDY Y. Boundary control of Timoshenko beam[J]. SIAM Journal on Control and Optimization， 1987， 25（6）： 1417?1429. DOI： 10.1137/0325078.

[28]?REDDY J N. On dynamic behaviour of Timoshenko beam finite elements[J]. Sadhana， 1999， 24（3）： 175?198. DOI： 10.1007/BF02745800.

[29]?CHANDRUPATLA T R. Introduction to finite elements in engineering[M]. Gewerbestrasse： Pearson Schweiz Ag， 2010： 69.

[30]?HA Y D， BOBARU F. Studies of dynamic crack propagation and crack branching with peridynamics[J]. International Journal of Fracture， 2010， 162（1/2）： 229?244. DOI： 10.1007/s10704?010?9442?4.

[31]?PINEDA E， WAAS A. A Thermodynamically?based mesh objective work potential theory for predicting intralaminar progressive damage and failure in fiber?reinforced laminates [C]// Proceedings of 53rd AIAA/ASME/ASCE/AHS/ASC Structures， Structural Dynamics and Materials Conference. Honolulu： AIAA， 2012. DOI： 10.2514/6.2012?1612.

[32]?劉貞谷. 基于近場動(dòng)力學(xué)的復(fù)合材料單向板漸進(jìn)損傷模擬與實(shí)驗(yàn)研究[D]. 哈爾濱： 哈爾濱工業(yè)大學(xué)， 2017.

（編輯?武曉英）