基于木材彈塑性損傷本構(gòu)的古建木結(jié)構(gòu)殘損梁柱構(gòu)件損傷非線性分析

張利朋，謝啟芳，吳亞杰，劉伊津

(西安建筑科技大學 土木工程學院；結(jié)構(gòu)工程與抗震教育部重點實驗室，西安 710055)

古建筑木結(jié)構(gòu)是中國重要的文化遺產(chǎn)。木梁和木柱組成古建筑木結(jié)構(gòu)的承重骨架，并維持其整體穩(wěn)定。然而，由于數(shù)百年來環(huán)境與荷載的共同作用，木梁和木柱發(fā)生了不同程度的殘損，如整體老化、局部腐朽、蟲蛀和開裂等，不僅導致其受力性能退化，也使古建筑木結(jié)構(gòu)的整體力學性能受到不同程度的削弱。因此，研究殘損梁柱構(gòu)件的受力行為可為現(xiàn)存殘損古建筑木結(jié)構(gòu)的保護提供科學依據(jù)[1]。

近年來，學者們圍繞殘損梁柱構(gòu)件的調(diào)查實測和受力性能分析開展了深入研究，并取得了一定的研究成果。在殘損梁柱構(gòu)件的調(diào)查實測方面，周乾[2]通過對故宮古建筑木結(jié)構(gòu)柱子的勘察實測，指出木柱的典型殘損問題主要為局部糟朽和開裂，并進行了成因分析；王曉麗[3]利用快速普查法調(diào)查了飛云樓的承重木柱和梁枋的殘損情況，確定了各構(gòu)件的殘損等級，為古建筑木結(jié)構(gòu)的保護及修繕提供了依據(jù)；喬冠峰[4]從承重木柱、承重梁枋、結(jié)構(gòu)整體變形三個方面入手，對古建筑木結(jié)構(gòu)的殘損點評估界限進行了探討，并對殘損界限評定的適用條件提出了建議。在殘損梁柱構(gòu)件的受力分析方面，Guan等[5-6]基于彈塑性本構(gòu)模型，采用單元生死技術研究了跨中開方孔木梁的受力性能和方孔角部的損傷演化；湯永紅[7]基于ANSYS模擬了縱向開裂木柱的軸壓性能，并將一階屈曲特征位移向量的2%作為初始幾何變形進行了縱向開裂木柱的非線性分析；Zhang等[8]通過人工預制貫通縫的方式模擬了木柱的干縮裂縫，并通過軸心和偏心受壓試驗研究了失效模式和殘余極限承載力的變化規(guī)律；朱忠漫[9]基于ABAQUS對縱向開裂木柱的受力性能進行了參數(shù)分析，首先對木柱開展線性屈曲模擬，進一步采用弧長法進行非線性屈曲分析，將第一、二階位移特征向量按照一定比例加權(quán)后，作為初始幾何缺陷輸入有限元模型；Mosallam[10]對人工預制殘損木梁進行了抗彎試驗，分析了其力學性能退化規(guī)律；陳立濤[11]針對實際結(jié)構(gòu)中存在的局部腐朽和蟲蛀殘損，通過預制局部缺口受彎木梁的受彎試驗和數(shù)值模擬，研究了殘損對木梁受力性能的影響規(guī)律；路鵬[12]采用材料性能折減法在ABAQUS中建立了老化木柱的有限元模型，并分析了不同老化程度的影響；Li等[13]通過人工開槽方式開展了柱腳糟朽木構(gòu)架的試驗，并采用有限單元法對其進行了數(shù)值模擬。謝啟芳等[14]基于ABAQUS模擬了局部帶缺口木柱的軸壓性能，通過分析缺口的大小和位置，考察了初始缺口缺陷對木柱承載力和非線性的影響。

可見，有關古建筑木結(jié)構(gòu)梁柱構(gòu)件的研究主要集中在開裂殘損方面，而關于木梁、木柱局部糟朽和蟲蛀影響下受力性能方面的研究較少，并且該方面的研究大多基于試驗方法開展，僅少數(shù)研究采用了有限單元法，且其采用的本構(gòu)模型均為ABAQUS自帶的彈塑性模型，導致其結(jié)果無法合理反映局部殘損部位與完好部位交界處的材料損傷行為。通常情況下，局部殘損梁柱構(gòu)件在局部交界面的破壞涉及損傷的演化、應力的軟化，而這是ABAQUS自帶彈塑性模型所無法直接反映的，所建有限元模型便無法較好地反映整個梁柱構(gòu)件的受力非線性和損傷演化過程。筆者通過建立木材的彈塑性損傷本構(gòu)模型，分析帶局部缺陷梁柱構(gòu)件的損傷演化、發(fā)展過程，并基于課題組試驗進行驗證。研究成果可為殘損梁柱構(gòu)件的受力分析及性能評價提供借鑒。

