陸曉敏，袁 濤

(河海大學(xué)工程力學(xué)系，南京 210098)

運(yùn)用塊體單元法的疊層梁應(yīng)力分析及破壞模擬

陸曉敏，袁 濤

(河海大學(xué)工程力學(xué)系，南京 210098)

應(yīng)用平面塊體單元法建立由剛體位移、常應(yīng)變、曲率組成的變形模式，推導(dǎo)了系統(tǒng)的平衡方程，開發(fā)了相應(yīng)的數(shù)值模擬程序，計(jì)算了簡支梁在均布荷載作用下的位移和應(yīng)力。結(jié)果表明:在較少單元情況下位移和應(yīng)力非常接近于理論解;相對于有限單元法，塊體單元法的位移精度更高，一般是理論值的上限;此外，接觸面中的應(yīng)力精度很高，確保了模擬接觸面破壞的可靠性。為模擬疊層結(jié)構(gòu)變形的非線性，在塊體之間的接觸面內(nèi)引入材料的非線性，分別考慮剪切屈服和開裂2種破壞形式。屈服判據(jù)采用M-C準(zhǔn)則，開裂采用最大拉應(yīng)力強(qiáng)度準(zhǔn)則，并編制了計(jì)算機(jī)程序。最后，計(jì)算模擬了一疊層式、受均布荷載懸臂梁的應(yīng)力、變形及層間剪切的破壞。結(jié)果表明:此方法能較好地模擬結(jié)構(gòu)的不連續(xù)變形及破壞。

疊層梁;塊體單元法;數(shù)值模擬

目前，疊層結(jié)構(gòu)的復(fù)合材料在生活和工程中逐漸得到廣泛應(yīng)用，如家具材料中的復(fù)合板、木工板、工程結(jié)構(gòu)中的復(fù)合梁等。這些材料具有層狀結(jié)構(gòu)，其破壞的特點(diǎn)常呈現(xiàn)層間錯(cuò)動(dòng)和開裂、層內(nèi)彎曲斷裂等，與均勻連續(xù)各向同性的材料有明顯的不同。所以，建立起能反映層狀結(jié)構(gòu)變形及破壞特點(diǎn)的數(shù)值模擬方法具有一定的工程應(yīng)用價(jià)值。

數(shù)值模擬是力學(xué)中研究結(jié)構(gòu)變形和安全度的一種有效方法。其中有限單元法［1］是最常用的一種數(shù)值方法，但其在處理層狀結(jié)構(gòu)時(shí)常需進(jìn)行特殊處理，異常麻煩，尤其是其常用的位移模式不能很好地反映層體材料彎曲變形的特點(diǎn)。在塊體單元法中，由于其以塊體的變形為求解未知量，塊體與塊體之間有單獨(dú)處理的接觸面單元，因而在模擬非連續(xù)變形時(shí)具有一定的優(yōu)勢。本文據(jù)此采用塊體單元法來研究疊層結(jié)構(gòu)復(fù)合材料的變形、應(yīng)力及破壞形態(tài)。

塊體單元法的基本思想:將固體分割成塊體系統(tǒng)，假設(shè)塊體的變形模式，在塊體之間人為引入接觸面單元;接觸面單元的變形由相鄰塊體的相對變形而決定;系統(tǒng)的平衡方程由最小勢能原理導(dǎo)出。早期的塊體單元法采用剛體位移模式，塊體的變形常作簡化處理。剛體位移模式的塊體單元法增加了系統(tǒng)的剛度，求得的變形量往往偏小，因而也影響了計(jì)算的精度。有文獻(xiàn)提出了一些簡化處理的方法［2-3］，但過于簡單，尤其是沒有考慮材料變形的泊松比效應(yīng)，給實(shí)際工程結(jié)構(gòu)分析帶來不少困難。

