本帖最后由 ME! 于 2014-3-13 16:50 编辑
你是用时间积分法做的吧,它是无条件稳定的,表示在0.55的时候趋于稳定值了;轴心轨迹也是对的
我想请教下楼主轴心轨迹的图是怎么画出来的,我的画出来就是一个圆
dt=pi/100; %步长
t0=0; %初始时刻
t=t0:dt:100*pi; %时间向量
zeta=ones(Sys_dof,1)*0.1; % 主振型阻尼比
C=diag((2*wn).*zeta); %对角化阻尼矩阵
alpha=0.3; beta=0.6; % 稳定条件
ff2=zeros(Sys_dof,1);
P2=[10,1,2]; %加载力
ff2(No_dof*(P2(2)-1)+P2(3))=P2(1); %加载力在某一节点上
u2=cos(t); %u1=sin(t);
ff2=ff2*u2;
dsp=zeros(length(K),length(t)); % 位移
vel=zeros(length(K),length(t)); % 速度
acc=zeros(length(K),length(t)); % 加速度
acc(:,1)=inv(M)*(ff2(:,1)-K*dsp(:,1)-C*vel(:,1)); % 计算初始加速度 (t=0)
ekk=K+M/(alpha*dt^2)+C*beta/(alpha*dt); % 计算有效刚度矩阵
for it=1:length(t)-1 % 时间步循环
cfm=dsp(:,it)/(alpha*dt^2)+vel(:,it)/(alpha*dt)+acc(:,it)*(0.5/alpha-1);
cfc=dsp(:,it)*beta/(alpha*dt)+vel(:,it)*(beta/alpha-1)...
+acc(:,it)*(0.5*beta/alpha-1)*dt;
efd=ff2(:,it)+M*cfm+C*cfc; % 计算有效力矢量
dsp(:,it+1)=inv(ekk)*efd; % t+dt时刻的位移
acc(:,it+1)=(dsp(:,it+1)-dsp(:,it))/(alpha*dt^2)-vel(:,it)/(alpha*dt)...
-acc(:,it)*(0.5/alpha-1); % t+dt时刻的加速度
vel(:,it+1)=vel(:,it)+acc(:,it)*(1-beta)*dt+acc(:,it+1)*beta*dt; % t+dt时刻的速度
end
加载部分:u2=sin(t);u2=cos(t);提取位移画图,不知道是不是这样画,请指教下,万分感谢
|