1 木材彈塑性損傷本構(gòu)模型

1.1 理論模型

σ=Cd(d):(ε-εp)

(1)

式中：Cd(d)為木材的彈性損傷剛度張量；ε、εp和σ分別為木材的總應變張量、塑性應變張量和應力張量；d為木材的損傷變量集，如式(2)所示。

d={d1t,d1c,d2t,d2c,d3t,d3c,d12s,d13s,,drs}

(2)

式中：dit和dic(i=1,2,3)為木材在各主軸方向的受拉和受壓損傷變量；d12s、d13s和d23s分別為各剪切面內(nèi)的剪切損傷變量。

由式(1)可知，彈塑性損傷本構(gòu)模型的關鍵在于確定塑性應變和損傷變量的演化規(guī)律。

1.1.2 塑性部分 木材的塑性屈服函數(shù)為

(3)

木材強度準則選用Hoffman準則[15]。

(4)

式中：Pα和pα可展開為式(5)和式(6)。

(5)

pα=[α11α22α330 0 0]T

(6)

式中：αij為與木材單軸受拉和受壓強度有關的參數(shù)，其表達形式和確定方法參見文獻[16]。

(7)

采用相關性流動法則，則木材的塑性應變演化方程可表達為

(8)

1.1.3 損傷部分 木材彈塑性損傷本構(gòu)模型的損傷部分需給出損傷演化起始條件和損傷變量演化方程。

參考復合材料的損傷建模理論，木材損傷準則可以考慮失效機制并表達為[17]

fI(φI,rI)=φI-rI≤0

(I=1t,1c,2t,2c,2v,3t,3c,3v)

(9)

式中：φI為與木材損傷機制有關的函數(shù)，采用Sandhaas[18]提出的木材損傷準則，如式(10)所示。

(10)

式(9)中，rI為損傷閾值函數(shù)，其表達式為

(11)

本文模型采用的損傷演化方程為[15-16]

J={1t;1c;2t;2c;3t;3c;12s;13s;rs}

(12)

gJ=GJ/Lc

(13)

式中：L為單元特征長度，由ABAQUS主程序計算后賦值給變量CELENT，具體計算方法參見ABAQUS幫助文檔[19]；GL為沿單元特征長度的應變能釋放，N/mm。

(14)

對于增量過程，式(14)可顯化表達為

(15)

1.2 本構(gòu)模型的數(shù)值算法

建立的彈塑性損傷本構(gòu)模型由隱式的后退歐拉算法進行求解，該方法通過在每個增量步結(jié)束時更新內(nèi)變量以避免求解結(jié)果從屈服面漂移。

有效應力空間中的積分算法可以表達為

εn+1=εn+Δεn+1

(16)

1.2.1 彈性預測 根據(jù)式(16)，可得

(17)

(18)

增量步開始時，總應變增量全部被視作彈性應變，據(jù)此計算試探應力，如果試探應力在屈服面內(nèi)，則說明試算過程為彈性過程，有效應力就等于試探應力，否則應進行塑性修正。

1.2.2 塑性修正 經(jīng)過初始彈性預測，可得塑性應變增量

(19)

將式(19)代入式(16)，可得

(20)

通過求解式(20)方程組，可得

(21)

式中:

(22)

將式(20)代入式(19)，可得

(23)

因此，可求得塑性乘數(shù)增量

(24)

進一步可得有效應力張量、塑性應變張量和等效塑性應變的表達式。

(25)

1.2.3 損傷修正 基于塑性修正過程求得的有效應力狀態(tài)，可計算損傷內(nèi)變量，進而對有效應力進行折減，即可求得名義應力。由于該過程不涉及迭代計算，計算過程簡單，不再展開。

1.3 本構(gòu)模型驗證

通過將所建立的木材本構(gòu)模型基于ABAQUS的材料二次開發(fā)接口UMAT進行程序?qū)崿F(xiàn)，可借助其計算模塊實現(xiàn)二次開發(fā)。文獻[16]給出了詳細的本構(gòu)算法和程序開發(fā)過程，且已基于多種加載工況下的材性試驗(順紋及橫紋方向的單軸受拉和受壓、順紋及橫紋重復受壓等[20-21])進行了模型驗證，模擬結(jié)果與試驗結(jié)果吻合較好?；谌c受彎木梁在構(gòu)件層面對所建本構(gòu)模型和所開發(fā)的本構(gòu)程序進行模型驗證。