本文借鑒DDA的思想［4-5］，把塊體的應(yīng)變作為獨(dú)立未知量，并在塊體的變形模式中進(jìn)一步引入曲率作為獨(dú)立變量，建立起具有彎曲變形特點(diǎn)的塊體單元法理論，并在接觸面內(nèi)引入材料剪切屈服和開裂的非線性模型，使之能很好地模擬層狀結(jié)構(gòu)的彎曲變形。此法的優(yōu)點(diǎn)一方面能提高計(jì)算精度;另一方面，由于考慮了塊體的變形，從而可以模擬材料非線性［6］的本構(gòu)關(guān)系，尤其在塊體單元的接觸面中可以較方便地引入各種材料性質(zhì)，如彈塑性、黏性流變等。本文算例表明:塊體單元法理論具有較高的精度，能很方便地模擬層狀結(jié)構(gòu)的層間屈服、開裂的破壞過程。

1 平面塊體單元法基本理論

1.1 塊體的變形模式

對于平面的塊體元，塊體的位移和變形模式可以用3個(gè)剛體位移、3個(gè)平均應(yīng)變和2個(gè)曲率來表示，用矩陣表示為

式中:ui，vi為i號塊體在整體坐標(biāo)系下x，y軸向的位移;θi為繞形心C的轉(zhuǎn)角，在小變形的前提下這些位移是微小量;εxi，εyi，γxyi為i號塊體在整體坐標(biāo)下的平均線應(yīng)變和切應(yīng)變;κx'i，κy'i為i號塊體在單元局部坐標(biāo)系Cx'y'下的平均彎曲曲率，為使用方便，曲率定義在局部坐標(biāo)系中。約定局部坐標(biāo)系中的x'軸與單元某邊平行，此方向常定義為層狀結(jié)構(gòu)的接觸面方位。曲率與曲率半徑的關(guān)系為

式中ρx'i，ρy'i分別為i號塊體在x'和y'軸向的平均曲率半徑。

如圖1所示:i號塊體的形心坐標(biāo)用(xCi，yCi)表示;塊體上某一點(diǎn)A的整體坐標(biāo)為(xi，yi)，局部坐標(biāo)為(x'i，y'i);變形前局部坐標(biāo)軸x'與x軸的夾角用φi表示。在小變形假設(shè)的前提下，A點(diǎn)的位移可表示為

式中Ni稱為塊體的形函數(shù)矩陣，其各元素為:

1．2 塊體的應(yīng)變及應(yīng)力

根據(jù)式(3)可以將塊體內(nèi)任一點(diǎn)A的位移用

式中Bi稱為i號塊體單元應(yīng)變矩陣，各元素為:

假設(shè)塊體是各向同性的材料，彈性模量為E，泊松比為υ，則考慮平面應(yīng)力時(shí)的彈性矩陣為

任一點(diǎn)A的應(yīng)力為

圖1 塊體單元與局部坐標(biāo)示意圖

1．3 接觸面單元的應(yīng)力和應(yīng)變

在塊體單元法中，人為在塊體接觸面之間引入很薄的接觸面單元，見圖2。設(shè)在塊體i和j之間有一薄層單元(簡稱為接觸面單元，編號為k)，厚度為h(很小)，C為形心，β為接觸面單元的局部坐標(biāo)x'軸的方位角。接觸面上的初始接觸點(diǎn)B和B'的位移可以用式(3)分別進(jìn)行計(jì)算。

將位移轉(zhuǎn)換至接觸面的局部坐標(biāo)系中，則

式中Ti為局部坐標(biāo)與整體坐標(biāo)間的轉(zhuǎn)換矩陣。假設(shè)沿接觸面法向(y'軸)的切應(yīng)變和線應(yīng)變均勻變化，x'向的線應(yīng)變?yōu)閴K體i，j在此線應(yīng)變的平均值，即

圖2 接觸面單元與局部坐標(biāo)示意圖

據(jù)此可以推導(dǎo)出接觸面單元內(nèi)的應(yīng)變矩陣ε'k為

式中:(x'B'，y'B')為B'點(diǎn)在塊體i局部坐標(biāo)系中的坐標(biāo);(x'B，y'B)為B點(diǎn)在塊體j局部坐標(biāo)系中的坐標(biāo)。式(14)中B'k稱為k號縫面單元的應(yīng)變矩陣。

