《最新地震工程作业哈工大MATLAB程序资料.doc》由会员分享，可在线阅读，更多相关《最新地震工程作业哈工大MATLAB程序资料.doc（14页珍藏版）》请在三一文库上搜索。
1、精品文档10 203040时间（s）506070戶卩度速加-0.250地震工程大作业哈工大李金平弹性反应谱原创性声明，主程序由本人独立编写完成，支持各种检验。所选地震动 RSN2345CHICHI.AT2 peer 索取号：23450.250.20.150.12)0.050-0.05-0.1-0.15-0.2精品文档位移反应谱mEL移位局部放大9.559.5)mm9.45移9.4位9.359.39.79.659.69.259.211.121.021.041.061.081.1周期T (s)求比值1.0411LL1Seismosignal 位移 /MATLAB 位移比值1.031.021.010
2、.990.980.970.960.950.94 匚0910周期T (s)绝对加速度反应谱周期T (s)局部放大-0.25-0.2-0.15-0.1-0.0500.050.10.150.20.25周期T (s)精品文档精品文档求比值进行比较1.041.031INI,LIL1Seismosignal力廿速度/MATLAB加速度比值1.021.0110.990.980.970.960.950.940 15678910周期T (s)速度反应谱局部放大观察差异卩度速X 102015105000.020.040.060.080.10.12周期T (s)求比值进行比较。1.3Seismosignal 速度
3、/MATLAB 速度 比值1.21.11值比0.90.80.7efr亍rr012345678910周期T (s)结论：对比：拟合效果非常好，短周期 0s<T<0.2s 误差较大，其中相对速度谱的误差较大，接近 30%，加速度和位移误差都在 3%以内，周期 0.2s<T<10s 各反应谱的误差都保持在 1%以内， 精度相同，猜想，两种计算程序应该用的同一种计算方法( newmark 法)。5%)不同周期的单自理解： 本次计算的反应谱是在给定的地震动下，具有相同阻尼比( 由度结构的线弹性反应幅值，得到的速度、加速度、移幅值随周期变化的三条曲线。 a1(i)=max(abs(
4、a+xg');v1(i)=max(abs(v);x1(i)=max(abs(x);用法： 今后给定一个结构，我们可以计算其周期T1,然后在相应的反应谱图表里面找到对应的谱值如 Sa、Sv Sd，即为我们所要求得的该结构在指定地震动作用下最大动力反 应，将动力分析简化成为静力计算，简单方面。 (弹性范围内)主程序 DZZY.mclear clcM=1;%fid=fope n( 'RSN2345CHICHI.txt');% 读取地震动加速度记录 xg=9.8*fscanf(fid,'%f');fclose(fid); n=length(xg);精品文档F=-
5、M*xg'%生成地震力for i=1:1000tn(i)=0.01*i;K=2*pi*2*pi/tn(i)/tn(i);lamda=2*pi/tn(i);dt=0.005;t=(0:0.005:(n-1)*0.005)'x0=0;v0=0;a0=0;ksi=0.01*5;DAMPER=2*ksi*lamda*M; x,v,a=newmarkb(M,K,DAMPER,1,F,x0,v0,a0,dt,n); a1(i)=max(abs(a+xg');v1(i)=max(abs(v);x1(i)=max(abs(x);endtest=importdata('a2.tx
6、t');%读取 Seismosignal 绝对加速度a2=test(:,2)*9.8;testv=importdata('v2.txt');%读取 Seismosignal 速度v2=testv(:,2);testd=importdata('x2.txt');%读取 Seismosignal 相对位移x2=testd(:,2);tn=0.01:0.01:10;figure (1)plot(tn,x1*1000,'r-')xlabel('周期 T (s)')ylabel(' 位移 (mm)')hold onp
7、lot(tn,x2*1000,'linewidth',2) legend('matlab 相对位移 ','Seismosignal 相对位移 ')%figure (2)plot(tn,a1,'r-','linewidth',1)xlabel('周期 T (s)')ylabel('Sa (m/sA2)')hold on plot(tn,a2,'k','linewidth',2) legend('matlab-Sa','Seismo
8、signal-Sa') figure (3)plot(tn,v1)xlabel('周期 T (s)')ylabel(' 速度 (m/s)')hold onplot(tn,v2,'linewidth',2)legend('matlab-Sv','Seismosignal-Sv')hold off figure (4) hold on plot(t,xg)legend('地震动加速度') xlabel(' 时间 (s)') ylabel('加速度(m/sA2)')
9、for j=1:1000%求比值ca(j)=a2(j)/a1(j);cv(j)=v2(j)/v1(j);cx(j)=x2(j)/x1(j);endfigure (6) plot(tn,ca)legend('Seismosignal 加速度 /MATLAB 加速度比值 ') xlabel('周期 T (s)')ylabel(' 比值 ')figure (7) plot(tn,cv) legend('Seismosignal 速度 /MATLAB 速度 比值 ') xlabel('周期 T (s)')ylabel(&#
10、39; 比值 ')figure (8) plot(tn,cx) legend('Seismosignal 位移 /MATLAB 位移 比值 ') xlabel('周期 T (s)')ylabel(' 比值 ')子程序 newmarkb.mfunction x,v,a=newmarkb(M,K,C,N,P,x0,v0,a0,dt,RecordLength)% newmark-beta method%obtain the response of the dynamic system%x,v,a=newmarkb(M,K,C,N,P,x0,v0
11、,a0,dt,RecordLength)%M - mass matrix%K - stiffness matrix%C - damping matrix%N - DOF%P - loads%x0 - initial displacement%v0 - initial velocity%a0 - initial acceleration%dt - interval% RecordLength - number of sampling pointsx=zeros(N,RecordLength);v=zeros(N,RecordLength);a=zeros(N,RecordLength);x(:,
12、1)=x0;v(:,1)=v0;a(:,1)=a0;deta=0.50;alpha=0.25;a0=1/alpha/dtA2;a1=deta/alpha/dt;a2=1/alpha/dt;a3=1/2/alpha-1;a4=deta/alpha-1;a5=dt*(deta/alpha-2)/2;a6=dt*(1-deta);a7=deta*dt;K_=K+a0*M+a1*C;iK=inv(K_);for i=1:RecordLength-1P_(:,i+1)=P(:,i+1)+M*(a0*x(:,i)+a2*v(:,i)+a3*a(:,i)+C*(a1*x(:,i)+a4* v(:,i)+a5*a(:,i);x(:,i+1)=iK*P_(:,i+1);a(:,i+1)=a0*(x(:,i+1)-x(:,i)-a2*v(:,i)-a3*a(:,i);v(:,i+1)=v(:,i)+a6*a(:,i)+a7*a(:,i+1);end精品文档