www.gusucode.com > matlab非线性混凝土模型计算分析源码程序 > nonlinear.m

    function [yu,d,v,a,ep]=nonlinear(para,sa,np,per,xg,t,dt,eb)
% 能量分析则分别有ei、ek、ed、ep、ee；
st=length(xg);
m=1;%先假定为单位质量；
ke=m*(2*pi/per)^2;%初始弹性刚度，计算时频反应谱是一个系列；
%此为恢复力模型输入参数；
eta=0.05;
afa=0.5;bta=0.5;p=0.05;%%钢混结构，假定uc/uy=afa;kc/ke=bta；
ebuxi=eb;
fy=ebuxi*sa*m;%一般工程结构,双线性屈服强度系数ebuxi为0.25~0.4左右，三线型适当降低；
kc=ke*bta;ky=kc*p;%初始刚度为ke，开裂刚度为kc，屈服后刚度为ky；
uc=fy/(ke+(1/afa-1)*kc);uy=uc/afa;

f1=0;
kk1=ke;%迭代前的刚度为ke，此处将初始弹性刚度赋予它；
txg=zeros(1,2);
x1=zeros(1,3);
x2=zeros(1,3);   
x1(3)=-xg(1);
up=uy;fp=fy;%%%%初始化；
un=-uy;fn=-fy;
ucom=uc;
pd=0;
for nt=1:st    
	if (nt==1)%初始时刻位移和速度为零，加速度为负初始时刻地震加速度；
        d(1)=0;
        v(1)=0;
        a(1)=-xg(1);
        ff(1)=0;        
        ep=0;%能量分析值；
	else
		txg(1)=xg(nt-1);
		txg(2)=xg(nt);
        [x2]=newmark_single(dt,kk1,m,eta,txg,x1);
        [kk2,x2,f2,pd,ucom,up,fp,un,fn]=concrete_trilinear(pd,per,nt,m,eta,kk1,ke,kc,ky,uc,uy,dt,txg,x1,x2,f1,ucom,up,fp,un,fn);
        %能量分析值，采用近似梯形公式，最好是考虑拐点分割为精确；
        if (f1*f2>=0)
            ep=ep+(f2+f1)*(x2(1)-x1(1));  
        else
            ep=ep-f1^2/kk1+f2^2/kk2;
        end
        kk(nt)=kk2;
	    kk1=kk2;
	    f1=f2;
	    x1=x2;
		d(nt)=x2(1);
		v(nt)=x2(2);
		a(nt)=x2(3);
    end
    ff(nt)=f1;
end
yu=abs(d)/uy;
ee=f2^2/kk2/2;
ep=ep-ee;

% if (per==0.2)
%     for ni=1:500
% 	figure(31)
% 	scatter(d(ni),ff(ni));
% 	hold on;
% 	end
%     for ni=501:1000
% 	figure(32)
% 	scatter(d(ni),ff(ni));
% 	hold on;
% 	end
%     for ni=1001:st
% 	figure(33)
% 	scatter(d(ni),ff(ni));
% 	hold on;
% 	end
% end
% if (per==0.2)
%     for ni=1:300
% 	figure(41)
% 	scatter(d(ni),ff(ni));
% 	hold on;
% 	end
%     for ni=301:400
% 	figure(42)
% 	scatter(d(ni),ff(ni));
% 	hold on;
% 	end
%     for ni=401:500
% 	figure(43)
% 	scatter(d(ni),ff(ni));
% 	hold on;
% 	end
% end