假設(shè)接觸面單元材料的彈性模量為E'，泊松比為υ'，引入彈性矩陣

則接觸面內(nèi)的應(yīng)力為

1．4 系統(tǒng)的平衡方程

系統(tǒng)的平衡方程可利用最小勢能原理導(dǎo)出。系統(tǒng)的內(nèi)力勢能由塊體的變形勢能和接觸面單元的變形勢能兩部分組成。塊體i單元的變形勢能為:

式中kBi稱為塊體i的單元?jiǎng)哦染仃?，具體計(jì)算式為

接觸面單元k的變形勢能為

式中kCk為接觸面單元k的勁度矩陣，其計(jì)算式為

系統(tǒng)總的彈性變形勢能為:

式中:Nc，NB分別為系統(tǒng)塊體單元和接觸縫面單元的個(gè)數(shù);d為系統(tǒng)的未知量(位移和應(yīng)變)列陣;K稱為系統(tǒng)的勁度矩陣，由所有塊體單元和接觸面單元的勁度矩陣組裝而成。

系統(tǒng)的外力勢能為

式中:Fi，pi，qi分別為作用在i塊體上的集中力、面力和體力矢量(矩陣);為相應(yīng)的等效荷載列陣

式中:s，v分別為面力作用的面積和體力作用的體積。式(21)中R為系統(tǒng)的等效荷載列陣，由各塊體單元的等效荷載列陣組裝而成

1．5 彈塑性塊體單元法

在塊體元模型中，既可以考慮塊體材料的塑性，又可以考慮接觸縫面的塑性?？紤]到對于疊層梁的屈服和開裂破壞主要發(fā)生在結(jié)構(gòu)面內(nèi)，故本文暫只考慮縫面材料的非線性。如接觸面單元屈服，則在計(jì)算其勁度矩陣和縫面應(yīng)力時(shí)，需用彈塑性矩陣Dep代替彈性矩陣D'。設(shè):縫面材料的抗拉強(qiáng)度為σt;抗壓強(qiáng)度為σc;粘結(jié)強(qiáng)度為c;摩擦因數(shù)為f。縫面單元的破壞形式主要是剪切屈服和法向開裂，剪切屈服采用M-C準(zhǔn)則判別，開裂由法向拉應(yīng)力判別。

如果接觸面單元法向應(yīng)力σy'＞σt，接觸面單元開裂，則認(rèn)為其剛度喪失。如果

則縫面發(fā)生剪切屈服。根據(jù)相關(guān)聯(lián)的流動(dòng)法則可導(dǎo)出其塑性矩陣［7］:

式中:λ=E+f2G，f=tanφ，δ=sign(τx'y')，E，G，υ分別為接觸面材料的彈性模量、剪切模量和泊松比。

當(dāng)接觸縫面單元進(jìn)入屈服或開裂狀態(tài)后，支配方程將成為非線性方程。此時(shí)，勁度矩陣K將是位移δ的函數(shù)。為了模擬加載過程，并加快收斂速度，本文對非線性方程的求解采用子增量變剛度迭代法，在每級荷載增量中剛度矩陣保持不變，并進(jìn)行迭代計(jì)算。

2 數(shù)值模擬結(jié)果

2.1 算例1:受均布荷載的簡支梁的線彈性計(jì)算分析

如圖3所示，取一單位厚度的矩形截面簡支梁，頂部受均勻分布荷載作用。分布荷載的集度為10 kN/m2;梁的跨度為l=8 m;截面高度為h= 0.8 m?？缰薪孛鏋閥軸(向下為正);中性層為x軸。梁材料的彈性模量E=2.0×104MPa，泊松比υ=0.20。計(jì)算時(shí)接觸面單元厚度人為取1.0× 10-3m，網(wǎng)格水平方向?yàn)?0等分，豎直為8層等分，塊體單元的邊長為0.2 m×0.1 m。

