火灾下钢筋混凝土柱非线性有限元分析

1、火灾下钢筋混凝土柱非线性有限元分析1丁发兴*，余志武(中南大学 土木建筑学院，湖南 长沙 410075)*E-mail:摘 要：本文建立了火灾下钢筋混凝土柱非线性有限元分析理论，编制了非线性有限元程序NFEMRCLF。首先提出了高温下混凝土单轴应力-应变滞回关系，采用了合理的钢材和混凝土的热-力耦合本构模型；然后基于火灾下U.L.列式虚功增量方程，采用非线性分层梁单元理论，给出了火灾下钢筋混凝土柱非线性有限元方程组的求解方法；最后对典型火灾下钢筋混凝土柱的试验资料进行双重非线性有限元分析，计算结果与试验结果较吻合。分析结果表明：火灾下钢筋混凝土柱的轴向变形-火灾时间曲线的计算结果基本上反映了钢

2、筋混凝土柱的变形特性，计算的耐火极限基本上是试验结果的上限，且采用材料应力-应变滞回关系的变形性能计算结果比采用材料应力-应变骨架关系的相应结果更合理。关键词：钢筋混凝土柱，火灾，应力-应变滞回关系，耐火极限，非线性分析，有限元1. 引 言随着社会经济的发展和城市化的进程，现代建筑的火灾发生频率和火灾危险性不断增加，建筑结构的抗火分析成为当前研究的热门课题。钢筋混凝土结构在我国应用最广泛，研究钢筋混凝土结构的抗火性能具有重要的现实意义。火灾下建筑结构的材料模型极为复杂，高温下钢材和混凝土单轴力学试验研究表明1：钢材和混凝土的高温本构关系随温度、应力和时间而有不同的变化规律，构成了复杂的热-力耦

3、合本构关系。由于温度场不均匀的影响，高温下钢筋混凝土梁跨中受拉区混凝土和钢筋混凝土柱形心处混凝土容易出现卸载和反向再加载情形，目前，笔者还没有见到反映这种情形的高温下混凝土应力-应变滞回关系的研究报道。已有钢筋混凝土柱抗火非线性分析中，分析方法包括模型柱法(也叫纤维模型法)24和有限元法58。采用的模型柱法，以升温加载的方式来计算恒载升温途径下钢筋混凝土柱的抗火性能，不能反映合理其加载途径；采用大型商业有限元程序三维建模5,7或自编的非线性梁单元模型1,6,8的有限元法，这些方法各有特点，但较少有学者将其进一步应用于局部火灾下钢筋混凝土足尺框架结构抗火分析中。笔者提出火灾下混凝土应力-应变滞回

4、关系，采用合理的钢材和混凝土的热-力耦合本构模型，并基于火灾下结构的U.L.列式虚功增量方程，采用平面非线性梁单元理论，给出火灾下钢筋混凝土柱非线性有限元方程组的求解方法，编制非线性有限元程序NFEMRCLF (non-linear finite element method of reinforced concrete columns under loading and fire)，比较充1 本课题得到国家自然科学基金项目(50438020；50578162)资助。- 1 -分地反映钢筋和混凝土的高温力学性能与构件的高温受力特点，实现钢筋混凝土柱在温度和荷载共同作用下的受力全过程分析，可望进

5、一步应用于局部火灾下钢筋混凝土足尺框架结构抗火分析中。2. 基本假设和材料热-力耦合本构关系2.1 基本假设火灾下钢筋混凝土柱受力分析采用两节点梁单元理论，节点截面采用分块计算模型，见图1。为简化分析，做如下假设：(1)平截面假定，即任一截面沿高度方向的轴向应变呈线性分布。(2)无滑移假定，即不考虑钢筋与核心混凝土之间粘结滑移的影响。(3)无剪切假定，即对于钢筋混凝土柱，其构件变形以弯曲变形为主，而忽略剪切变形的影响。(4)(5)蠕变仍是单调增长。(6)2.2 欧洲规范9在往复循环荷载作用下，一般假设常温下混凝土的应力-应变滞回关系骨架曲线与单向加载应力-应变曲线一致，高温下也采用这种假设。在

