# suanli1.zip

• 谁谁谁0311
了解作者
• matlab
开发工具
• 1KB
文件大小
• zip
文件格式
• 0
收藏次数
• 10 积分
下载积分
• 17
下载次数
• 2017-12-05 14:15
上传日期

suanli1.zip
• 新建文本文档.txt
2.5KB

fid=fopen('Ei-centro.txt','r');%读取地震波 dzhbo=fscanf(fid,'%f',1501); xg=dzhbo/max(abs(dzhbo))*2; dt=0.02; t=(0:length(dzhbo)-1)*dt; syms w; K=1.0e8*[4 -2 0;-2 4 -2;0 -2 2];%刚度矩阵 M=1.0e5*[4 0 0;0 4 0;0 0 4];%质量矩阵 a=K-w^2*M;%频率方程 b=solve(det(a));%求解频率方程行列式等于0 z=real(double(b)); w1=z(3);w2=z(5);zeta1=0.05;zeta2=0.05;%解得频率w1，w2 arphc=(2*w1*w2*(zeta1*w2-zeta2*w1))/(w2^2-w1^2); betac=(2*(zeta2*w2-zeta1*w1))/(w2^2-w1^2); C=arphc*M+betac*K;%阻尼矩阵 Bs=[1 -1 0;0 1 -1;0 0 1];%地震作用位置矩阵 Ds=-M*ones(3,1);%环境位置干扰矩阵 A=[zeros(3) eye(3);inv(-M)*K inv(-M)*C]; B=[zeros(3);inv(M)*Bs]; D=[zeros(3,1);inv(M)*Ds]; I=eye(3); arph=100; beta=[0.01 0.1 0.5 1 3 5 8 12 15 20 40 60]*1.0e-6; for i=1:12; betai=beta(i); Q=arph*[K zeros(3);zeros(3) M]; R=betai*I; G=lqr(A,B,Q,R);%连续状态方程控制力状态反馈增益矩阵 [y,V]=lsim(A,D,eye(6),zeros(6,1),xg,t);%无控时输入响应 [Y,Z]=lsim((A-B*G),D,eye(6),zeros(6,1),xg,t);%有控时输入响应 V1=A*V'+D*xg';%无控结构状态方程 Z1=(A-B*G)*Z'+D*xg';%有空结构状态方程 U=-G*Z'/1000;%最优控制力 f01=max(abs(V(:,1))); %各层无控位移：m max(abs(V(:,2))); max(abs(V(:,3))); f02=max(abs(V(:,2)))-max(abs(V(:,1)));%第二层无控层间位移：m f03=max(abs(V(:,3)))-max(abs(V(:,2)));%第三层无控层间位移：m max(abs(V(:,4)));%各层无控速度：m max(abs(V(:,5))); max(abs(V(:,6))); f11=max(abs(Z(:,1))); %各层有控位移：m max(abs(Z(:,2))); max(abs(Z(:,3))); f12=max(abs(Z(:,2)))-max(abs(Z(:,1)));%第二层有控层间位移：m f13=max(abs(Z(:,3)))-max(abs(Z(:,2)));%第三层有控层间位移：m f=[f01 f02 f03 f11 f12 f13]; max(abs(Z(:,4)));%各层有控速度：m max(abs(Z(:,5))); max(abs(Z(:,6))); max(abs(V1(1,:))); %各层无控速度：m max(abs(V1(2,:))); max(abs(V1(3,:))); s01=max(abs(V1(4,:)));%各层无控加速度：m s02=max(abs(V1(5,:))); s03=max(abs(V1(6,:))); max(abs(Z1(1,:))); %各层有控速度：m max(abs(Z1(2,:))); max(abs(Z1(3,:))); s11=max(abs(Z1(4,:)));%各层有控加速度：m s12=max(abs(Z1(5,:))); s13=max(abs(Z1(6,:))); s=[s01 s02 s03 s11 s12 s13]; n(i)=max(abs(U(1,:))); %各层最大控制力：KN max(abs(U(2,:))); max(abs(U(3,:))); a1=figure; plot(t,V(:,1),'--k',t,Z(:,1),'-r'); xlabel('时间/s'); ylabel('第一层位移/m'); a2=figure; plot(t,U(1,:)); xlabel('时间/s'); ylabel('第一层控制力/m'); weiyi(i)=f11; end c=figure; plot(beta,weiyi*100,'-s',beta,n/500); xlabel('β'); ylabel('第一层位移/cm和第一层控制力/500KN'); figure(30) hold on o1=0:0.02:30; Y1=o1; o2=1:1:1501; y1=xg(o2,1); plot(Y1,y1,'k'); grid on xlabel('时间') ylabel('地震加速度') title('地震动')

相关推荐