利用材料力學(xué)［8-14］的方法計(jì)算梁跨中撓度的值，本文計(jì)算的最大撓度發(fā)生在跨中，其值為6.411 1×10-3m，與材料力學(xué)法結(jié)果的誤差為0.549%。塊體元變形計(jì)算結(jié)果略偏大的原因主要是人為引入接觸面單元降低了梁的剛度。

圖3 簡支梁網(wǎng)格與荷載約束圖

表1給出了靠1/4跨截面處(x=-1.1 m)塊體單元形心的x，y向正應(yīng)力和切應(yīng)力的計(jì)算值及對應(yīng)的理論解［15］。比較塊體元結(jié)果與理論解可以發(fā)現(xiàn)兩者非常接近，誤差不超過1.0%，說明塊體元計(jì)算結(jié)果正確且精度足夠。

表2給出了靠1/4跨截面處(x=-1.0 m)接觸面單元形心的x，y向正應(yīng)力和切應(yīng)力的計(jì)算值及對應(yīng)的理論解［10］。比較兩者可以發(fā)現(xiàn):接觸面單元內(nèi)的切應(yīng)力及法向應(yīng)力與理論解非常接近，誤差不超過1.0%;但接觸面方向(y向)的正應(yīng)力與理論解有一定的誤差，變化規(guī)律相同，最大誤差為14.7%。這說明:塊體元法中接觸面內(nèi)的切應(yīng)力和法向應(yīng)力計(jì)算結(jié)果可靠且精度足夠，接觸面內(nèi)切向正應(yīng)力誤差略大。

表1 截面x=-1.1 m塊體形心應(yīng)力分量

表2 截面x=-1.0 m接觸縫面形心應(yīng)力分量

2.2 算例2:受均布荷載的懸臂疊層梁屈服非線性分析

如圖4所示，單位厚度的矩形截面懸臂梁長6.0 m，截面高2.0 m，由4層等厚的板通過粘結(jié)膠合疊加而成。網(wǎng)格圖見圖10，板材料的彈性模量均取8.0 MPa，泊松比取0.22，接觸面材料彈性模量和泊松比的取值與板材一致。模擬接觸面2種特性:①粘結(jié)力取0 kPa，抗拉強(qiáng)度取2.0 MPa，摩擦角取25°;②粘結(jié)力取300 kPa，抗拉強(qiáng)度取2.0 MPa，摩擦角取25°。計(jì)算時(shí)僅考慮接觸面的非線性。

圖4 懸臂梁網(wǎng)格荷載約束圖

①屬于無粘結(jié)的疊層梁問題，接觸面上有法向約束力及摩擦產(chǎn)生的切應(yīng)力。由于相對滑動(dòng)，粘結(jié)面兩側(cè)的軸向正應(yīng)力及位移不再連續(xù)。計(jì)算了5種荷載情況:P=30，80，120，150，200 kN/m2。圖5給出了懸臂端最大撓度與荷載的關(guān)系(撓度向下定義為負(fù)值)其中:“”為整體梁計(jì)算結(jié)果(不考慮接觸面的非線性);“”為疊層梁的計(jì)算結(jié)果。由變形結(jié)果可見:接觸面由于摩擦屈服，梁軸向線應(yīng)變沿截面高度不再連續(xù)，梁的剛度明顯降低。當(dāng)取均布荷載集度為150 kN/m2時(shí)，懸臂端的最大撓度分別為 -0.737 8 cm(疊梁)和-0.675 45 cm(整體梁)，疊梁的撓度增加了10.0%。圖6給出梁變形后的網(wǎng)格圖。從圖6可以明顯看出:粘結(jié)層處變形不連續(xù)，層與層之間有相互錯(cuò)動(dòng)。

圖5 懸臂梁端部結(jié)點(diǎn)位移分量與荷載圖

圖6 懸臂梁網(wǎng)格變形圖(P=150 kN/m2)

圖7、8給出了截面A(截面位置x=4.20 m)處結(jié)點(diǎn)在5種不同荷載情況下的軸向正應(yīng)力和切應(yīng)力，數(shù)據(jù)見表3。軸向正應(yīng)力圖表明:由于粘結(jié)層的滑移，軸向應(yīng)力在粘結(jié)層兩側(cè)不再連續(xù);隨著荷載增加，突變量也隨之增加，在層內(nèi)正應(yīng)力基本呈線性變化。切應(yīng)力圖表明:在粘結(jié)層摩擦屈服后，切應(yīng)力降低，不再呈拋物線變化。

