地震工程学反应谱和地震时程波的相互转化matlab编程Word文档格式.docx
《地震工程学反应谱和地震时程波的相互转化matlab编程Word文档格式.docx》由会员分享，可在线阅读，更多相关《地震工程学反应谱和地震时程波的相互转化matlab编程Word文档格式.docx（17页珍藏版）》请在冰豆网上搜索。
%N读入的记录的量
time=0:
0.005:
（N-1）*0.005;
%单位s
%初始化各储存向量
Displace=zeros（1,N）;
%相对位移
Velocity=zeros（1,N）;
%相对速度
AbsAcce=zeros（1,N）;
%绝对加速度
%***********A,B矩阵***********
Damp=0.02;
%阻尼比0.02
TA=0.0:
0.05:
6;
%TA=0.000001:
0.02:
%结构周期
Dt=0.005;
%地震记录的步长
%记录计算得到的反应，MaxD为某阻尼时最大相对位移，MaxV为某阻尼最大相对速度，MaxA某阻尼时最大绝对加速度，用于画图
MaxD=zeros（3,length（TA））;
MaxV=zeros（3,length（TA））;
MaxA=zeros（3,length（TA））;
t=1;
forT=0.0:
6
NatualFrequency=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-Damp*Damp）;
A（2,2）=e_t*（-s*Damp/sqrt（1-Damp*Damp）+c）;
d_f=（2*Damp^2-1）/（NatualFrequency^2*Dt）;
d_3t=Damp/（NatualFrequency^3*Dt）;
B=zeros（2,2）;
B（1,1）=e_t*（（d_f+Damp/NatualFrequency）*s/DampFrequency+（2*d_3t+1/NatualFrequency^2）*c）-2*d_3t;
B（1,2）=-e_t*（d_f*s/DampFrequency+2*d_3t*c）-1/NatualFrequency^2+2*d_3t;
B（2,1）=e_t*（（d_f+Damp/NatualFrequency）*（c-Damp/sqrt（1-Damp^2）*s）-（2*d_3t+1/NatualFrequency^2）*（DampFrequency*s+Damp*NatualFrequency*c））+1/（NatualFrequency^2*Dt）;
B（2,2）=e_t*（1/（NatualFrequency^2*Dt）*c+s*Damp/（NatualFrequency*DampFrequency*Dt））-1/（NatualFrequency^2*Dt）;
fori=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）*Displace（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）-NatualFrequency^2*Displace（i+1）;
end
MaxD（1,t）=max（abs（Displace））;
MaxV（1,t）=max（abs（Velocity））;
ifT==0.0
MaxA（1,t）=max（abs（Accelerate））;
else
MaxA（1,t）=max（abs（AbsAcce））;
%初始化各储存向量，避免下次不同周期计算时引用到前一个周期的结果
t=t+1;
End
%***********PLOT***********
closeall
figure%绘制地震记录图
plot（time（:
）,Accelerate（:
））
title（'
PEERSTRONGMOTIONDATABASERECORD'
）
xlabel（'
time（s）'
ylabel（'
acceleration（g）'
）
grid
figure%绘制位移反应谱
plot（TA,MaxD（1,:
）,'
-.b'
TA,MaxD（2,:
-r'
TA,MaxD（3,:
:
k'
Displacement'
Tn（s）'
Displacement（m）'
legend（'
ζ=0.02'
Grid
figure%绘制速度反应谱
plot（TA,MaxV（1,:
TA,MaxV（2,:
TA,MaxV（3,:
Velocity'
velocity（m/s）'
figure%绘制绝对加速度反应谱
plot（TA,MaxA（1,:
TA,MaxA（2,:
TA,MaxA（3,:
AbsoluteAcceleration'
absoluteacceleration（m/s^2）'
3运行的结果得到的反应谱
图2位移反应谱
图3速度反应谱
图4加速度反应谱
一、反应谱生成地震波
1所取的反应谱为上海市设计反应谱
图5上海市设计反应谱
2反应谱取值程序为：
%%规范反应谱取值程序参照01年抗震规范
functionrs_z=r_s_1（pl,zn,ld,cd,fz）%%%pl圆频率,zn阻尼比,ld烈度,cd场地类型,场地分组fz
%%%%烈度选择
ifld==6
arfmax=0.11;
ifld==7
arfmax=0.23;
ifld==8
arfmax=0.45;
ifld==9
arfmax=0.90;
%%%%场地类别，设计地震分组选择
ifcd==1
iffz==1
Tg=0.25;
end
iffz==2
Tg=0.30;
iffz==3
Tg=0.35;
ifcd==2
Tg=0.40;
Tg=0.45;
ifcd==3
Tg=0.55;
Tg=0.65;
ifcd==4
Tg=0.75;
Tg=0.90;
%%%%%%%%%
ceita=zn;
%%%%%阻尼比
lmt1=0.02+（0.05-ceita）/8;
iflmt1<
lmt1=0;
lmt2=1+（0.05-ceita）/（0.06+1.7*ceita）;
iflmt2<
0.55
lmt2=0.55;
sjzs=0.9+（0.05-ceita）/（0.5+5*ceita）;
%%%%%分段位置T1T2T3
T1=0.1;
T2=Tg;
T3=5*Tg;
T_jg=2*pi./pl;
%%%%第一段0～T1
ifT_jg<
=T1
arf_jg=0.45*arfmax+（lmt2*arfmax-0.45*arfmax）/0.1*T_jg;
%%%%第二段T1～T2
ifT1<
T_jg&
T_jg<
=T2
arf_jg=lmt2*arfmax;
%%%%第三段T2~T3
ifT2<
=T3
arf_jg=（（Tg/T_jg）^sjzs）*lmt2*arfmax;
%%%%第四段T3～6.0
ifT3<
=6.0
arf_jg=（lmt2*0.2^sjzs-lmt1*（T_jg-5*Tg））*arfmax;
%%%%第五段6.0～
if6.0<
T_jg
arf_jg=（lmt2*0.2^sjzs-lmt1*（6.0-5*Tg））*arfmax;
%%%%%%反应谱值拟加速度值
rs_z=arf_jg*9.8;
3生成人造地震波主程序：
%%%主程序%%%%
%%%%确定需要控制的反应谱Sa（T）（T=T1,...,TM）的坐标点数M，反应谱控制容差rc
Tyz=[0.04:
0.016:
0.1,0.15:
3.0,3.2:
5.0];
rc=0.06;
nTyz=length（Tyz）;
ceita=0.035;
%%%阻尼比：
0.035
nTyz
Syz（i）=r_s_1（2*pi/Tyz（i）,ceita,8,2,1）;
%%%%8度，2类场地，第1地震分组
%%%%%%变换的频率差：
2*pi*0.005（可以保证长周期项5s附近有5项三角级数）；
%%%%频率变化范围N1=30,30*0.005*2*pi;
N2=3000,5000*0.005*2*pi
plc=2*pi*0.005;
pl=30*0.005*2*pi:
0.005*2*pi:
10000*0.005*2*pi;
npl=length（pl）;
P=0.9;
%%%保证率
%%%%%%人造地震动持续时间40s，时间间隔：
0.02s
Td=40;
dt=0.02;
t=0:
40;
nt=length
升级会员 