威尔逊-theta法MATLAB程序及算例说明.docx

单自由度：function u,u1,u2=wiltheta(m,k,c,d,e,f,p,theta,t)n=length(p);u=zeros(n,1);u1=zeros(n,1);u2=zeros(n,1);uu=zeros(n,1);pp=zeros(n,1);u(1)=d;u1(1)=e;u2(1)=f;l=6/(theta*t)2;q=3/(theta*t);r=6/(theta*t);s=theta*t/2;kk=k+l*m+q*c;for i=1:n-1 pp(i)=p(i)+theta*(p(i+1)-p(i)+(l*u(i)+r*u1(i)+2*u2(i)*m+(q*u(i)+2*u1(i)+s*u2(i)*c; uu(i)=pp(i)/kk; u2(i+1)=u2(i)+l/theta*(uu(i)-u(i)-l*t*u1(i)-3/theta*u2(i); u1(i+1)=u1(i)+t/2*(u2(i+1)+u2(i); u(i+1)=u(i)+t*u1(i)+t2/6*(u2(i+1)+2*u2(i+1);endend实例：，初始时刻位移、速度、加速度均为零。以下为0.02-1时间段时程曲线。其中p=-1*0.000508665 0.000906218 0.00134557 0.00169929 0.00212334 0.00248216 0.00289602 0.00326809 0.00366972 0.0040581 0.00443833 0.00485321 0.00519776 0.0056473 0.00488073 -0.000121305 -0.00759837 -0.00272681 0.00588379 0.00300815 0.00207034 0.00169215 0.001187971 0.00298879 0.00417125 -0.000565749 -0.00590112 -0.00725178 0.00413965 0.00858614 -0.00128848 -0.00422426 0.00104893 0.00502956 0.00237105 -0.00238736 -0.00554638 -0.00149032 0.00601325 0.00408767 0.000469929 -0.00323547 -0.00552192 2.34455E-05 0.0111213 0.000141692 -0.0116616 -0.0109857 0.00715902 0.00928542；多自由度function u u1 u2=multi_wiltheta(M,K,C,d,e,f,P,theta,t)n=length(P(:,1);m=length(P(1,:);u=zeros(n,m);u1=zeros(n,m);u2=zeros(n,m);uu=zeros(n,m);PP=zeros(n,m);u(:,1)=d;u1(:,1)=e;u2(:,1)=f;l=6/(theta*t)2;q=3/(theta*t);r=6/(theta*t);s=theta*t/2;KK=K+l*M+q*C;for i=1:m-1 pp(:,i)=P(:,i)+theta*(P(:,i+1)-P(:,i)+M*(l*u(:,i)+r*u1(:,i)+2*u2(:,i)+C*(q*u(:,i)+2*u1(:,i)+s*u2(:,i); uu(:,i)=KKpp(:,i); u2(:,i+1)=u2(:,i)+l/theta*(uu(:,i)-u(:,i)-l*t*u1(:,i)-3/theta*u2(:,i); u1(:,i+1)=u1(:,i)+t/2*(u2(:,i+1)+u2(:,i); u(:,i+1)=u(:,i