地震工程学反应谱和地震时程波的相互转化matlab编程Word文档格式.docx

1、 %绝对加速度% *A,B矩阵*Damp=; %阻尼比TA=:6; %TA=: %结构周期 Dt=; %地震记录的步长 %记录计算得到的反应，MaxD为某阻尼时最大相对位移，MaxV为某阻尼最大相对速度，MaxA某阻尼时最大绝对加速度，用于画图MaxD=zeros（3,length（TA）;MaxV=zeros（3,length（TA）;MaxA=zeros（3,length（TA）;t=1;for T=:6NatualFrequency=2*pi/T ; %结构自振频率DampFrequency=NatualFrequency*sqrt（1-Damp*Damp）; %计算公式化简e_t=ex

2、p（-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-Damp*Damp）;A（2,2）=e_t*（-s*Damp/sqrt（1-Damp*Damp）+c）;d_f=（2*Damp2-1）/（NatualFrequency2*Dt）;d_3t=Damp/（Natual

3、Frequency3*Dt）;B=zeros（2,2）;B（1,1）=e_t*（d_f+Damp/NatualFrequency）*s/DampFrequency+（2*d_3t+1/NatualFrequency2）*c）-2*d_3t;B（1,2）=-e_t*（d_f*s/DampFrequency+2*d_3t*c）-1/NatualFrequency2+2*d_3t;B（2,1）=e_t*（d_f+Damp/NatualFrequency）*（c-Damp/sqrt（1-Damp2）*s）-（2*d_3t+1/NatualFrequency2）*（DampFrequency*s+Dam

4、p*NatualFrequency*c）+1/（NatualFrequency2*Dt）;B（2,2）=e_t*（1/（NatualFrequency2*Dt）*c+s*Damp/（NatualFrequency*DampFrequency*Dt）-1/（NatualFrequency2*Dt）;for i=1:（N-1） %根据地震记录,计算不同的反应Displace（i+1）=A（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）*Displ

5、ace（i）+A（2,2）*Velocity（i）+B（2,1）*Accelerate（i）+B（2,2）*Accelerate（i+1）;AbsAcce（i+1）=-2*Damp*NatualFrequency*Velocity（i+1）-NatualFrequency2*Displace（i+1）;endMaxD（1,t）=max（abs（Displace）;MaxV（1,t）=max（abs（Velocity）;if T=MaxA（1,t）=max（abs（Accelerate）;elseMaxA（1,t）=max（abs（AbsAcce）;%初始化各储存向量，避免下次不同周期计算时引用

6、到前一个周期的结果t=t+1;End% *PLOT*close all figure %绘制地震记录图plot（time（:）,Accelerate（:） title（PEER STRONG MOTION DATABASE RECORD）xlabel（time（s）ylabel（acceleration（g） gridfigure %绘制位移反应谱plot（TA,MaxD（1,:）,TA,MaxD（2,:-r,TA,MaxD（3,:kDisplacementTn（s）Displacement（m）legend（=Gridfigure %绘制速度反应谱plot（TA,MaxV（1,:,TA,Ma

7、xV（2,:,TA,MaxV（3,:Velocityvelocity（m/s）figure %绘制绝对加速度反应谱plot（TA,MaxA（1,:,TA,MaxA（2,:,TA,MaxA（3,:Absolute Accelerationabsolute acceleration（m/s2）3 运行的结果得到的反应谱 图2 位移反应谱图3 速度反应谱图4 加速度反应谱一、 反应谱生成地震波1 所取的反应谱为上海市设计反应谱图5 上海市设计反应谱2 反应谱取值程序为：%规范反应谱取值程序 参照01年抗震规范function rs_z=r_s_1（pl,zn,ld,cd,fz） %pl 圆频率,zn

8、阻尼比,ld烈度,cd场地类型,场地分组fz%烈度选择if ld=6 arfmax=;if ld=7if ld=8if ld=9%场地类别，设计地震分组选择if cd=1 if fz=1 Tg=; end if fz=2 if fz=3if cd=2if cd=3if cd=4%ceita=zn; %阻尼比lmt1=+/8;if lmt1 lmt1=0;lmt2=1+/+*ceita）;if lmt2 lmt2=;sjzs=+/+5*ceita）;%分段位置 T1 T2 T3T1=;T2=Tg;T3=5*Tg;T_jg=2*pi./pl;% 第一段 0T1if T_jg=T1 arf_jg=*

9、arfmax+（lmt2*arfmax）/*T_jg;% 第二段 T1T2if T1T_jg&T_jg=T2 arf_jg=lmt2*arfmax;% 第三段 T2T3if T2=T3 arf_jg=（Tg/T_jg）sjzs）*lmt2*arfmax;% 第四段 T3if T3= arf_jg=（lmt2*sjzs-lmt1*（T_jg-5*Tg）*arfmax;% 第五段 if T_jg arf_jg=（lmt2*sjzs-lmt1*Tg）*arfmax;%反应谱值 拟加速度值rs_z=arf_jg*;3 生成人造地震波主程序：%主程序%确定需要控制的反应谱Sa（T）（T=T1,.,TM）

10、的坐标点数M， 反应谱控制容差rcTyz=:,:;rc=;nTyz=length（Tyz）;ceita=;%阻尼比：nTyz Syz（i）=r_s_1（2*pi/Tyz（i）,ceita,8,2,1）; %8度，2类场地，第1地震分组% 变换的频率差：2*pi*（可以保证长周期项5s附近有5项三角级数）；%频率变化范围 N1=30, 30*2*pi ;N2=3000, 5000*2*pi plc=2*pi*;pl=30*2*pi:*2*pi:10000*2*pi;npl=length（pl）;P=; %保证率%人造地震动持续时间40s，时间间隔：Td=40;dt=;t=0:40;nt=leng

11、th（t）;% 衰减包络函数t1=8; %上升段t2=8+24; %平稳段; 下降段则为40328sc=; %衰减段参数nt if t（i）t1 & t（i）rc %循环体函数 blxs=Syz./rspa; for xsxs=1: if 2*pi/pl（xsxs）=Tyz（sxsx） & （2*pi/pl（xsxs）Tyz（nTyz） blxs1（xsxs）=blxs（nTyz）; Aw=Aw.*blxs1; % 合成地震动 at=zeros（nt,1）; atj=zeros（nt,1）; % 计算反应谱验证是否满足rc在5%的要求 % response spectra of callida

12、r % parameter g=; jsnum=jsnum+1 max（rcrsp（:）% 最终的反应谱 与规范谱% Tjs=: % nTjs=length（Tjs）; m=1; x0=0; v0=0; ww=2*pi./Tyz; % load ag=at; % solution for y=1: rspa（i）=maxas（i）/g; rspa_S（i）=r_s_1（2*pi/Tyz（i）,ceita,8,2,1）/g;subplot（2,1,1）;plot（t,at）;subplot（2,1,2）;plot（Tyz,rspa）;hold on;plot（Tyz,rspa_S）;4 生成的人造地震波如图所示。图6 人造地震波和初始反应谱