三點受彎木梁模型采用的材料和試驗數(shù)據(jù)取自Khennane等[22]的試驗，木材為云杉，木梁橫截面尺寸為60 mm×60 mm，長1 600 mm，如圖1所示。

圖1 三點受彎木梁尺寸及試驗加載裝置[22]Fig.1 Size and test setup of three-point bending timber

為了防止木梁局部壓潰，試驗過程中，在施加集中力處放置了60 mm × 60 mm鋼墊塊，采用位移控制加載，直到試件破壞。試驗結(jié)束時，試件底部跨中發(fā)生了脆性受拉斷裂。

1.3.1 幾何模型、網(wǎng)格劃分、接觸設置 三點受彎木梁有限元模型采用與試驗相同的尺寸，梁頂面的加載墊塊和梁底面兩端支座處采用截面尺寸為60 mm×60 mm，厚20 mm的墊塊進行模擬。選取C3D8R單元類型，網(wǎng)格長度設置為10 mm。確保墊塊與梁的網(wǎng)格節(jié)點重合以提高收斂性。如果在鋼墊塊與梁頂部接觸面之間采用綁定約束(Tie)，則會在鋼墊塊邊緣與梁頂接觸處引起應力集中，且損傷將首先在此產(chǎn)生，并進一步引起收斂性問題，而采用考慮摩擦系數(shù)的方法則可以避免這一問題[23]。因此，墊塊與梁之間的接觸面相互作用采用法向硬接觸和切向摩擦系數(shù)法?？紤]到墊塊與梁頂面之間發(fā)生相對位移的可能性很小，取摩擦系數(shù)為0.45。

表1 本構(gòu)模型參數(shù)確定Table 1 Determination of wood properties

1.3.3 邊界條件及加載方式 有限元模型的邊界條件按試驗情況設為簡支，左端約束梁的水平位移和豎向位移，右端僅約束豎向位移。采用位移加載模式加載至計算過程終止。

1.3.4 模擬結(jié)果分析 由三點受彎木梁有限元模擬結(jié)果的失效模式與試驗結(jié)果對比(順紋受拉損傷變量SDV20)，如圖2所示。由圖2可見，順紋受拉損傷分布與試驗結(jié)果較為一致，均為梁跨中截面(圖3中的A點)發(fā)生。此時，SDV20的最大值為0.81，計算過程中斷。通過所建立的木材彈塑性損傷本構(gòu)模型可以較好地預測三點受彎木梁的脆性斷裂失效。

圖2 三點受彎木梁有限元模擬與試驗結(jié)果對比Fig.2 Comparison of finite element simulation and experimental results of three-point bending wood

圖3 三點受彎木梁底部損傷單元的順紋受拉應力應變曲線與應力云圖分析Fig.3 Analysis of stress-strain curve and stress nephogram of damage element along grain at the bottom of three-point

1.4 黏性調(diào)整對收斂性的影響

以簡單受拉試件為例，考察黏性規(guī)則化對脆性斷裂損傷引起收斂困難的改善，主要討論兩個問題：一是有無黏性調(diào)整對收斂性的影響；二是黏性系數(shù)取值范圍的問題。所選試件模型如圖4所示，構(gòu)件尺寸為210 mm × 10 mm × 10 mm，所采用的主要材料參數(shù)選自文獻[16]，順紋抗拉強度為24 MPa，順紋彈性模量為11 000 MPa，黏性規(guī)則化系數(shù)為0.000 1。具體操作方法為：沿構(gòu)件長度方向選定部分單元為初始損傷單元，并通過人工降低損傷單元材料性能(強度設置為23 MPa)的方式形成斷裂帶，實現(xiàn)初始損傷位置的預設，從而排除軟化段對網(wǎng)格敏感性的影響而單獨研究黏性調(diào)整對收斂性的影響。

圖4 黏性調(diào)整對單向受拉木構(gòu)件收斂性的影響Fig.4 Influence of viscosity adjustment on convergence of unidirectional tension wood

