閤東東， 張 純， 苗啟松

(1. 北京工業(yè)大學(xué) 建筑工程學(xué)院， 北京 100124； 2. 北京市建筑設(shè)計研究院有限公司， 北京 100045；3. 碧桂園集團北京區(qū)域公司， 北京 100070)

動力彈塑性時程分析是罕遇地震作用下結(jié)構(gòu)變形驗算的重要方法，也是結(jié)構(gòu)性能化設(shè)計的重要手段。在鋼筋混凝土結(jié)構(gòu)的動力彈塑性時程分析中，應(yīng)用最為廣泛的通用有限元軟件主要有ABAQUS和MSC.MARC(以下簡稱MARC)，由于這兩款軟件并非結(jié)構(gòu)工程專用分析軟件，采用其進行動力彈塑性時程分析需要進行纖維梁子程序二次開發(fā)。

已有部分學(xué)者對MARC在結(jié)構(gòu)動力彈塑性分析中的應(yīng)用進行了研究。陸新征等[1,2]基于MARC平臺開發(fā)的THUFIBER程序可以準確模擬地震作用下結(jié)構(gòu)的動力響應(yīng)。Han等[3]將MARC用于組合結(jié)構(gòu)的動力彈塑性分析。汪訓(xùn)流等[4]編制的混凝土和鋼筋單軸本構(gòu)分別采用Légeron-Paultre模型[5]和修正的Légeron模型[6]，對鋼筋混凝土柱的受力狀態(tài)進行了預(yù)測，與試驗結(jié)果吻合良好。聶建國等[7]也在MARC的基礎(chǔ)上，開發(fā)了纖維梁單元模型，用于分析鋼-混凝土組合結(jié)構(gòu)。柳國環(huán)等[8]也在MARC的基礎(chǔ)上，開發(fā)了可用于模擬鋼材和混凝土纖維梁單元的用戶子程序。Lu等[9]運用MARC成功的模擬了632 m高的上海中心大廈在強震作用下的連續(xù)倒塌過程。

為采用MARC對鋼筋混凝土結(jié)構(gòu)進行動力彈塑性時程分析和倒塌分析奠定基礎(chǔ)，本文開發(fā)了基于剛度法的纖維梁單元子程序，并與2011年中國建筑學(xué)會抗震防災(zāi)分會建筑結(jié)構(gòu)抗倒塌專業(yè)委員會進行了鋼筋混凝土框架結(jié)構(gòu)擬靜力倒塌試驗進行了對比，驗證了程序的正確性和可靠性。

1 MARC中的纖維梁單元

纖維模型的主要原理是將桿件截面離散化為若干個纖維，并假定構(gòu)件在各受力階段，整個截面滿足平截面假定，各纖維受力分別采用各自的單軸應(yīng)變-應(yīng)力曲線來描述，并根據(jù)相應(yīng)纖維材料的單軸應(yīng)變-應(yīng)力關(guān)系來計算截面內(nèi)力與變形的非線性關(guān)系。纖維梁單元通過截面方向和長度方向上兩次積分得到單元的剛度矩陣，通過材料級別的應(yīng)變-應(yīng)力關(guān)系實現(xiàn)精細化的有限元建模。纖維梁模型包含基于剛度法和柔度法的兩種模型，剛度法采用基于位移的插值函數(shù)，并忽略構(gòu)件橫截面剪力對構(gòu)件正截面受力的影響。目前在通用有限元軟件中，都采用基于剛度法的模型。

MARC提供的纖維梁子程序接口UBEAM，包含有98號可以考慮剪切變形的兩節(jié)點單積分點線性Timoshenco梁單元，和52號不考慮剪切變形的兩節(jié)點三積分點Euler-Bernoulli梁單元。在接口UBEAM中，只需要給出單元截面層次的剛度矩陣，即可實現(xiàn)粱式結(jié)構(gòu)的復(fù)雜非線性行為模擬。因此，在UBEAM中，可以采用纖維梁模型通過截面方向上的積分得到截面的剛度矩陣，也可以提供特征截面的塑性鉸模型。相比較而言，前者更容易實現(xiàn)，且更容易考慮軸力與雙向彎矩的耦合作用。下面主要對基于剛度法的纖維梁單元截面內(nèi)力與變形的關(guān)系進行介紹。

圖1所示為矩形構(gòu)件截面示意圖，截面被離散為混凝土纖維和鋼筋纖維。對于第i根纖維，其形心坐標為(yi,zi)，面積為Ai。首先進行截面分析，根據(jù)平截面假定，截面各點的軸向應(yīng)變εi可以表達為：

