用于ABAQUS梁單元的混凝土單軸本構模型

王 強，常 凱，侯康康，路 炯，劉 琳

(沈陽建筑大學土木工程學院，遼寧沈陽 110168)

0 引 言

在強震作用下，復雜高層建筑結構會進入塑性狀態(tài)，因此進行結構的動力彈塑性時程分析對于評估結構的抗震安全性具有重要意義。整體結構分析時，采用梁單元能夠以較少的計算量合理描述梁柱構件的力學性態(tài)。梁單元模型通常有塑性鉸和纖維模型2類，纖維梁單元模型直接從材料單軸應力-應變本構關系和截面構造配筋出發(fā)，相對于塑性鉸模型，其能夠更準確地描述梁柱構件的軸力-雙向彎曲耦合效應，在工程界已得到了廣泛認可[1-2]。

通用有限元軟件ABAQUS擁有強大的非線性分析能力，在分析大型復雜結構時得到了廣泛應用。然而目前該軟件缺乏用于三維纖維梁單元的混凝土本構模型，無法采用梁單元描述鋼筋混凝土梁柱桿件的非線性力學行為，因此采用該軟件進行建筑結構的動力彈塑性時程分析存在困難。為此，已有學者通過ABAQUS的用戶子程序接口，基于Scott-Kent-Park的混凝土本構模型對用于三維梁單元的混凝土單軸本構模型進行二次開發(fā)[3-4]。

2008年汶川地震后，中國修訂了《混凝上結構設計規(guī)范》(GB 50010—2010，以下簡稱《規(guī)范》)，給出了考慮損傷的混凝土本構關系。采用該本構關系進行結構的彈塑性分析對于中國工程界具有重要的現(xiàn)實意義。目前，已有部分學者基于該本構模型進行了二次開發(fā)，如胡曉斌等[5]、齊虎等[6]，但前者并未給出受拉段的滯回規(guī)則且模擬結果與試驗結果差距較大，而后者缺少實例數(shù)據(jù)對比。因此，本文對《規(guī)范》附錄C混凝土本構關系進行了深入研究，首先補充了受拉加卸載滯回規(guī)則與考慮箍筋對混凝土的約束作用，并修正了受壓損傷恢復的問題，隨后利用ABAQUS用戶子程序接口UMAT，VUMAT進行二次開發(fā)，編寫了用于顯、隱式分析的梁單元混凝土單軸本構模型子程序，用于鋼筋混凝土梁柱構件在多維受力狀態(tài)下的彈塑性分析，滿足結構動力彈塑性分析的需求。

1 纖維梁單元

基于材料單軸本構關系的纖維模型是將梁柱構件沿縱向劃分為若干子段，再沿構件橫截面劃分成纖維束。ABAQUS中的B31纖維梁單元構型如圖1所示。每個纖維(以截面積分點表征)只考慮它的軸向本構關系，且可定義不同的本構關系。構件橫截面變形符合平截面假定。對截面纖維的當前狀態(tài)積分就可以得到截面的雙向抗彎剛度、雙向抵抗矩以及軸力，進而沿構件長進行積分，可以得到精確的構件單元剛度矩陣。纖維模型可以自然、簡單地描述構件的雙向彎曲-軸力耦合效應[7-8]。

2 材料的本構模型

纖維模型的準確度在很大程度上依賴其使用的材料本構模型與材料參數(shù)的合理性。為準確描述梁柱構件中鋼筋與混凝土的材料性能，本文采用2個共節(jié)點的梁單元來分別模擬構件中的縱向鋼筋與混凝土材料，且縱向鋼筋通過等面積、等位置原則代換成箱形截面單元，從而分別定義混凝土和鋼筋的本構模型與材料參數(shù)。

2.1 鋼筋的本構模型

