地震工程学作业
课程名称: 地震工程学______
指导老师:_______翟永梅_________
姓 名: 史先飞________
学 号: 1232627________
一、地震波生成反应谱
1 所取的地震波为Elcentro地震波加速度曲线,如图1所示。
图1 Elcentro地震波加速度曲线
2 所调用的Matlab程序为:
% ***********读入地震记录***********
ElCentro;
Accelerate= ElCentro(:,1)*9.8067;%单位统一为m和s
N=length(Accelerate);%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: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=0.0:0.05: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);
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)*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));
if T==0.0
MaxA(1,t)=max(abs(Accelerate));
else
MaxA(1,t)=max(abs(AbsAcce));
end
Displace=zeros(1,N);%初始化各储存向量,避免下次不同周期计算时引用到前一个周期的结果
Velocity=zeros(1,N);
AbsAcce=zeros(1,N);
t=t+1;
End
% ***********PLOT***********
close all
figure %绘制地震记录图
plot(time(:),Accelerate(:))
title('PEER STRONG MOTION DATABASE RECORD')
xlabel('time(s)')
ylabel('acceleration(g)')
grid
figure %绘制位移反应谱
plot(TA,MaxD(1,:),'-.b',TA,MaxD(2,:),'-r',TA,MaxD(3,:),':k')
title('Displacement')
xlabel('Tn(s)')
ylabel('Displacement(m)')
legend('ζ=0.02')
Grid
figure %绘制速度反应谱
plot(TA,MaxV(1,:),'-.b',TA,MaxV(2,:),'-r',TA,MaxV(3,:),':k')
title('Velocity')
xlabel('Tn(s)')
ylabel('velocity(m/s)')
legend('ζ=0.02')
Grid
figure %绘制绝对加速度反应谱
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('ζ=0.02')
Grid
3 运行的结果得到的反应谱
图2 位移反应谱
图3 速度反应谱
图4 加速度反应谱
一、 反应谱生成地震波
1 所取的反应谱为上海市
反应谱
图5 上海市设计反应谱
2 反应谱取值程序为:
%%
反应谱取值程序 参照01年抗震规范
function rs_z=r_s_1(pl,zn,ld,cd,fz) %%%pl 圆频率,zn阻尼比,ld烈度,cd场地类型,场地分组fz
%%%%烈度选择
if ld==6
arfmax=0.11;
end
if ld==7
arfmax=0.23;
end
if ld==8
arfmax=0.45;
end
if ld==9
arfmax=0.90;
end
%%%%场地类别,设计地震分组选择
if cd==1
if fz==1
Tg=0.25;
end
if fz==2
Tg=0.30;
end
if fz==3
Tg=0.35;
end
end
if cd==2
if fz==1
Tg=0.35;
end
if fz==2
Tg=0.40;
end
if fz==3
Tg=0.45;
end
end
if cd==3
if fz==1
Tg=0.45;
end
if fz==2
Tg=0.55;
end
if fz==3
Tg=0.65;
end
end
if cd==4
if fz==1
Tg=0.65;
end
if fz==2
Tg=0.75;
end
if fz==3
Tg=0.90;
end
end
%%%%%%%%%
ceita=zn; %%%%%阻尼比
lmt1=0.02+(0.05-ceita)/8;
if lmt1<0
lmt1=0;
end
lmt2=1+(0.05-ceita)/(0.06+1.7*ceita);
if lmt2<0.55
lmt2=0.55;
end
sjzs=0.9+(0.05-ceita)/(0.5+5*ceita);
%%%%%分段位置 T1 T2 T3
T1=0.1;
T2=Tg;
T3=5*Tg;
T_jg=2*pi./pl;
%%%% 第一段 0~T1
if T_jg<=T1
arf_jg=0.45*arfmax+(lmt2*arfmax-0.45*arfmax)/0.1*T_jg;
end
%%%% 第二段 T1~T2
if T1
参数