第 2 章function VTB2(m,c,k,x0,v0,tf,w,f0)%单自由度系统得谐迫振动clcwn=sqrt(k/m);z=c/2/m/wn;lan=w/wnwd=wn*sqrt(1-z^2);A=sqrt(((v0+z*wn*x0)^2+(x0*wd)^2)/wd^2);t=0:tf/1000:tf;phi1=atan2(x0*wd,v0+z*wn*x0);phi2=atan2(2*z*lan,1-lan^2);B=wn^2*f0/k/sqrt((wn^2-w^2)^2+(2*z*wn*w)^2);x=A*exp(-z*wn*t)、*sin(sqrt(1-z^2)*wn*t+phi1)+B*sin(w*t-phi2);plot(t,x),gridxlabel('时间(s)')ylabel('位移')title('位移与时间得关系')function VTB1(m,c,k,x0,v0,tf)%VTB1 用来计算单自由度有阻尼自由振动系统得响应%VTB1 绘出单自由度有阻尼自由振动系统得响应图%m 为质量;c 为阻尼;k 为刚度;x0 为初始位移;v0 为初始速度;tf 为仿真时间%VTB1(zeta,w,x0,v0,tf)绘出单自由度有阻尼自由振动系统得响应图%程序中 z 为阻尼系数 ξ;wn 为固有频率 ωn;A 为振动幅度;phi 为初相位 θclcwn=sqrt(k/m);z=c/2/m/wn;wd=wn*sqrt(1-z^2);fprintf('固有频率为%、3g、rad/s、\n',wd);fprintf('阻尼系数为%、3g、\n',z);fprintf('有阻尼得固有频率为%、3g、rad/s、\n',wd);t=0:tf/1000:tf;if z<1 A=sqrt(((v0+z*wn*x0)^2+(x0*wd)^2)/wd^2); phi=atan2(x0*wd,v0+z*wn*x0) x=A*exp(-z*wn*t)、*sin(wd*t+phi); fprintf('A=%、3g\n',A); elseif z==1 a1=x0; a2=v0+wn*x0; fprintf('a1=%、3g\n',a1); fprintf('a2=%、3g\n',a2); x=(a1+a2*t)、*exp(-wn*t);else a1=(-v0+(-z+sqrt(z^2-1))*wn*x0)/2/wn/sqrt(z^2-1); a2=(v0+(z+sqrt(z^2-1))*wn*x0)/2/wn/sqrt(z^2-1); fprintf('a1=%、3g\n',a1); fprintf('a2=%、3g\n',a2); x=exp(-z*wn*t)、*(a1*exp(-wn*sqrt(z^2-1)*t)+a2*exp(wn*sqrt(z^2-1)*t));endplot(t,x),gridxlabel('时间(s)')ylabel('位移')title('位移与时间得关系')function jzdd%矩阵迭代法求系统得三阶固有频率与主阵型clear allclose allfid1=fopen('A-1','wt'); %建立主振型文件fid2=fopen('B-1','wt'); %建立固有频率文件%输入质量矩阵M(1,1)=2;M(2,2)=1、5;M(3,3)=1;%输入刚度矩阵K(1,1)=5;K(1,2)=-2;K(2,1)=-2;K(2,2)=3;K(2,3)=-1;K(3,2)=-1;K(3,3)=1%计算特征值与特征向量D=inv(K)*M; %原始动力矩阵A=ones(3,1); %初始振型for i=1:3 %计算三阶固有频率与主振型 pp0=0; i B=D*A; pp=1、0/B(3); %B(3)为 B 中得最后一个元素 A=B/B(3); while abs((pp-pp0)/pp)>1、e-6 pp0=pp; B=D*A; pp=...