大部分學者在模擬鋼筋混凝土構件受力性能時，鋼筋本構通常采用ABAQUS自帶本構，但該本構對于包興格效應考慮不足。因此，本文鋼筋本構采用自編單軸本構，可較為充分地考慮鋼材的包興格效應以及屈服、硬化與軟化特性，其中本構的骨架曲線采用Esmaeily-Xiao模型[9]，加卸載曲線采用Légeron模型[10]，具體圖形如圖2所示。圖2中，σ為應力，ε為應變。該本構模型的準確性驗證和模型的參數(shù)取值方法可參考文獻[11]。

2.2 混凝土的本構模型

2.2.1 受壓區(qū)骨架曲線與滯回規(guī)則

《規(guī)范》給出了素混凝土的受壓區(qū)骨架曲線及滯回規(guī)則，如圖3所示，其中，fc,r為混凝土抗壓強度代表值，εc,r為混凝土受壓峰值應力對應的應變，σun，εun分別為受壓混凝土骨架曲線開始卸載時的應力和應變，Er為受壓混凝土卸載/再加載的變形模量，εz為受壓混凝土卸載至零應力點時的殘余應變。骨架曲線具體公式參見文獻[12]。

(1)受壓區(qū)骨架曲線上的應力-應變關系

受壓區(qū)骨架曲線上的應力-應變關系采用混凝土單軸受壓損傷演化參數(shù)dc表征，即

σ=(1-dc)E0ε

(1)

(2)

(3)

(4)

(5)

式中：E0為素混凝土的彈性模量；αc為單軸受壓應力-應變曲線下降段參數(shù)。

文獻[13]指出，鋼筋混凝土梁柱構件中由于箍筋對核心區(qū)混凝土的約束作用，使得混凝土的受壓應力-應變曲線隨著約束效應的增大由明顯的陡峰向平緩飽滿的曲線過渡。因此在構件滯回性能模擬時有必要考慮箍筋的約束作用。本文參照文獻[11]，采用約束指標λt來考慮鋼筋混凝土梁柱構件中箍筋對混凝土的約束作用，對《規(guī)范》中的混凝土應力-應變曲線方程進行修正。λt計算公式如下

(6)

式中：fyv為箍筋屈服強度；ρv為混凝土梁柱構件中箍筋的體積配箍率。

此時，約束混凝土的抗壓強度fcc、峰值應變εpc及曲線下降段參數(shù)αdc為

fcc=(1+0.5λt)fc,r

(7)

εpc=(1+2.5λt)εc,r

(8)

αdc=(1-1.75λt)αc

(9)

將fcc，εpc，αdc代替式(2)～(5)中相應的fc,r，εc,r，αc，則可得梁柱構件中考慮箍筋約束作用的混凝土應力-應變骨架曲線。

考慮箍筋約束作用的混凝土與非約束混凝土的應力-應變曲線對比如圖4所示。由圖4可以看出，當柱構件中混凝土考慮箍筋約束作用時，混凝土抗壓強度、峰值應變等都有所提高，尤其是曲線下降段形狀趨于平緩。

(2)受壓區(qū)滯回規(guī)則的修正

《規(guī)范》已給出了重復荷載作用下受壓混凝土卸載與再加載路徑方程，即

σ=Er(ε-εz)

(10)

(11)

(12)

(13)

式中：εca為附加應變。

已有研究表明，對于混凝土材料，重復加卸載條件下的材料損傷是不可恢復的，即受壓混凝土卸載/再加載的變形模量單調降低[13-14]。當卸載應變較大時，采用式(11)～(13)確定的受壓混凝土卸載/再加載的變形模量Er，會出現(xiàn)隨卸載應變增大而增大的現(xiàn)象，即混凝土受壓損傷恢復。

由圖5可以看出，采用《規(guī)范》公式，在卸載應變超過A點時出現(xiàn)了變形模量Er增大的情況。A點應變?yōu)?.87εc，即為εc/(εc+εun)=0.09εun/εc時對應的卸載應變εun。因此，本文對式(13)進行了修正，采用式(14)計算附加應變εca，即

(14)

計算結果表明，按公式(14)確定的附加應變得到的變形模量Er單調下降，能夠合理表征循環(huán)加載下混凝土材料的受壓損傷不可恢復特性。