εi=ε0+κyyi+κzzi

(1)

式中：ε0為形心處的軸向應(yīng)變，以受拉為正；κy和κz分別為截面主軸的曲率。

圖1 矩形截面纖維模型示意

以增量形式表達的混凝土纖維和鋼纖維的應(yīng)力-應(yīng)變關(guān)系為：

Δσi=EiΔεi

(2)

式中：Ei為第i根纖維的切線模量。

截面上所有纖維的合力增量為：

(3)

(4)

式(4)可簡化為：

ΔF=DTΔε

(5)

式中:ΔF為與正應(yīng)力相關(guān)的截面內(nèi)力增量向量；DT為截面切線剛度矩陣；Δε為截面應(yīng)變增量向量。

動力時程分析時， MARC采用隱式積分算法和Newton-Raphson迭代法求解。迭代過程中，在每個增量步的開始，各單元上一步的截面應(yīng)變、應(yīng)力和應(yīng)變增量由主程序傳遞到UBEAM中，在增量步的結(jié)束時，UBEAM完成截面應(yīng)力和切線剛度矩陣的更新并返回主程序。在應(yīng)力更新的過程中，需要借助于一些狀態(tài)變量，并開辟數(shù)組空間記錄所有單元的纖維關(guān)鍵歷史應(yīng)力、應(yīng)變信息，并在每一步迭代收斂時進行更新。

2 混凝土和鋼筋材料應(yīng)變-應(yīng)力關(guān)系

2.1 混凝土材料本構(gòu)模型

混凝土采用了兩種單軸本構(gòu)模型，分別命名為Concrete01，Concrete02。Concrete01模型采用了OpenSees中的Concrete02模型。Concrete02為GB 50010-2010《混凝土結(jié)構(gòu)設(shè)計規(guī)范》中提供的模型[10]。兩種混凝土材料模型的循環(huán)加載曲線如圖2所示。

圖2 Concrete01與Concrete02混凝土材料本構(gòu)模型

2.2 鋼筋材料本構(gòu)模型

鋼筋單軸本構(gòu)共包括Steel01和Steel02兩種模型，其中Steel01為運用廣泛的雙線性隨動強化模型。Steel02采用OpenSees中的修正Menegotto-Pinto模型。

3 纖維梁子程序測試

3.1 試驗概況

2011年，中國建筑學(xué)會抗震防災(zāi)分會建筑結(jié)構(gòu)抗倒塌專業(yè)委員會進行了鋼筋混凝土框架結(jié)構(gòu)擬靜力倒塌試驗，其中對一組鋼筋混凝土框架柱進行了滯回試驗(圖3)，并公布了其試驗數(shù)據(jù)和材料參數(shù)[11]。試驗過程中，水平加載位置位于距柱底750 mm處，軸壓千斤頂球鉸距柱底1030 mm。底層中柱C的配筋如圖4所示。其中混凝土立方體強度為30.1 MPa。Φ10鋼筋屈服強度和極限強度分別為481和745 MPa，彈性模量為265433 MPa，屈服應(yīng)變?yōu)?.002，伸長率為23.6%；Φ8鋼筋屈服強度和極限強度分別為582和855 MPa，彈性模量為289850 MPa，屈服應(yīng)變?yōu)?.002，伸長率為28.8%。

圖3 鋼筋混凝土框架柱滯回試驗

圖4 底層中柱C配筋/mm

3.2 C柱測試結(jié)果

將柱單元沿長度方向上劃分為5個單元，柱截面纖維劃分示意圖如圖1所示，其中混凝土部分一共劃分了8×8=64個纖維，鋼筋根據(jù)實際鋼筋根數(shù)劃分了8根纖維。采用上述混凝土和鋼筋材料單軸本構(gòu)模型，基于UBEAM子程序接口，分別編制了52號Euler-Bernoulli梁和98號Timoshenko梁單元的截面本構(gòu)關(guān)系。采用MARC對C柱進行時程分析，持時與擬靜力加載試驗時間相同。對采用Euler-Bernoulli梁和Timoshenco梁單元，及不同混凝土和鋼筋單軸本構(gòu)時的結(jié)果進行了歸納總結(jié)，結(jié)果定義如下：E和T分別代表52號和98號梁的計算結(jié)果，Concrete01Steel01代表混凝土和鋼筋分別采用Concrete01和Steel01單軸本構(gòu)模型時的結(jié)果，其余工況依此類推。