6、参考已有常温下混凝土单轴应力-应变滞回关系10的基础上，根据常温下混凝土单轴损伤本构模型11的卸载规则和高温下单轴受压应力-应变关系12，笔者构造高温下混凝土单轴应力-应变滞回关系，图3为其示意图。(1)包络线采用高温下混凝土单轴受压应力-应变关系式12：A1xx21+(A2)x x1 (1) 1y=x x121（）+xx1式中，y=/fcT、x=/cT，fcT为T下混凝土轴心抗压强度，cT为T下混凝土受压峰值应变， - 2 -cucT/c=1+0.23(T20)/1001.5,4/937/67/18, c=383fcu，1=2.5105fcu。 fc=0.4fcu106, A1=9.1fcu

7、高温下混凝土单轴受拉应力-应变关系式采用与常温下混凝土单轴受拉应力-应变关系式13一致的形式：A2xx21+(A2)x x1 (2) 2y=x x11.71xx（）+2式中y=/ftT，x=/tT，ftT为T下混凝土轴心抗拉强度，tT为T下混凝土受拉峰值应变，ftT1,tT/t=1+0.23(T20)/1001.5, =ft1+19(T20)/900b122/31/3, t=33fcu104。 ft=0.24fcu106, A2=1.306，2=1+3.4fcu需要指出，恒高温下混凝土轴心抗拉强度(ftT)研究很少，笔者暂用高温下轴心抗压强度折减系数来代替；由于未曾见到高温下混凝土受拉峰值应变

8、(tT)的研究报道，笔者暂用高温下受压峰值应变放大系数来代替。(2)受压卸载、再加载法则表1 不同温度下钢材屈服强度和弹性模量折减系数Es/ Es温度/fs/ fs Es/ Es温度/ fs/ fs1.0000.2300.1301.0000.0900.9000.06750.8000.0400.04500.7000.0200.02250.6000.0000.0000.310- 3 -当0.5cT时，按原曲线卸载、再加载；当0.5cT时，卸载、再加载为同一直线，其刚度为(1DcT)EcT，如从B点卸载至C点，C为自卸载点B卸载至C=0时的残余应变：B (3) C=BTT(1Dc)Ec式中DcT为高

9、温下混凝土受压区损伤变量，在文献11的基础上，引入不同温度下混凝土的受压峰值应变，其他参数不变，于是1(1x)c1D0 x1T (4) Dc=1D01c(1D)(x1)3+1x120式中x=/cT，D02.10.4ln(fcu+41)，c10.560.004fcu，c21.174.3410-5fcu2.8，c30.320.3ln(fcu10)。EcT为高温下混凝土弹性模量，计算如下12TEc1 (5) =b21.5Ec1+120(T20)/900+0.23(T20)/100式中7.65 20 fcu40MPa b2=3.054.60+1+0.001(f40)3 fcu40MPacu反向加载时，

10、当应变历史上出现的最大拉应变maxtT。即受拉混凝土尚未发生开裂时，则应力应变将沿直线CD发展，D(tT, ftT)为骨架线上峰值拉应力的对应点；当应变历史上出则应力应变将沿直线CE发展，E(E, E)为骨架线上最大拉应变现的最大拉应变maxtpT时，的对应点。(3)受拉卸载、再加载法则当0.5tT时，按原曲线卸载、再加载；当0.5tT时，卸载、再加载为同一直线，其刚度为(1DtT)EcT，如从E点卸载至F点，F为自卸载点E卸载至F=0时的残余应变：F=EE(1D)ETtTc (6)式中DtT为高温下混凝土受拉区损伤变量，在文献11的基础上，引入不同温度下混凝土的受拉峰值应变，其他参数不变，于

