matlab求解模态,[转载]MATLAB与模态叠加法求解高等动力学问题，两者误差比较.pdfVIP

matlab求解模态,[转载]MATLAB与模态叠加法求解⾼等动⼒学 问题，两者误差⽐较 图⽰表⽰四层楼的抗剪模型，其剪切刚度系数及楼板质量均表⽰在图中，在顶层受⼀⽔平⼒的简谐激振⼒pcos(Ωt)，仅考虑其稳态响应。 求顶层的响应函数。 %%%%%%%%%% 运动⽅程求解 %%%%%% 运动⽅程求解 %%%% T=input(Time is ); tb=input(Buchang is ); %中⼼差分法 t=0:tb:T; M=[1 0 0 0; 0 2 0 0; 0 0 2 0; 0 0 0 3]; K=[800 -800 0 0; -800 2400 -1600 0; 0 -1600 4000 -2400; 0 0 -2400 5600]; a0=(tb^2)1; a1=(2*tb)1; a2=2*a0; u0=[0.00 0.00 0.00 0.00];v0=[0.00 0.00 0.00 0.00];ac0=[100.00 0.00 0.00 0.00]; u00=u0-tb*v0+a2ac0; %(-tb时刻的位移) DM=a0*M; %等效质量阵 n=fix(tbT); u=zeros(4,n); %定义４x N的位移矩阵 u(:,1)=u00;u(:,2)=u0; t=zeros(1,n); %定义１x N的时间矩阵 t(1,1)=tb;t(1,2)=2*tb; for i=3:n ti=i*tb; t(1,i)=ti; Qt=[100.00*cos(0.5*13.294*(i-1)*tb) 0.00 0.00 0.00]; %载荷 DQt=Qt-(K-a2*M)*u(:,i-1)-a0*M*u(:,i-2); %等效载荷 u(:,i)=DMDQt; dis=u(1,i); %顶层楼的位移 end %模态叠加法 w0=[13.294 29.660 41.079 55.882]; %固有频率 w1=w0.^2; Q=[1.00000 1.00000 -0.90145 0.15436; 0.77910 -0.09963 1.00000 -0.44817; 0.49655 -0.53989 -0.15859 1.00000; 0.23506 -0.43761 -0.70797 -0.63688]; %模态矩阵 M1=Q(:,1)*M*Q(:,1);M2=Q(:,2)*M*Q(:,2);M3=Q(:,3)*M*Q(:,3);M4=Q(:,4)*M*Q(:,4); K1=w1(1)*M1;K2=w1(2)*M2;K3=w1(3)*M3;K4=w1(4)*M4; %模态刚度 F=[100.00 0.00 0.00 0.00]*cos(0.5*13.294*t); F1=Q(:,1)*F;F2=Q(:,2)*F;F3=Q(:,3)*F;F4=Q(:,4)*F; %模态⼒ p=0.5*w0(1); %激振⼒频率 w3=p^2; u1=1.000*F1/(K1*(1-w3/w1(1)))+1.000*F2/( K2*(1-w3/w1(2)))+ (-0.90145)*F3/(K3*(1-w3/w1(3)))+0.15436*F4/(K4*(1-w3/w1(4))); %--顶层楼稳的态响应幅值 s=u1; %计算均⽅差 l1=s; l2=dis; y0=zeros(1,n); y0=l1-l2; wc=nnorm(y0) plot(t,l1,t(1,:),u(1,:),k--); title(模态叠加法数值和中⼼差分法数值); xlabel(时间); ylabel(位移); legend((- 模态叠加法数值),(-- 中⼼差分法数值)); grid on; 分别取T=10，100，1000 Time is 10 Buchang is 0.01 wc = 0.0156 Time is 100 Buchang is 0.01 wc = 0.0030 Time is 1000 Buchang is 0.01 wc = 0.0011