中柱C的配筋結(jié)果詳細數(shù)據(jù)可見網(wǎng)站。其中位移加載時程和軸壓時程分別如圖5a，5b所示，加載持時265.427 s。

圖5 中柱C往復(fù)加載方案

圖6所示為52號Euler-Bernoulli梁數(shù)值模擬結(jié)果與試驗結(jié)果對比，從圖中可以看出，數(shù)值模擬結(jié)果與試驗結(jié)果較為一致。分別從圖6a，6b可以看出，采用相同的混凝土和不同的鋼筋單軸本構(gòu)模型時，得到的計算結(jié)果較為接近，說明混凝土單軸本構(gòu)對計算結(jié)果不會產(chǎn)生嚴重影響。圖7所示為98號Timoshenko梁數(shù)值模擬結(jié)果與試驗結(jié)果對比，可見數(shù)值模擬結(jié)果與試驗結(jié)果吻合較好，仔細對比圖6，7，可以發(fā)現(xiàn)當(dāng)計算參數(shù)相同時，52號梁和98號梁得到的計算結(jié)果略有區(qū)別。數(shù)值模擬結(jié)果與試驗結(jié)果較為一致，說明計算結(jié)果具有一定的可靠性。

圖6 采用Euler-Bernoulli梁時中柱C滯回曲線

圖7 采用Timoshenko梁時中柱C滯回曲線

為詳細和直觀地考察子程序的正確性，利用UBEAM接口輸出了混凝土和鋼筋纖維應(yīng)變、應(yīng)力以及切線剛度時程曲線進行分析。圖8所示為52號Euler-Bernoulli梁沿長度方向第一個積分點處，角部第一根混凝土纖維滯回曲線，從圖中可以看出，混凝土材料滯回曲線與預(yù)期的本構(gòu)模型一致，當(dāng)鋼筋本構(gòu)模型不同時，得到的混凝土滯回曲線略有區(qū)別。圖9所示為角部第一根鋼筋纖維應(yīng)變-應(yīng)力關(guān)系曲線，鋼筋材料滯回曲線與采用的模型一致，對比圖8可以看出，相對而言，因混凝土材料模型不同引起的鋼筋滯回曲線差別較因鋼筋材料模型不同引起的混凝土滯回曲線差別大。

圖8 Euler-Bernoulli梁角部第一根混凝土纖維滯回曲線

圖9 Euler-Bernoulli梁長度方向上第一個積分點滯回曲線

圖10輸出了沿長度方向第一個積分點處，角部第一根混凝土纖維的切線剛度時程曲線，從圖中可以看出不同鋼筋材料本構(gòu)得到的時程曲線較為一致，采用Concrete01模型時局部情況下出現(xiàn)了負剛度，采用Concrete02模型時，可以看出隨著構(gòu)件進入屈服，混凝土材料剛度呈下降趨勢。圖11所示為角部第一根鋼筋纖維的切線剛度時程曲線，從圖中可以看出，初始剛度約為2.8×105N/mm2。

圖10 采用Euler-Bernoulli梁時角部第一根混凝土纖維切線剛度時程曲線

圖11 采用Euler-Bernoulli梁時角部第一根鋼筋纖維切線剛度時程曲線

以上計算結(jié)果說明編制的UBEAM子程序是正確的，能較好地模擬往復(fù)荷載作用下鋼筋混凝土柱的承載力和變形性能的變化過程，在極限承載能力和極限變形能力計算上具有一定精度，得到的計算結(jié)果較為可靠。

4 結(jié) 論

(1)本文基于MSC.MARC提供的接口UBEAM，編制了52號三積分點Euler-Bernoulli梁單元，和98號單積分點線性Timoshenco梁單元纖維梁子程序，并可以分別采用兩種不同的鋼筋和混凝土材料單軸本構(gòu)。

(2)經(jīng)過鋼筋混凝土柱擬靜力往復(fù)加載試驗結(jié)果對比，表明編制的程序可以較好地模擬其承載力和變形的變化過程，得到的計算結(jié)果較為可靠。對子程序得到的結(jié)果進行了詳細討論，驗證了程序的正確性，并可進一步用于鋼筋混凝土結(jié)構(gòu)的動力彈塑性時程分析和倒塌數(shù)值模擬研究中。

(3)對于同時考慮壓、彎、剪、扭的彈塑性行為及應(yīng)變率的影響，這部分工作還需進一步研究。