11、是1(1x)c1D0 x1 (7) DtT=1D01c(1D)(x1)c3+1x120式中x=/tT，D00.19，c10.31，c21.561.8310-4fcu2.08，c31.13.5410-3fcu。反向加载时，当应变历史上出现的最大压应变cT。即受压混凝土尚未发生开裂时，则应力应变将沿曲线FGA发展，A(cT, fcT)为骨架线上峰值压应力的对应点；当应变历史上出现的最大拉应变cT时，则应力应变将沿直线FB发展，B(B, B)为骨架线上最大压应变的对应点。- 4 -曲线FGA反映了裂面效应使刚度增大，由于裂面效应应力-应变关系将沿下列曲线到达受压峰值应力(cT)：A4xx2 (8)y

12、=1+(A42)xTT式中y=,x=F,A4=Ec(cF)。fcTfcTcTF当F0.3(cTF)时，曲线FGA上卸载时按原曲线卸载、再加载；当F0.3(cTF)时，卸载、再加载为同一直线，其刚度为(1Dc1T)EcT，从G点卸载至H点，H为自卸载点G卸载至H=0时的残余应变：H=GGT(1Dc1)EcT(9)式中Dc1T为高温下混凝土曲线FGA上的损伤变量，计算如下：TDc1=1(1x)c1D0 (10)式中x=(F)/(cTF)，D02.10.4ln(fcu+41)，c10.560.004fcu。当卸载至H点后反向加载，则沿直线HE返回骨架曲线。2.3 高温下混凝土的纵向热-力耦合本构模型

13、核心混凝土的纵向总应变(c)由4部分组成，即应力作用产生的应变(c,)、自由膨胀应变(c,th)、瞬态热应变(c,tr)和短期高温徐变(c,cr)。于是可建立cc,c,thc,trc,cr (11)则对于火灾时间从tf到tf+tf时c,cc,T(dc,thdc,trdc,cr) (12)式中，混凝土的应力以拉为正，压为负，应变以伸长为正，缩短为负，dc,th=cdT，dc,tr=c,tr(c,T)Tt=nt0dT，dc,cr=c,cr(c,T,tf)tfdtf，c,T为火灾时间从0tf时刻的dc,th、dc,tr和dc,cr的总和，即c,T=(dc,th+dc,tr+dc,cr)混凝土纵向热膨