圖7 截面A結(jié)點(diǎn)正應(yīng)力與荷載圖

圖8 截面A結(jié)點(diǎn)切應(yīng)力與荷載圖

表3 截面A(x=4.02 m)結(jié)點(diǎn)應(yīng)力分量

②屬于有粘結(jié)的疊合梁問題，接觸面上有法向約束力及摩擦和粘結(jié)產(chǎn)生的切應(yīng)力。同樣計(jì)算了5種荷載情況:P=150，200，300，500，800 kN/m2。圖9給出了P=150 kN/m2時(shí)疊層梁的變形及粘結(jié)層的屈服區(qū)范圍。由圖9可見:粘結(jié)屈服區(qū)范圍及疊層梁變形比無粘結(jié)疊層梁小。由此得出:有粘結(jié)的疊層梁結(jié)構(gòu)相對更穩(wěn)定。

圖10、11給出了截面A(截面位置x=4.20 m)處單元形心在5種荷載情況下的軸向正應(yīng)力和切應(yīng)力。軸向正應(yīng)力圖表明:粘結(jié)層一旦屈服，軸向應(yīng)力在粘結(jié)層兩側(cè)不再連續(xù)，切應(yīng)力也有所降低。

圖9 懸臂梁網(wǎng)格變形圖(P=150 kN/m2)

圖10 截面A塊體單元形心正應(yīng)力圖

圖11 截面A塊體單元形心切應(yīng)力圖

2.3 算例3:受集中荷載兩端固定疊層梁開裂非線性分析

將上述梁的兩端改為固定端約束，在第2層(從底到高編號)作用2個(gè)向下對稱的集中荷載F(見圖12)，網(wǎng)格示意圖見圖12。粘結(jié)面的抗拉強(qiáng)度取2.0 MPa，其他參數(shù)不變。隨著荷載的增加，簡支疊層梁中間的粘結(jié)層的最大法向應(yīng)力接近材料抗拉強(qiáng)度，接觸面層開裂。圖12～14給出了3種不同荷載集度網(wǎng)格變形及粘結(jié)層的開裂情況。由圖可知:隨著荷載的逐漸增大，粘結(jié)層的開裂范圍相繼擴(kuò)大。圖15給出梁跨中底部的豎向位移分量與荷載的關(guān)系。結(jié)果表明:粘結(jié)層一旦發(fā)生開裂，梁的擾度就會(huì)發(fā)生突變。

圖12 荷載作用位置與網(wǎng)格變形圖(P=2 000 kN)

圖13 荷載作用位置與網(wǎng)格變形圖(P=3 000 kN)

圖14 荷載作用位置與網(wǎng)格變形圖(P=4 000 kN)

圖15 懸臂梁跨中底部豎向位移分量與荷載關(guān)系圖

3 結(jié)束語

本文研究表明:在塊體的變形模式中引入應(yīng)變和曲率后能顯著提高計(jì)算位移和應(yīng)力的精度，且接觸面的應(yīng)力也有足夠的精度，確保了模擬接觸面材料非線性力學(xué)行為的可靠性。

考慮接觸面材料的屈服和開裂后，非線性塊體單元法能較好地模擬有結(jié)構(gòu)面材料的非線性變形及破壞規(guī)律，如疊層復(fù)合板、層狀巖體等。本文算例驗(yàn)證了非線性塊體單元法的可靠性和合理性。如在塊體內(nèi)引入材料的非線性，使塊體單元法模擬的結(jié)果更加接近實(shí)際。

［1］陳國榮.有限單元法原理及應(yīng)用［M］.北京:科學(xué)出版社，2009.

［2］任青文，余天堂.塊體單元法的理論和計(jì)算模型［J］.工程力學(xué)，1999，16(1):67-77.