由圖4可知，未進行黏性調(diào)整時，模型僅能預測木材順紋受拉峰值強度，而無法反映其斷裂過程；黏性調(diào)整后，則可以較好地預測斷裂路徑?？梢姡ば哉{(diào)整有利于在使用彈塑性損傷本構(gòu)模型程序計算構(gòu)件損傷破壞過程時提高其收斂性。

圖5 黏性規(guī)則化系數(shù)的優(yōu)化取值分析Fig.5 Optimization value analysis of viscosity

2 局部殘損梁柱構(gòu)件的數(shù)值模型

2.1 局部殘損木梁的數(shù)值模型

2.1.1 殘損木梁幾何模型 建立的殘損木梁有限元模型的幾何模型選自陳立濤進行的殘損木梁試驗[11]，其幾何尺寸如圖6所示。

圖6 殘損木梁幾何模型Fig.6 Geometric model of damaged wooden

圖7 底部缺口殘損木梁的破壞模式Fig.7 Failure mode of timber beam with bottom

2.1.2 殘損木梁有限元模型 根據(jù)上述帶缺口木梁的幾何條件建立有限元模型，其邊界條件、加載位置及支座布置、材料模型、分析步設置等與試驗基本等效。有限元模型的邊界條件按試驗情況設為簡支，左端約束梁的水平位移和豎向位移，右端僅約束豎向位移。采用位移加載模式加載至計算過程終止。區(qū)別之處僅在于網(wǎng)格劃分模塊，殘損木梁有限元模型在缺口附近的網(wǎng)格尺寸劃分較小。

2.1.3 殘損木梁的破壞模式 底部缺口殘損木梁有限元模型的順紋剪切損傷(SDV24)發(fā)展和演化情況如圖8所示。由圖8可見，使用所開發(fā)的木材彈塑性損傷本構(gòu)模型可較好地預測其損傷的發(fā)展和演化過程，但由于使用隱式材料程序(UMAT)，該模型不具有顯式有限元程序的網(wǎng)格刪除功能，因此，裂縫的開展情況只能通過損傷進行分析而無法直觀地在有限元模型中表達。

圖8 底部缺口殘損木梁有限元模型的剪切損傷發(fā)展和演化Fig.8 Shear damage development and evolution of fintte element model of damaged wood beam with bottom

圖9 殘損木梁荷載位移曲線Fig.9 Load-displacement curves of damaged wooden

由圖9可見，建立的木材彈塑性損傷本構(gòu)模型可以較好地預測帶有局部缺口的四點受彎木梁的承載能力，但由于模型假定在材料各個方向達到屈服前均為線性，導致構(gòu)件模擬的剛度偏大。

2.2 局部殘損木柱的數(shù)值模型

建立的殘損木柱有限元模型所采用的幾何模型選自王玄[24]進行的殘損木柱試驗。邊界條件、加載位置及支座布置、材料模型、分析步設置等與試驗基本等效。采用位移加載模式加載至計算過程終止。帶缺口木柱有限元模型的順紋受壓損傷(SDV21)發(fā)展和演化情況如圖10所示。由圖10可見，使用建立的木材彈塑性損傷本構(gòu)模型可較好地預測其損傷的發(fā)展和演化過程，但由于所使用的是隱式材料程序，不具有網(wǎng)格刪除功能，因此，局部壓屈的發(fā)展情況只能通過損傷進行分析而無法直觀地在有限元模型中進行表達。

圖10 局部缺口木柱破壞模式的有限元與試驗結(jié)果對比Fig.10 Comparison of finite element and experimental results of failure modes of wood columns with local

殘損木柱豎向荷載位移曲線的有限元分析結(jié)果與試驗結(jié)果的對比情況如圖11所示。由圖11可見，建立的彈塑性損傷本構(gòu)模型可以較好地預測其承載能力，且在局部損傷出現(xiàn)之前與試驗結(jié)果吻合較好。但在木柱損傷逐漸趨近于缺口附近時，有限元結(jié)果出現(xiàn)較為明顯的軟化段，總體帶有脆性失效特征，因此無法反映出試驗結(jié)果的延性失效過程。

圖11 殘損木柱豎向荷載位移曲線Fig.11 Vertical load-displacement curve of

3 結(jié)論

1)在連續(xù)介質(zhì)力學和損傷力學框架內(nèi)建立了木材的彈塑性損傷本構(gòu)模型，可用于模擬三點受彎木梁的受彎損傷模式及分布規(guī)律。

2)所開發(fā)的彈塑性損傷本構(gòu)模型可用于反映局部殘損梁柱構(gòu)件在殘損界面與完好界面處的損傷行為，可較好地反映殘損梁柱構(gòu)件受力過程中的非線性行為和損傷演化。

3)從文中梁柱構(gòu)件的3個模擬實例來看，建立的木材彈塑性損傷本構(gòu)模型在應用于殘損梁柱構(gòu)件的數(shù)值模擬時，無法完全使分析過程收斂到最終，這是由于建立的木材彈塑性損傷本構(gòu)模型并未完全解決塑性與損傷的耦合問題：文中采用的強度準則和屈服準則都是有效應力的函數(shù)，二者在本質(zhì)上相同，即損傷準則某種程度上也是強度準則，強度準則某種程度上也是損傷準則，都是對線性區(qū)域的限定，但由于其表達式不同，所限定的區(qū)域并不相同，從而出現(xiàn)了塑性與損傷的“打架”問題。該問題需在后續(xù)研究中進一步解決。