非线性大作业

非线性大作业-标准化文件发布号：（9456-EUATWK-MWUB-WUNN-INNUL-DDQTY-KII工程结构非线性分析学院：姓名：学号：指导教师：目录1、预应力混凝土梁截面非线性 01.1 材料的本构关系 01.2 平截面假定 (2)1.3 预应力筋作用下截面初应变的求解 (3)2、预应力混凝土梁构件的非线性 (4)2.1 构件弯曲的一般理论 (4)2.2 共轭梁分析法 (4)2.3 预应力钢筋混凝土梁非线性分析的数值法 (5)3、算例分析 (7)3.1 试验梁简介 (7)3.2 截面非线性与构件非线性分析程序编制 (8)3.3 试验结果验证 (9)3.4 结果分析 (11)参考文献 (12)附录 (13)作业2：预应力混凝土梁的非线性全过程分析要求：1.阐述预应力混凝土梁截面和构件非线性全过程分析的理论背景；2.编制相应的截面和构件非线性分析程序，给出具体算例分析结果，方法及程序的适用性必须有试验结果的验证。

1、预应力混凝土梁截面非线性1.1 材料的本构关系 1.1.1 混凝土本构关系混凝土受压采用Rush 建议的应力—应变曲线，如图1-1所示。

0cuf cσ图1-1 混凝土受压应力-应变曲线000[1(1)]n c cc c cc cu f f εεεεεεε⎧-- 0≤≤⎪σ=⎨⎪ <≤⎩ 式中c σ——对应于混凝土应变为c ε时的混凝土压应力；c f ——混凝土抗压强度标准值；cu ε——正截面处于非均匀受压时的混凝土极限压应变，50.0033(50)10cu cu f ε-=--⨯，当0.0033cu ε>时，取为0.0033；0ε——受压峰值应变，500.0020.5(50)10cu f ε-=+-⨯，当00.002ε<时，取为0.002；n ——系数，,12(50)60cu k n f =--，当 2.0n >时，取为2.0。

为计算方便，混凝土受拉应力-应变曲线采用线性式，如图1-2所示。

γtσtuf图1-2 混凝土受拉应力-应变曲线t t tt t tE f f εγγ σ≤⎧σ=⎨0 σ>⎩式中t σ—对应于混凝土应变为ε时的混凝土拉应力；t f —混凝土轴心抗拉强度值；对于高强混凝土，取2/30.21t cu f f =；γ—混凝土受拉塑性系数，取为γ=1.5。

t E —混凝土受拉弹性模量，计算式如下：42(1.450.628)10/t t E f N mm =+⨯1.1.2 普通钢筋的本构关系认为钢筋承受拉力和承受压力的本构关系类同，如图1-3所示。

ss,maxf yσ0图1-3 普通钢筋应力-应变曲线0,max s s s s s ys s s E f εεεεεε 0≤≤⎧⎪σ=⎨<≤⎪⎩ 式中y f ——钢筋的强度设计值；0s ε——对应于钢筋应力为y f 时的钢筋应变；s E ——对应于钢筋弹性阶段的弹性模量；,max s ε——钢筋极限拉应变。

1.1.3 预应力钢筋的本构关系预应力钢筋应力-应变关系取理想模型，如图1-4所示。

pp,maxf yσ0图1-4 预应力钢筋应力-应变曲线0p,max p p p p p yp p E f εεεεεε 0≤≤⎧⎪σ=⎨ <≤⎪⎩式中y f －为预应力钢筋的强度设计值；p E －对应于钢筋弹性阶段的弹性模量。

1.2 平截面假定p A 总量（应变量）不在平截面上，增量在平截面上。

在任一受力P 下，pε与c ε、s ε'、s ε不成线性关系，p ε∆与c ε∆、sε'∆、s ε∆成线性关系，即p p p εεε∆+=0。

1.3 预应力筋作用下截面初应变的求解（1）将截面划分为条带；（2）将预应力筋引起的偏心压力作为外力加在截面上，分解为轴心压力p N 和预弯矩p M 作用在截面上；（3）根据轴向力平衡二分法迭代求解在轴向力单独作用下的初应变（如图1-5-b ）；（4）在预弯矩单独作用下，假定初始曲率，受拉区高度，求得弯矩作用下的初应变分布，再加上轴向力引起初应变，得到总的截面初应变分布，由轴向力平衡方程（11csN N p ci ci sk sk i k N A A σσ===∆+∆∑∑）二分法迭代求解初曲率对应的受拉区高度；（5）由弯矩平衡方程（11csN N p i ci ci k sk sk i k M y A y A σσ===∆+∆∑∑）二分法迭代求解初曲率，直到两个平衡方程均满足为止，求得初曲率和其纯弯时对应的受拉区高度；（6）由第（5）步求得和，算得预弯矩p M 作用下的各条带应变分布（如图1-5-a ）；（7）将（3）、（6）步求得的应变叠加起来得到在预应力筋作用下的截面各条带的初应变分布。