［3］任青文，余天堂.邊坡穩(wěn)定的塊體單元法分析［J］.巖石力學(xué)和工程學(xué)報(bào)，2001，20(1):20-24.

［4］石根華.數(shù)值流形方法與非連續(xù)變形分析［M］.北京:清華大學(xué)出版社，1997.

［5］陸曉敏.邊坡的彈粘塑變形及穩(wěn)定性研究［J］.巖石力學(xué)與工程學(xué)報(bào)，2002，21(4):493-496.

［6］殷有泉.非線性有限元基礎(chǔ)［M］.北京:北京大學(xué)出版社，2007.

［7］李詠偕，施澤華.塑性力學(xué)［M］.北京:水利電力出版社，1987.

［8］孫訓(xùn)芳.材料力學(xué)［M］.北京:高等教育出版社，2003.

［9］TIMOSHENKO S P.高等材料力學(xué)［M］.北京:科學(xué)出版社，1979.

［10］陳杰.考慮層間接觸變形時(shí)疊層梁層間接觸壓力分析［J］.力學(xué)與實(shí)踐，2001(4):45-47.

［11］黃傳躍，諸德超.基于剖面翹曲修正理論的復(fù)合材料疊層梁動(dòng)態(tài)特性分析［J］.航空學(xué)報(bào)，1996(5):96 -101.

［12］周一勤.夾層疊層梁的內(nèi)力計(jì)算及其應(yīng)用［J］.華東公路，1991(5):36-39.

［13］鄧梁波.地基上復(fù)合材料疊層梁的穩(wěn)定性［J］.華南理工大學(xué)學(xué)報(bào):自然科學(xué)版，1995(4):77-80.

［14］VIJAY K G.非守恒荷載作用下疊層梁的動(dòng)力穩(wěn)定性［J］.鋼結(jié)構(gòu)，2009(2):74.

［15］徐芝倫.彈性力學(xué)［M］.北京:高等教育出版社，1988.

(責(zé)任編輯陳 艷)

Analysis of Stress and Deformation of Layered Beam Based on Block Element Method

LU Xiao-min，YUAN Tao
(Department of Engineering Mechanics，Hohai University，Nanjing 210098，China)

Based on block element method，deformation mode included rigid displacement，constant strains and curvatures was constructed.The equilibrium equations were conducted by principle of minimum potential energy.The displacement and stress under blanced load of simple beam were calculated.The results of block element method show that stress and displacement are very close to theoretical values only with fewer elements，and the stresses in interface elements have more accuracy yet alway it will be the upper limit of theoretical value.Besides，the stress accuraly is very high，which can ensule the reliability in simulating the damage of the contact surface.To analyze nonlinear material mechanics，yielding and cracking in interface elements were simulated.Yielded criterion usedM-Cformulation principle and cracked condition used the principle that the max-tensile stress in interfaceelement exceeds maximum tensile strength.Using the computer program of upper method，the stress and deformation of a layered cantilever beam were studied under uniform pressure.The results are correct and show that the method of the paper is feasible.

layered beam;block element method;numerical analysis

TB301

A

1674-8425(2015)11-0051-09

10.3969/j.issn.1674-8425(z).2015.11.009

2015-06-22

國家自然科學(xué)基金重點(diǎn)資助項(xiàng)目(11132003)

陸曉敏(1963—)，男，江蘇無錫人，博士，教授，主要從事水工結(jié)構(gòu)、巖土工程等的力學(xué)分析及穩(wěn)定性研究;袁濤(1989—)，男，江蘇鎮(zhèn)江人，碩士研究生，主要從事計(jì)算力學(xué)與工程仿真研究。

陸曉敏，袁濤.運(yùn)用塊體單元法的疊層梁應(yīng)力分析及破壞模擬［J］.重慶理工大學(xué)學(xué)報(bào):自然科學(xué)版，2015 (11):51-59.

format:LU Xiao-min，YUAN Tao.Analysis of Stress and Deformation of Layered Beam Based on Block Element Method［J］.Journal of Chongqing University of Technology:Natural Science，2015(11):51-59.