1、工程结构非线性分析学 院：姓 名：学 号：指导教师：目录1、预应力混凝土梁截面非线性11.1 材料的本构关系11.2 平截面假定31.3 预应力筋作用下截面初应变的求解42、预应力混凝土梁构件的非线性52.1 构件弯曲的一般理论52.2 共轭梁分析法52.3 预应力钢筋混凝土梁非线性分析的数值法63、算例分析83.1 试验梁简介83.2 截面非线性与构件非线性分析程序编制93.3 试验结果验证93.4 结果分析11参考文献12附录13作业2：预应力混凝土梁的非线性全过程分析要求：1. 阐述预应力混凝土梁截面和构件非线性全过程分析的理论背景；2.编制相应的截面和构件非线性分析程序，给出具体算例分
2、析结果，方法及程序的适用性必须有试验结果的验证。1、预应力混凝土梁截面非线性1.1 材料的本构关系1.1.1 混凝土本构关系混凝土受压采用Rush建议的应力应变曲线，如图1-1所示。图1-1 混凝土受压应力-应变曲线式中对应于混凝土应变为时的混凝土压应力； 混凝土抗压强度标准值；正截面处于非均匀受压时的混凝土极限压应变，当时，取为0.0033； 受压峰值应变，当时，取为0.002； 系数，当时，取为2.0。为计算方便，混凝土受拉应力-应变曲线采用线性式，如图1-2所示。图1-2 混凝土受拉应力-应变曲线式中对应于混凝土应变为时的混凝土拉应力；混凝土轴心抗拉强度值；对于高强混凝土，取； 混凝土受
3、拉塑性系数，取为=1.5。 混凝土受拉弹性模量，计算式如下：1.1.2 普通钢筋的本构关系认为钢筋承受拉力和承受压力的本构关系类同，如图1-3所示。图1-3 普通钢筋应力-应变曲线式中钢筋的强度设计值；对应于钢筋应力为时的钢筋应变；对应于钢筋弹性阶段的弹性模量； 钢筋极限拉应变。1.1.3 预应力钢筋的本构关系预应力钢筋应力-应变关系取理想模型，如图1-4所示。图1-4 预应力钢筋应力-应变曲线式中为预应力钢筋的强度设计值；对应于钢筋弹性阶段的弹性模量。1.2 平截面假定总量（应变量）不在平截面上，增量在平截面上。在任一受力下，与、不成线性关系，与、成线性关系，即。1.3 预应力筋作用下截面初
4、应变的求解（1）将截面划分为条带；（2）将预应力筋引起的偏心压力作为外力加在截面上，分解为轴心压力和预弯矩作用在截面上；（3）根据轴向力平衡二分法迭代求解在轴向力单独作用下的初应变（如图1-5-b）； （4）在预弯矩单独作用下，假定初始曲率，受拉区高度，求得弯矩作用下的初应变分布，再加上轴向力引起初应变，得到总的截面初应变分布，由轴向力平衡方程（）二分法迭代求解初曲率对应的受拉区高度；（5）由弯矩平衡方程（）二分法迭代求解初曲率，直到两个平衡方程均满足为止，求得初曲率和其纯弯时对应的受拉区高度；（6）由第（5）步求得和，算得预弯矩作用下的各条带应变分布（如图1-5-a）；（7）将（3）、（6）
5、步求得的应变叠加起来得到在预应力筋作用下的截面各条带的初应变分布。图1-5 截面初应变几何示意2、预应力混凝土梁构件的非线性2.1 构件弯曲的一般理论截面曲率： 对于初等梁有： 即有： 则有：处的截面转角： 处的位移： 因此当构件截面曲率的大小及分布确定以后，则构件任一截面的转角和位移即可确定。上述分析公式与材料特性无关，故对线性材料和非线性材料组成的构件均可适用。2.2 共轭梁分析法对于初等梁有： 根据上述微分方程的相似性，如果将构件每一截面的曲率视为沿构件分布作用的虚拟“荷载”，则每一截面相应的“剪力”和“弯矩”即为该截面的转角和位移。但应注意内力和位移所对应的边界条件不同。图2-1 实际
6、梁与共轭梁转换关系图2.3 预应力钢筋混凝土梁非线性分析的数值法2.3.1 基本假定 截面的弯矩曲率关系为数值型； 构件的剪切变形忽略不计。2.3.2 基本思路图2-2 截面网格划分及应变分布 将截面划分成个混凝土条带、个钢筋条带，求出每一截面离散型的弯矩曲率关系； 将梁划分为m段，借以刻划截面刚度沿梁长的分布； 根据荷载P作用下的截面弯矩和截面的弯矩曲率关系，确定此时每一截面的曲率； 根据每一截面的曲率积分或采用共轭梁法求制定截面的转角和挠度。2.3.3 分析步骤l 分级加变形分级加控制截面的曲率（1）分级施加控制截面的曲率：（2）由根据截面的关系求相应的截面弯矩；（3）由计算相应的；（4）
7、由计算其他截面的弯矩；（5）由计算其他截面的曲率；（6）由计算挠度；（7）满足构件的破坏条件？yes结束；No增加继续计算。分级加变形也可采用分级加控制截面的挠度或分级加受压混凝土极限压应变。l 分级加荷载（1）（2）由计算相应的；（3）由计算截面的曲率；（4）由计算挠度；（5）满足构件的破坏条件？yes结束；No增加继续计算。3、算例分析（文献1中PPB07试验梁）3.1 试验梁简介试验梁为矩形截面简支梁，构件基本尺寸为BHL=2504005700mm，其中支点距梁端150mm，计算跨径。截面下缘普通受拉钢筋为3根12的HRB335钢筋，钢筋距梁底40mm，横向间距85mm；截面上缘配有2根
8、10的HRB335架立钢筋；箍筋采用10的HRB335钢筋，箍筋间距为200mm。试验梁配置有2根高强低松驰钢铰线，钢铰线距梁底80mm，钢饺线之间横向间距采用100mm。试验梁尺寸如图3-1所示，梁截面如图3-2所示。图3-1 试验梁结构示意图（单位：cm）图3-2 试验梁截面尺寸（单位：cm，钢筋尺寸采用mm）试验总体布置采用简支梁3分点加载方式。在试验梁两端分别设置专用固定铰支座与滚动铰支座，支座中心距梁端15cm。在试验梁上端处分别精确设置专用固定铰支座与活动铰支座，其上设置一分配梁。加载受力示意如图3-3所示。图3-3 加载示意简图（单位：cm）3.2 截面非线性与构件非线性分析程序
9、编制使用Matlab计算，通过控制跨中截面梁顶应变分级加载。源程序见附录。该程序适用于有粘结直线预应力钢筋混凝土梁。3.3 试验结果验证程序计算的截面非线性曲线关系如图3-4所示，构件非线性曲线关系如图3-5所示，关系的计算值与试验值的对比结果如图3-6所示。图3-4 截面非线性曲线图3-5 构件非线性曲线图3-6 计算值与试验值的对比3.4 结果分析曲线基本呈现四折线，体现了预应力混凝土梁截面破坏的四个阶段，第一个折点为外加弯矩克服预压弯矩后，下缘混凝土达到极限拉应变后，受拉区混凝土开裂点；第二个折点表示在外加荷载增加到一定程度后，受拉钢筋达到屈服点；第三个折点表示随着荷载进一步加大，预应力
10、钢筋达到屈服点，梁的抗弯刚度趋于零。计算结果的精确程度主要取决于计算所采用的本构关系的准确性，本程序误差来源主要体现在以下三个方面：混凝土、普通钢筋、预应力钢筋本构关系的选取与实际结构间必然存在一定偏差，这是影响计算精度的主要因素；本文认为混凝土与预应力筋粘结良好，不考虑两者之间的滑移，这样的处理在与预应力筋相邻的混凝土开裂之前是符合实际情况的，在混凝土开裂以后会存在一定的误差；在程序中，混凝土和钢筋的强度取的都是标准值，而非设计值。混凝土受拉的塑性放大系数取为2.0，这些都将影响到弯矩曲率曲线，如果取值偏大，将使得全梁的极限承载力计算值略大于试验值，而且刚度也会变大，这也就使得与文献中试验梁
11、数值有一定的误差。总体来说，本程序得到的截面和梁关系曲线是一个四折线，与试验梁数值存在一定误差，但基本符合，这说明了本程序的正确性，故可以利用本程序对有粘结直线预应力钢筋混凝土梁进行非线性分析。参考文献1黄文雄.基于新型弯起器的折线配筋先张梁力学性能研究D.武汉：华中科技大学，20122唐国金.预应力混凝土梁的非线性有限元分析与参数评估D.长沙：国防科学技术大学，2006附录1.混凝土应力-应变关系子函数conc.mfunction y=conc(s) %混凝土应力-应变关系（开裂点前用）fc=29.52;s0=0.002;su=0.0038;ft=-2.51;Ec=3.36*104; if
12、s=0 if s=s0 y=fc*(1-(1-s/s0)2); elseif s=ft y=Ec*s;else y=0;endend2.普通钢筋应力-应变关系子函数reinf.mfunction y=reinf(s) %普通钢筋应力-应变关系Es=2*105; %MPafy=335;if abs(s)=fy/Es y=Es*s;else y=s/abs(s)*fy;end3.预应力钢筋应力-应变关系子函数pres.mfunction y=pres(s) %预应力钢筋应力-应变关系Es=1.977*105; %MPafy=1910.5;if abs(s)=fy/Es y=Es*s;else y=
13、s/abs(s)*fy;end4.截面非线性初始点求解子函数find_ini.mfunction y=find_ini(x) %初始点求解器str_t=x(1);phi=x(2);global h b as h0 As1 As2 As3 h1 ps sigma Es Np str_t0 phi0n=1000; %混凝土截面划分为n层hs=h/n; %混凝土分层厚度c=0;mcc=0;for k=1:n strainc(k)=str_t-(k-1/2)*hs*phi; %混凝土第 k 层应变（压为正）曲率以下缘受拉为正 stressc(k)=conc(strainc(k); c=c+stress
14、c(k)*hs*b; %混凝土合力计算（以受压为正） mcc=mcc+stressc(k)*hs*b*(k-1/2)*hs; %对顶部取矩endstrain1=str_t-phi*as; %上部钢筋应变（压为正）stress1=reinf(strain1); %上部钢筋应力strain2=str_t-phi*h0; %下部钢筋应变（压为正）stress2=reinf(strain2); %下部钢筋应力c=c+stress1*As1+stress2*As2-Np; %平衡条件1：合力（以受压为正）m1=stress1*As1*as; %对顶部取矩m2=stress2*As2*h0;mp=Np*h
15、1;m=mcc+m1+m2-mp; %平衡条件2：合力矩y(1)=c;y(2)=m;5.截面弯矩计算子函数find_m.mfunction y=find_m(s,phi)global h b as h0 As1 As2 As3 h1 ps sigma Es Np str_t0 phi0n=1000;hs=h/n;c=0;mc=0;for k=1:n strainc(k)=-(k-1/2)*hs*phi+s; %混凝土第 k 层应变（压为正） stressc(k)=conc(strainc(k); c=c+stressc(k)*hs*b; mc=mc-stressc(k)*hs*b*(k-1/2
16、)*hs; %对顶部取矩endstrain1=-as*phi+s; %上部钢筋应变（压为正）stress1=reinf(strain1); %上部钢筋应力strain2=s-h0*phi; %下部钢筋应变（压为正）stress2=reinf(strain2); %下部钢筋应力strain3=s-h1*phi-(str_t0-h1*phi0)-sigma/Es; %预应力钢筋应变（压为正）stress3=pres(strain3); %预应力钢筋应力m1=-stress1*As1*as; %对顶部取矩m2=-stress2*As2*h0; m3=-stress3*As3*h1;m=+mc+m1+
17、m2+m3;y=m;6.截面曲率求解子函数find_phi.mfunction y=find_phi(s,a)global h b as h0 As1 As2 As3 h1 ps sigma Es Np str_t0 phi0n=1000;hs=h/n;c=0;mc=0;for k=1:n strainc(k)=-(k-1/2)*hs*a+s; %混凝土第 k 层应变（压为正） stressc(k)=conc(strainc(k); c=c+stressc(k)*hs*b; mc=mc+stressc(k)*hs*b*(k-1/2)*hs; %对顶部取矩endstrain1=-as*a+s;
18、%上部钢筋应变（压为正）stress1=reinf(strain1); %上部钢筋应力strain2=s-h0*a; %下部钢筋应变（压为正）stress2=reinf(strain2); %下部钢筋应力strain3=s-h1*a-(str_t0-h1*phi0)-sigma/Es; %预应力钢筋应变（压为正）stress3=pres(strain3); %预应力钢筋应力c=c+stress1*As1+stress2*As2+stress3*As3; %合力y=c;7.主程序main.mclcclearglobal h b as h0 As1 As2 As3 h1 ps sigma Es N
19、p str_t0 phi0ntd=1000; %截面条带划分数1000options=optimset(Algorithm,levenberg-marquardt,0.005,ScaleProblem,Jacobian,Display,off,TolFun,1e-10,TolX,1e-16); %迭代选项h=400;b=250;as=40;h0=h-as;ps=80;h1=h-ps;L=5400;As1=157;As2=339;As3=278;sigma=1019.73;Es=1.977*105;Np=As3*sigma;%初始点确定，初始曲率,预应力筋所在处混凝土应变x0=fsolve(x)
20、find_ini(x),0,0,options);str_t0=x0(1);phi0=x0(2); %初始曲率和初始上缘应变（以压为正）str_pc0=str_t0-phi0*h1;%求截面M_phi关系采用分级加上缘应变str_t=str_t0:5e-5:0.0038; %上缘应变列表m(1)=0; %初始弯矩phi(1)=phi0; %初始曲率for i=2:length(str_t) phi(i)=fsolve(a)find_phi(str_t(i),a),phi(i-1),options); %上缘应变对应的弯矩 m(i)=find_m(str_t(i),phi(i);endfigur
21、e;plot(phi,m/1000000);grid;xlabel(itphi (m-1);ylabel(itM (kNcdotm);saveas(gcf,M_phi关系图,emf);%求构件P_f关系，由对称性取半跨进行结构计算dL=10; %两端划分密度（每段10mm）共540段。xx=(L/2:-dL:0); %节点距一段支座坐标，保证第一个元素为跨中P=m/(L/3)*2; %荷载列表mm=zeros(length(str_t),length(xx);ss=zeros(length(str_t),length(xx);aa=zeros(length(str_t),length(xx);
22、f=zeros(length(str_t),1);for i=1:length(str_t) %跨中混凝土上缘应变作为控制变量 for j=1:length(xx) %求节点弯矩 if xx(j)L/3 mm(i,j)=3*m(i)/L*xx(j); else mm(i,j)=m(i); end end aa(i,1)=phi(i);ss(i,1)=str_t(i); for j=2:length(xx) aa(i,j)=interp1(m,phi,mm(i,j); %插值求得每个节点的曲率 f(i)=f(i)+(aa(i,j)+aa(i,j-1)/2*dL*(xx(j)+xx(j-1)/2;
23、 %共轭梁法求跨中挠度 endendfigure;plot(f,P/1000);grid;xlabel(itf (mm);ylabel(itP (kN);saveas(gcf,P_f关系图,emf);P_exp=0 72.3 182.6 239.12 242.68;f_exp=-3.85 3.92 26.96 78.05 96.70;figure;plot(f,P/1000);hold all;plot(f_exp,P_exp,-,Marker,*);hold off;grid;xlabel(itf (mm);ylabel(itP (kN);legend(计算值,试验值,Location,SouthEast);saveas(gcf,P_f对比图,emf);15