2.2.2 受拉區(qū)骨架曲線與滯回規(guī)則

《規(guī)范》中給出了受拉區(qū)骨架曲線，如圖6所示，其中，ft,r，εt,r分別為受拉區(qū)應力、應變，具體公式參見文獻[12]，但是《規(guī)范》中并沒有給出明確的受拉區(qū)滯回規(guī)則。目前王瑾等[15-16]在OpenSEES中采用Concrete02模型在模擬混凝土構件滯回性能的分析中取得了較好的結果。本文參照OpenSEES中Concrete02模型，采用原點指向規(guī)則，對《規(guī)范》混凝土本構補充混凝土受拉區(qū)滯回規(guī)則。

混凝土在拉-壓與壓-拉往復加載條件下，本文采用拉區(qū)應力-應變曲線隨受壓殘余應變點遷移的滯回規(guī)則。若混凝土初始受力即為受拉，則沿受拉區(qū)骨架曲線加載。若受拉卸載，卸載曲線由受拉卸載點B指向初始受拉點A，如圖7所示。隨后受壓加載，則沿受壓區(qū)骨架曲線加載，加載至C點時卸載至受壓殘余應變點D點。繼續(xù)受拉加載時其路徑為DE(DE∥BA)，至F點進入受拉卸載段時，依據(jù)原點指向規(guī)則，其路線為指向受壓殘余應變點D的FD直線。受壓再加載按DC路線回到受壓曲線上的C點。

3 用戶自定義材料子程序的實現(xiàn)

ABAQUS軟件提供了基于隱式算法的UMAT材料子程序接口和基于顯式算法的VUMAT材料子程序接口，分別應用于隱式求解器(Standard)和顯式求解器(Explicit)[17]。

采用前述的混凝土與鋼筋本構關系，筆者基于UMAT接口進行二次開發(fā)，編制了材料子程序，并嵌入ABAQUS/Standard模塊中。由于隱式求解必須引入雅可比矩陣，而雅可比矩陣和N+1時刻的應力、應變以及狀態(tài)變量是實時更新的，這導致材料非線性發(fā)展較為強烈時會出現(xiàn)雅可比矩陣奇異而產(chǎn)生隱式求解的數(shù)值不收斂問題。因此，本文又基于VUMAT接口，編制了用于顯式分析的材料子程序。在進行梁柱構件乃至整體結構的非線性分析時，對于豎向加載工況下結構非線性發(fā)展較弱的工況，采用隱式求解方法；對于水平加載等結構非線性發(fā)展比較強烈的工況則采用顯式分析方法。通過顯、隱式結合來完成結構在不同荷載工況下的接力計算。

4 算例驗證

為驗證本文所采用混凝土本構的合理性，分別對文獻[18]中鋼筋混凝土柱試件邊柱A及文獻[19]中鋼筋混凝土柱試件TP74，TP77進行數(shù)值模擬。

4.1 試件邊柱A

本文首先模擬了文獻[18]中軸壓-單軸受彎鋼筋混凝土柱試件邊柱A在低周往復荷載下的滯回性能。

4.1.1 試驗概況

試件構造和配筋如圖8所示，加載圖示見圖9，加載規(guī)則如圖10所示，材料參數(shù)取值見表1。

4.1.2 結果分析

試驗與模擬結果對比如圖11所示?？梢钥闯?，當不考慮鋼筋混凝土柱中箍筋對混凝土的約束作用時，本文采用的混凝土本構模型已能夠較好地模擬混凝土柱的滯回性能，滯回曲線飽滿程度與試驗結果接近，但承載力略低。當考慮混凝土柱中箍筋對混凝土的約束作用時，計算所得承載能力有所提高，模擬結果與試驗結果更為吻合，表明本文所采用的混凝土本構模型較為合理。

表1 材料參數(shù)Tab.1 Material Parameters

圖12為在考慮箍筋約束情況下修正附加應變前后計算得出的滯回曲線對比?？梢钥闯?，2條曲線基本重合，證明本文對公式(13)的修正是可行的。