abcdεmφy*εm εnεm φy*εn图1-5 截面初应变几何示意2、预应力混凝土梁构件的非线性2.1 构件弯曲的一般理论截面曲率： =d dsθϕ 对于初等梁有： ''''3/2'2==1()d ds θνϕνν≈⎡⎤+⎣⎦即有： ''=ϕν 则有：x 处的截面转角： '0==xdx θνϕ⎰x 处的位移： 0=xx dx νϕ⎰因此当构件截面曲率的大小及分布确定以后，则构件任一截面的转角和位移即可确定。

上述分析公式与材料特性无关，故对线性材料和非线性材料组成的构件均可适用。

2.2 共轭梁分析法对于初等梁有：22dV q dx d M q dx ⎧=⎪⎪⎨⎪=⎪⎩ ⇒ 22d dxd v dx θϕϕ⎧=⎪⎪⎨⎪=⎪⎩根据上述微分方程的相似性，如果将构件每一截面的曲率ϕ视为沿构件分布作用的虚拟“荷载”，则每一截面相应的“剪力”和“弯矩”即为该截面的转角和位移。

但应注意内力和位移所对应的边界条件不同。

q实际梁共轭梁qφφ图2-1 实际梁与共轭梁转换关系图2.3 预应力钢筋混凝土梁非线性分析的数值法 2.3.1 基本假定① 截面的弯矩－曲率关系为数值型； ② 构件的剪切变形忽略不计。

2.3.2 基本思路bs's图2-2 截面网格划分及应变分布① 将截面划分成c N 个混凝土条带、s N 个钢筋条带，求出每一截面离散型的弯矩－曲率关系；② 将梁划分为m 段，借以刻划截面刚度沿梁长的分布；③ 根据荷载P 作用下的截面弯矩和截面的弯矩－曲率关系，确定此时每一截面的曲率；④ 根据每一截面的曲率积分或采用共轭梁法求制定截面的转角和挠度。

2.3.3 分析步骤分级加变形－分级加控制截面的曲率（1）分级施加控制截面的曲率ϕ：ϕϕϕ∆+=-1i i （2）由i ϕ根据截面的ϕ-M 关系求相应的截面弯矩i M ； （3）由i M 计算相应的i P ；（4）由i P 计算其他截面的弯矩mi M ； （5）由mi M 计算其他截面的曲率mi ϕ； （6）由mi ϕ计算挠度i f ；（7）满足构件的破坏条件？yes ⇒结束；No ⇒增加ϕ继续计算。

分级加变形也可采用分级加控制截面的挠度或分级加受压混凝土极限压应变。

分级加荷载（1）1i i P P P -=+∆（2）由i P 计算相应的mi M ； （3）由mi M 计算截面的曲率mi ϕ； （4）由mi ϕ计算挠度i f ； （5）满足构件的破坏条件？yes ⇒结束；No ⇒增加P 继续计算。

3、算例分析（文献1中PPB07试验梁）3.1 试验梁简介试验梁为矩形截面简支梁，构件基本尺寸为B×H×L=250×400×5700mm ，其中支点距梁端150mm ，计算跨径mm l 54000=。

截面下缘普通受拉钢筋为3根Φ12的HRB335钢筋，钢筋距梁底40mm ，横向间距85mm ；截面上缘配有2根Φ10的HRB335架立钢筋；箍筋采用Φ10的HRB335钢筋，箍筋间距为200mm 。

试验梁配置有2根)71(2.15⨯s φ高强低松驰钢铰线，钢铰线距梁底80mm ，钢饺线之间横向间距采用100mm 。

试验梁尺寸如图3-1所示，梁截面如图3-2所示。

图3-1 试验梁结构示意图（单位：cm ）图3-2 试验梁截面尺寸（单位：cm ，钢筋尺寸采用mm ）试验总体布置采用简支梁3分点加载方式。

在试验梁两端分别设置专用固定铰支座与滚动铰支座，支座中心距梁端15cm 。

在试验梁上端3/0l 处分别精确设置专用固定铰支座与活动铰支座，其上设置一分配梁。

加载受力示意如图3-3所示。

图3-3 加载示意简图（单位：cm ）3.2 截面非线性与构件非线性分析程序编制使用Matlab 计算，通过控制跨中截面梁顶应变分级加载。

源程序见附录。

该程序适用于有粘结直线预应力钢筋混凝土梁。

3.3 试验结果验证程序计算的截面非线性曲线ϕ-M 关系如图3-4所示，构件非线性曲线f P -关系如图3-5所示，f P -关系的计算值与试验值的对比结果如图3-6所示。

图3-4 截面非线性曲线ϕ-MP 图3-5 构件非线性曲线f图3-6 计算值与试验值的对比3.4 结果分析曲线基本呈现四折线，体现了预应力混凝土梁截面破坏的四个阶段，第一个折点为外加弯矩克服预压弯矩后，下缘混凝土达到极限拉应变后，受拉区混凝土开裂点；第二个折点表示在外加荷载增加到一定程度后，受拉钢筋达到屈服点；第三个折点表示随着荷载进一步加大，预应力钢筋达到屈服点，梁的抗弯刚度趋于零。

计算结果的精确程度主要取决于计算所采用的本构关系的准确性，本程序误差来源主要体现在以下三个方面：①混凝土、普通钢筋、预应力钢筋本构关系的选取与实际结构间必然存在一定偏差，这是影响计算精度的主要因素；②本文认为混凝土与预应力筋粘结良好，不考虑两者之间的滑移，这样的处理在与预应力筋相邻的混凝土开裂之前是符合实际情况的，在混凝土开裂以后会存在一定的误差；③在程序中，混凝土和钢筋的强度取的都是标准值，而非设计值。

混凝土受拉的塑性放大系数取为2.0，这些都将影响到弯矩曲率曲线，如果取值偏大，将使得全梁的极限承载力计算值略大于试验值，而且刚度也会变大，这也就使得与文献中试验梁数值有一定的误差。

总体来说，本程序得到的截面ϕ-P-关系曲线是一个四折线，与M和梁f试验梁数值存在一定误差，但基本符合，这说明了本程序的正确性，故可以利用本程序对有粘结直线预应力钢筋混凝土梁进行非线性分析。

参考文献[1]黄文雄.基于新型弯起器的折线配筋先张梁力学性能研究[D].武汉：华中科技大学，2012[2]唐国金.预应力混凝土梁的非线性有限元分析与参数评估[D].长沙：国防科学技术大学，2006附录1.混凝土应力-应变关系子函数[conc.m]function y=conc(s) %混凝土应力-应变关系（开裂点前用）fc=29.52;s0=0.002;su=0.0038;ft=-2.51;Ec=3.36*10^4;if s>=0if s<=s0y=fc*(1-(1-s/s0)^2);elseif s<=suy=fc;elsey=0;endelseif Ec*s>=fty=Ec*s;elsey=0;endend2.普通钢筋应力-应变关系子函数[reinf.m]function y=reinf(s) %普通钢筋应力-应变关系Es=2*10^5; %MPafy=335;if abs(s)<=fy/Esy=Es*s;elsey=s/abs(s)*fy;end3.预应力钢筋应力-应变关系子函数[pres.m]function y=pres(s) %预应力钢筋应力-应变关系Es=1.977*10^5; %MPafy=1910.5;if abs(s)<=fy/Esy=Es*s;elsey=s/abs(s)*fy;end4.截面非线性初始点求解子函数[find_ini.m]function 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:nstrainc(k)=str_t-(k-1/2)*hs*phi; %混凝土第 k 层应变（压为正）曲率以下缘受拉为正stressc(k)=conc(strainc(k));c=c+stressc(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*h1;m=mcc+m1+m2-mp; %平衡条件2：合力矩y(1)=c;y(2)=m;5.截面弯矩计算子函数[find_m.m]function 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:nstrainc(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)*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+m2+m3;y=m;6.截面曲率求解子函数[find_phi.m]function 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:nstrainc(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; %上部钢筋应变（压为正）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.m]clcclearglobal h b as h0 As1 As2 As3 h1 ps sigma Es Np 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*10^5;Np=As3*sigma;%初始点确定，初始曲率,预应力筋所在处混凝土应变x0=fsolve(@(x)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));endfigure;plot(phi,m/1000000);grid;xlabel('{\it\phi} (m^{-1})');ylabel('{\itM} (kN\cdotm)'); saveas(gcf,'M_phi关系图','emf');%求构件P_f关系，由对称性取半跨进行结构计算dL=10; %两端划分密度（每段10mm）共540段。