guaidian2.m

钢、混凝土结构时频反应谱程序

function [tm,xm]=guaidian2(dt,k1,m,eta,xg,x1);%第二类拐点处理； 
%%%%精确拐点非迭代法；纽马克法求解；
tm=-1;
%采用精确拐点非迭代法求出零速度拐点；
%gama/beta等于2则解二次方程，很方便的。
gama=0.5;beta=0.25;w=sqrt(k1/m);cc=2*eta*w*m;
pa=k1*x1(3)/2-k1*beta*x1(3)/gama;%%pa==0
pb=(xg(2)-xg(1))*m/dt+k1*x1(2)-k1*beta*x1(2)/gama;
pc=-m*x1(3)/gama-cc*x1(2);
pd=-m*x1(2)/gama;
if(pb==0)
    tm=-pd/pc;
else
    a=pc/pb;b=pd/pb;%%%%%%x=solve('x^2+a*x+b=0');
    det(1)=-1/2*a+1/2*(a^2-4*b)^(1/2);
    det(2)=-1/2*a-1/2*(a^2-4*b)^(1/2);
	for ni=1:2
        if det(ni)<dt&det(ni)>=0
            tm=det(ni);
        end
	end
end
if(length(tm)~=1)
    tm=min(tm);
end
if(tm>dt|tm<0)
    cuowu='222222?????'
end
txg(1)=xg(1);
txg(2)=xg(1)+(xg(2)-xg(1))*tm/dt;
if(tm==0)
    xm=x1;
else    
    [xm]=newmark_single(tm,k1,m,eta,txg,x1);
end
%%%不计误差则应该xm(2)=0;