圖13為試件邊柱A在加載幅值10，25，45 mm時的單圈滯回曲線對比。由圖13可知：加載幅值為

10 mm時，滯回環(huán)峰值點與實測曲線峰值點較為接近，相對誤差達到2%，但計算曲線的單圈耗能能力(單圈滯回環(huán)包含的面積)略高于實測曲線的耗能能力；加載幅值達到25 mm時，計算曲線的滯回環(huán)峰值點略高于實測曲線，計算曲線的耗能能力略高于實測曲線的耗能能力，相對誤差低于16%。當加載幅值達到45 mm時，滯回環(huán)峰值點較實測曲線峰值點高約20%，且計算曲線較“瘦”，其單圈耗能能力略低于實測曲線，相對誤差約為4%。這可能與模擬分析未考慮試件的剪切非線性、軸壓力波動、黏結滑移等因素有關。此外，計算曲線能夠反映實測曲線的捏縮現(xiàn)象。

圖14為計算所得試件邊柱A底部鋼筋、混凝土單元截面中的角部纖維應力-應變曲線。由圖14可知，在單軸往復加載過程中，混凝土纖維受壓已進入骨架曲線下降段，鋼筋纖維受拉已進入骨架曲線軟化段，且應力-應變曲線符合本文所構建的混凝土本構模型的滯回規(guī)則。

4.2 試件TP74，TP77

本文對文獻[19]中的試件TP74和TP77進行模擬分析。構件單元及截面劃分與試件邊柱A相同。試件及試驗概況詳見文獻[19]。

4.2.1 單軸加載試件TP74

圖15為TP74實測滯回曲線和采用本文構建的材料本構計算所得的滯回曲線對比，考慮箍筋約束和不考慮箍筋約束2種情況。可以看出，計算所得滯回曲線飽滿程度均與試驗結果較為接近，且考慮箍筋約束后計算所得承載力有所提高，與實測曲線更為吻合，證明了本文所構建的混凝土本構模型的合理性。

圖16為在考慮箍筋約束情況下計算得出的附加應變修正前后滯回曲線對比。通過對比發(fā)現(xiàn)，修正附加應變前后最大承載力一致，但修正后的承載力并無明顯退化現(xiàn)象，且滯回曲線與圖15(a)的實測滯回曲線更為吻合，這表明了本文修正附加應變的正確性。

4.2.2 雙軸加載試件TP77

圖17為試驗與數(shù)值模擬所得滯回曲線，由圖17可知，計算所得的承載力、強度退化程度、滯回曲線飽滿程度均與試驗結果較為接近。

此外，在試驗曲線與計算曲線中，在保持方向1位移恒定時，在方向2加載均會使方向1的承載力出現(xiàn)明顯的下降趨勢，反之亦然。此現(xiàn)象表明，本文提出的模型能夠合理地反映梁柱構件的雙向彎曲耦合性能。

5 結 語

(1)本文采用原點指向模型與拉區(qū)應力-應變曲線隨受壓殘余應變點遷移法，補充完善了《規(guī)范》混凝土本構模型的混凝土受拉段加卸載、拉壓受力狀態(tài)轉換的滯回規(guī)則。

(2)本文采用基于體積配箍率建立的約束指標λt考慮鋼筋混凝土梁柱構件中箍筋對混凝土的約束作用，對《規(guī)范》混凝土受壓應力-應變骨架曲線方程進行修正，有效提高了鋼筋混凝土梁柱構件滯回性能模擬結果的準確性。

(3)本文對《規(guī)范》混凝土本構模型中的附加應變算法進行了修正，所得的受壓混凝土卸載/再加載變形模量Er能夠合理地考慮循環(huán)加載下混凝土受壓損傷不可恢復特性。

(4)本文構建的混凝土本構模型及相應的用戶材料子程序能夠細致地描述梁柱構件受力過程中截面的內(nèi)力變化情況，準確地模擬鋼筋混凝土柱在低周往復荷載作用下的滯回性能，能進一步用于大型復雜結構的彈塑性反應分析。