《地震工程学反应谱和地震时程波的相互转化matlab编程》由会员分享，可在线阅读，更多相关《地震工程学反应谱和地震时程波的相互转化matlab编程（15页珍藏版）》请在金锄头文库上搜索。
1、TONGJI UNIVERSITY地震工程学作业名:课程名称：地震工程学史先飞号:1232627一、地震波生成反应谱1所取的地震波为Elcentro地震波加速度曲线，如图1所示。400BD02001003-100-200-2D0壬口 I图1 Elcentro地震波加速度曲线2所调用的Matlab程序为：% * 读入地震记录 *ElCentro;Accelerate= ElCentro(:,1)*9.8067;%单位统一为 m 和 sN=length(Accelerate);%N 读入的记录的量time=0:0.005:(N-1)*0.005; %单位 s%初始化各储存向量Displace=ze
2、ros(1,N); %相对位移Velocity=zeros(1,N); %相对速度AbsAcce=zeros(1,N); %绝对加速度% *A,B 矩阵*Damp=0.02; %阻尼比 0.02TA=0.0:0.05:6; %TA=0.000001:0.02:6; %结构周期Dt=0.005; %地震记录的步长%记录计算得到的反应，MaxD为某阻尼时最大相对位移，MaxV为某阻尼最大相对速度，MaxA某阻尼时最大 绝对加速度，用于画图MaxD=zeros(3,length(TA);MaxV=zeros(3,length(TA);MaxA=zeros(3,length(TA);t=1;for T
3、=0.0:0.05:6NatualFrequency=2*pi/T ; %结构自振频率DampFrequency=NatualFrequency*sqrt(1-Damp*Damp); %计算公式化简e_t=exp(-Damp*NatualFrequency*Dt);s=sin(DampFrequency*Dt);c=cos(DampFrequency*Dt);A=zeros(2,2);A(1,1)=e_t*(s*Damp/sqrt(1-Damp*Damp)+c);A(1,2)=e_t*s/DampFrequency;A(2,1)=-NatualFrequency*e_t*s/sqrt(1-Da
4、mp*Damp);A(2,2)=e_t*(-s*Damp/sqrt(1-Damp*Damp)+c);d_f=(2*Damp2T)/(NatualFrequency2*Dt); d_3t=Damp/(NatualFrequency3*Dt);B=zeros(2,2);B(l,l)=e_t*(d_f+Damp/NatualFrequency)*s/DampFrequency+(2*d_3t+l/NatualFrequency2)*c)-2*d_3 t;B(l,2)=-e_t*(d_f*s/DampFrequency+2*d_3t*c)T/NatualFrequency2+2*d_3t;B(2,l)
5、=e_t*(d_f+Damp/NatualFrequency)*(c-Damp/sqrt(l-Damp2)*s)-(2*d_3t+l/NatualFrequency2 )*(DampFrequency*s+Damp*NatualFrequency*c)+l/(NatualFrequency2*Dt);B(2,2)=e_t*(l/(NatualFrequency2*Dt)*c+s*Damp/(NatualFrequency*DampFrequency*Dt)T/(NatualF requency2*Dt);for i=1:(N-1) %根据地震记录,计算不同的反应Displace(i+1)=A(
6、1,1)*Displace(i)+A(1,2)*Velocity(i)+B(1,1)*Accelerate(i)+B(1,2)*Accelerate(i +1);Velocity(i+1)=A(2,1)*Displace(i)+A(2,2)*Velocity(i)+B(2,1)*Accelerate(i)+B(2,2)*Accelerate(i +1);AbsAcce(i+l)=-2*Damp*NatualFrequency*Velocity(i+l)-NatualFrequency2*Displace(i+l);endMaxD(1,t)=max(abs(Displace);MaxV(1,t)
7、=max(abs(Velocity);if T=0.0MaxA(1,t)=max(abs(Accelerate);elseMaxA(1,t)=max(abs(AbsAcce);endDisplace=zeros(l,N);%初始化各储存向量，避免下次不同周期计算时引用到前一个周期的结果Velocity=zeros(1,N);AbsAcce=zeros(1,N);t=t+1;End% *PLOT* close allfigure %绘制地震记录图plot(time(:),Accelerate(:)title(PEER STRONG MOTION DATABASE RECORD) xlabel(t
8、ime(s)ylabel(acceleration(g)gridfigure %绘制位移反应谱plot(TA,MaxD(1,:),-.b,TA,MaxD(2,:),-r,TA,MaxD(3,:),:k)title(Displacement)xlabel(Tn(s) ylabel(Displacement(m)legend( Z =0.02)Gridfigure %绘制速度反应谱plot(TA,MaxV(1,:),-.b,TA,MaxV(2,:),-r,TA,MaxV(3,:),:k) title(Velocity)xlabel(Tn(s)ylabel(velocity(m/s)legend(Z
9、=0.02)Gridfigure %绘制绝对加速度反应谱plot(TA,MaxA(1,:),-.b,TA,MaxA(2,:),-r,TA,MaxA(3,:),:k) title(Absolute Acceleration)xlabel(Tn(s)ylabel(absolute acceleration(m/s 2) legend(Z=0.02)Grid3 运行的结果得到的反应谱| c=ddq|；、r、i :JlrL-、.7.;i,|: jL/I3-*J-/LUiqslKimMii033O.IZ0 01&010M5Tn冊图2位移反应谱|&-+ 图3速度反应谱一、反应谱生成地震波1所取的反应谱为上
10、海市设计反应谱图5上海市设计反应谱2反应谱取值程序为：%规范反应谱取值程序参照01年抗震规范function rs_z=r_s_l(pl,zn,ld,cd,fz) %pl 圆频率,zn 阻尼比,ld 烈度,cd 场地类型，场地分组 fz %烈度选择if ld=6arfmax=0.11;endif ld=7arfmax=0.23;endif ld=8arfmax=0.45;endif ld=9arfmax=0.90;end%场地类别，设计地震分组选择if cd=1if fz=1Tg=0.25;endif fz=2Tg=0.30;endif fz=3Tg=0.35;endendif cd=2if
11、fz=1Tg=0.35;end if fz=2Tg=0.40;endif fz=3Tg=0.45;endendif cd=3if fz=1Tg=0.45;endif fz=2Tg=0.55;endif fz=3Tg=0.65;endendif cd=4if fz=1Tg=0.65;endif fz=2Tg=0.75;endif fz=3Tg=0.90;endend%ceita=zn; %阻尼比 lmt1=0.02+(0.05-ceita)/8;if lmt10lmt1=0;endlmt2=1+(0.05-ceita)/(0.06+1.7*ceita); if lmt20.55lmt2=0.55
12、;endsjzs=0.9+(0.05-ceita)/(0.5+5*ceita); %分段位置 T1 T2 T3T1=0.1;T2=Tg;T3=5*Tg;T_jg=2*pi./pl;%第一段0T1if T_jg=T1 arf_jg=0.45*arfmax+(lmt2*arfmax-0.45*arfmax)/0.1*T_jg;end% 第二段 T1T2if T1T_jg&T_jg=T2 arf_jg=lmt2*arfmax;end% 第三段 T2T3if T2T_jg&T_jg=T3 arf_jg=(Tg/T_jg)sjzs)*lmt2*arfmax;end% 第四段 T36.0if T3T_jg
13、&T_jg=6.0 arf_jg=(lmt2*0.2sjzs-lmtl*(T_jg-5*Tg)*arfmax;end% 第五段 6.0if 6.0T_jg arf_jg=(lmt2*0.2sjzs-lmtl*(6.0-5*Tg)*arfmax;end%反应谱值 拟加速度值rs_z=arf_jg*9.8;end3 生成人造地震波主程序：%主程序%确定需要控制的反应谱Sa (T)(T=T1,.,TM)的坐标点数M,反应谱控制容差rc Tyz=0.04:0.016:0.1,0.15:0.05:3.0,3.2:0.05:5.0;rc=0.06;nTyz=length(Tyz); ceita=0.035;%阻尼比：0.035for i=1:nTyzSyz(i)=r_s_1(2*pi/Tyz(i),ceita,8,2,1); %8度，2类场地，第1地震分组end%变换的频率差：2*pi*0.005(可以保证长周期项5s附近有5项三角级数)；%频率变化范围 N1=30, 30*0.005*2*pi ;N2=3000, 5000*0.005*2*piplc=2*pi*0.005;pl=30*0