14、胀系数c采用文献2建议的表达式c=(6+0.008T)106 (13)混凝土纵向受压时，在高温下混凝土单轴受压研究成果14基础上，针对文献2四面受火下矩形截面钢筋混凝土柱的变形性能试验结果，纵向瞬态热应变c,tr和纵向短期高温徐变c,cr的表达式修正为：c,tr=0.3fc0.17+0.73(T20T20103 (14) 100100c,cr=0.3fc(T20)1.25tf0.001106 (15)式中为混凝土应力，tf为火灾时间，单位为秒。而当混凝土在纵向受拉时，假设c,tr和c,cr- 5 -都为0。2.4 高温下钢材的纵向热-力耦合本构模型钢材在温度和应力共同作用下的总应变(s)由3部

15、分组成，即应力作用产生的应变(s,)、自由膨胀应变(s,th)和高温蠕变(s,cr)。于是可建立ss,s,ths,cr (16)则对于火灾时间从tf到tf+tf时s,ss,T(ds,thds,cr) (17)式中ds,th=sdT，ds,cr=s,cr(s,T,tf)tfdtf，s,T=t=nt0(ds,th+ds,cr)钢筋纵向热膨胀系数s计算式为15s=(11+0.0036T)106 (18)钢材单轴受拉、压时，按文献16建议，钢材的高温蠕变s,cr表达式为：s,cr=10a/(T+273)+b(/9.8)c/(T+273)+dtfe(T+273)+f (19)式中参数a8480，b2.5

16、0，c3060，d0.228，e0.002，f1.1；为钢材应力；tf为火灾时间，单位为min。3. 火灾下钢筋混凝土柱非线性方程组的求解火灾下钢筋混凝土梁单元的形心应变0和曲率的计算方法以及材料非线性的引入方法与常温的方法17一致，只是将常温的分层形式改为火灾下的分块形式。由火灾下增量格式 的U.L.列式虚功方程18，对于梁单元模型，火灾下梁单元的非线性方程组求解按下式进行：t+tt+tt(i1)(20) K0+K(i1)d(i)=t+tRtt+tF式中t+tt+tt(i1)t+t(i1)t+tt(i1)(21a) F(i1)=tt+t+tF(2i1)+tt+tK0+Ktd+tFt+tt+t

17、F(2i1)Mj+MiM+Mi,Mi,Nj,j,Mj (21b) =Ni,llTTt+tt+tmF(i1)nt+tM,j+M,iM+M,i,M,i,N,j,j,M,j (21c) =N,i,lltt+tN,i=(k=1l=1klklEklt+tkl,TAkl),M,i=ykl(t+tkltklt+tEklt+tkl,TAkl)k=1l=1mn(21d)tt+t式中tt+tK0为t+t位形的线性应变增量刚度矩阵，t+tK为非线性应变(几何或初始应力)t增量刚度矩阵，tt+tF为t+t位形单元应力引起的等效结点荷载向量，tF为t位形单元td为每次迭代的增量位移，温度、时间和应变不变温度变化引起的等

18、效结点荷载增量向量，t+t td为某荷载步下迭代后的增量总位移，右上标(i)为迭代次数，Ni、Nj为t位形加载终了- 6 -时单元杆端轴力，N,i、N,j为相应的t+t位形单元温度、时间和应变不变温度变化引起的等效杆端增量轴力，Mi、Mj为相应的杆端弯矩，M,i、M,j为相应的t+t位形的等效杆端增量弯矩，、E、和A分别为钢筋或混凝土块的纵向应力、纵向切线模量、纵向应变和面积。对于钢筋混凝土柱在常温加载时的求解方法，没有特殊之处，这里不再赘述。当实际钢筋混凝土柱在承受一定荷载下而发生火灾时(一般假设火灾环境温度随火灾时间的关系服从一定的火灾升温曲线，如ISO-834曲线4)，此时外荷载基本保持

19、不变，如对其进行非线性有限分析，当i1时，则按如下迭代：t+tt+tt(0) (22) K0+K(0)d(1)=tt+tF式中t+tF 是个未知量，不能直接计算，可按如下方法估算：在已知增量受火时间tf和增量温度T的情况下，计算得到温度、时间和应变不变温度变化引起的应变增量，假设每一块的应力和应力引起的应变不变，根据受火时间在tf+tf时刻的每一块材料应力-应变关系，计算得到该时刻该块材料的切线模量，按式(21c)计算即可。方程组求解采用Full Newton-Raphson迭代法，为确保求解精度，在迭代过程中每次迭代都更新等效结点荷载增量向量。根据本文思想编制了火灾下钢筋混凝土柱非线性有限元

20、程序NFEMRCLF。 t+t(0)4. 数值算例火灾下钢筋混凝土柱的试验研究由两部分组成，一为常温加载，模拟实际钢筋混凝土柱承受的正常使用荷载，二是恒载升温，荷载保持不变而钢筋混凝土柱周围环境温度不断升高，模拟钢筋混凝土柱在楼层发生火灾时的情形。文献2报道了3根两端固支的足尺钢筋混凝土轴压柱在恒载升温途径下的试验研究，试验结果合理，被国内外学者广泛应用。试验采用的火灾环境升温曲线按加拿大规程CAN-S1012进行：Tf=20+750(1e3.79553tf)+170.41tf0.5 (23)式中，Tf为火灾环境温度()，升温时间tf以小时计。一般认为，温度场的分析先独立于结构抗火计算，钢筋混

21、凝土柱在四面均匀受火下的温度场分布按文献19计算，图4给了钢筋混凝土柱截面温度场的计算结果与试验结果2的比较，两者符合较好，图中b为柱截面宽度，h为柱截面高度，d表示温度场测点沿短边离混凝土外表面的距离。结合本文提出的计算理论对其进行数值仿真分析，计算时取初始挠度v0L/1000，并假设挠曲线为正弦半波，L为柱的总长度。图5所示为典型的柱轴向变形(u)与火灾时间(tf)关系的计算曲线与试验曲线的比较，图中也给出了柱跨中侧向变形(v)与火灾时间(tf)关系的计算曲线(由于钢筋混凝土柱四面均匀受火，无法测得跨中侧向变形)，图中fs为钢筋屈服强度，fcu为混凝土立方体抗压强度，N0为柱子施加的轴向荷

22、载，C为柱混凝土保护层厚度。分析结果表明，该理论基本上可反映Lie T.T.等人2报道的试验结果的规律，数值计算的变形性能和耐火极限与试验结果相比基本合理，而计算的耐火极限基本上是试验结果的上限。分析表明，在火灾进程前期，构件截面温度场变化快，采用的时间步长较短(tf0.2分钟)；在火灾进程中期，构件截面温度场变化缓，采用的时间步长较长(tf0.52分钟)；在 - 7 -0.5火灾进程后期，此时的等效结点弯矩荷载增量较大，程序运算时较难收敛，采用的时间步长较短(tf0.2分钟)。1000800bh=305457mm900d=13mmd750bh=2036002=13mmdT/T/6004002

23、0004503001500=102mm计算曲线tf /mintf /min图4 钢筋混凝土柱温度场计算值与实测值的比较105105计算曲线-轴向u(v)/mm-5u(v)/mm00-5-102bhL=3053053810mm-15fs=414MPa fcu=44.1MPa425 C=40mmN=1067kN-2004080120160200240280-10tf /min2bhL=3053053810mm-15fs=414MPa fcu=44.1MPa425 C=40mmN=1067kN-200408012016020024028020100tf /min(a) 按应力-应变滞回关系计算 (b

24、) 按应力-应变骨架关系计算200u(v)/mm-20-40bhL=9142033810mmfs=414MPa fcu=50.1MPa 2819 C=40mmN=756kN-60060120180240300tf /min360u(v)/mm-10-20-30-40bhL=9142033810mmf=414MPa fcu=50.1MPa 2-50s819 C=40mmN=756kN-60060120180240300105tf /min360(a) 按应力-应变滞回关系计算 (b) 按应力-应变骨架关系计算105u(v)/mm-5u(v)/mm00-5-102bhL=4573053810mm

25、-15fs=414MPa fcu=50.5MPaC=40mm 822N=1413kN-20-10tf /min2bhL=4573053810mm -15fs=414MPa fcu=50.5MPaC=40mm 822N=1413kN-20tf /min(a) 按应力-应变滞回关系计算 (b) 按应力-应变骨架关系计算图5 钢筋混凝土柱计算曲线与试验曲线的比较图5中也给出了两种材料应力-应变关系对计算结果的影响，从中可以看出，采用应力-应变滞回关系(加载、卸载和再加载的路径不同)的变形计算结果比采用应力-应变骨架关系- 8 -(加载、卸载和再加载的路径相同)的变形计算结果更合理。此外，徐玉野5采用

26、ANSYS对图5中截面尺寸为457305mm的钢筋混凝土柱的抗火性能进行分析。计算时考虑结构的对称性，取1/4的柱子分析，混凝土采用Solid65单元，钢筋采用Link8单元，模型共分2716个单元，4606个节点，三维结构非线性分析耗时90h，可见计算量非常大。对于相同的算例，笔者编制的程序进行分析时，模型共分24个单元，25个节点，耗时在0.5h以内，可见采用梁单元的结构非线性分析，其计算速度是惊人的。5. 结论本文提出了高温下混凝土单轴应力-应变滞回关系，采用钢材和混凝土的热-力耦合本构模型，并基于火灾下结构的U.L.列式虚功增量方程，提出了较合理的火灾下钢筋混凝土柱双重非线性有限元分析

27、理论，编制了NFEMRCLF非线性有限元程序，计算结果与试验结果比较表明该理论可较好反映火灾下钢筋混凝土柱的变形性能和耐火极限性能，为建立火灾下钢筋混凝土柱的耐火极限计算方法和实现局部火灾下钢筋混凝土足尺框架结构抗火分析创造了良好的条件。参考文献1 过镇海,时旭东钢筋混凝土的高温性能及其计算M北京: 清华大学出版社,20032 Lie T.T., Irwin R.J. Method to calculate the fire resistance of reinforced concrete columns with rectangularcross sectionJ. ACI Structu

28、ral Journal, 1993, 90(1):52-60.3 Kodur V.K.R., Wang T.C., Cheng F.P.Predicting the fire resistance behaviour of high strength concretecolumnsJ. Cement and Concrete Composites, 2004,26(2):141-153.4 吴波,洪洲钢筋混凝土柱的耐火极限研究J自然灾害学报,2005,14(5):119-1265 徐玉野钢筋混凝土柱在高温下的数值模拟J建筑科学, 2005,21(6):41-436 陆洲导,朱伯龙,姚亚雄. 钢

29、筋混凝土框架火灾反应分析J. 土木工程学报,1995,28 (6):18-277 Zha X.X.Three-dimensional non-linear analysis of reinforced concrete members in fireJBuilding andEnviroment, 2003, 38(2):297-3078 Sidibe Karamoko, Duprat Frederic, Pinglot Michel, et al. Fire safety of reinforced concrete columns J. ACIStructural Journal, 200

30、0,97(4): 642-647.9 CEN(European Committee for Standardization), Eurocode 3Design of steel structures, draft part 1.4. fireresistanceS199210 滕智明,邹离湘. 反复荷载下钢筋混凝土构件的非线性有限元分析J. 土木工程学报,1996,29(2):19-17.11 丁发兴,余志武混凝土单轴损伤本构模型BD/OL中国科技论文在线( 2006-6-12.12 丁发兴,余志武恒高温下混凝土及钢管混凝土受压力学性能研究J铁道科学与工程学报,2005,2(6):9-14.

31、13 丁发兴,余志武混凝土受拉力学性能统一计算方法J华中科技大学学报(城市科学版),2004,21(3):293414 南建林温度-应力耦合作用下混凝土力学性能的试验研究D北京：清华大学,199415 ACI 216 R-94Guide for determining the fire endurance of concrete elementsS199416 孙金香,高伟(译)建筑物综合防火设计M天津: 天津科技翻译出版公司,199217 丁发兴,余志武,蒋丽忠圆钢管混凝土结构非线性有限元分析BD/OL中国科技论文在线(http:/ 2006-6-8.18 丁发兴,余志武,蒋丽忠火灾下圆钢管

32、混凝土柱非线性有限元分析BD/OL中国科技论文在线(http:/ 2006-6-14.19 余志武,唐国庆,丁发兴三面受火下钢筋混凝土梁温度场非线性分析J建筑科学与工程学报,2005,22(4):11-14.- 9 -Nonlinear Finite Element Analysis of Reinforced ConcreteColumns in FireDING Fa-xing, YU Zhi-wuSchool of Civil Engineering and Architecture, Central South University, Changsha 410075,ChinaAbst

33、ractAn analysis theory of nonlinear finite beam element for reinforced concrete (RC) columns in fire was presented and a FORTRAN program named NFEMRCLF was developed in this paper. First, the stress-strain hysteretic relations for concrete at high temperature was proposed and appropriate thermal-stress coupling constitutive model of both steel bars and concrete in fire were presented. Whats more, based on the U.L. formula of element increm