-
采用
EL-
CENTRO
地震波计算单自由体系的振动位移反应谱
前言
：本课程论文采用软件是
MATLAB R2010a
，所编程序应用的方法是线性加
速度法。
程序采用的是两层循环上实现位移谱的求取。
< br>内循环实现在地震波下由
上一步所得到的位移、速度和加速度得到下一时间间隔后
的位移速度和加速度；
外循环实现在不同结构自振周期
（频率）
下遭遇地震波的响应。
本文在求出相对
位移反应谱的同时，也求出了相对速度反应谱、绝对加速度反应谱。
一．线性加速度法的简述：
线性加速度法是直接数值积分法求解地震反应的方法之一，
本文所采用
的线性加
速度法参考大崎顺彦的
《地震动的谱分析入门》
第二版。
具体计算公式详见大崎
顺彦的《地震
动的谱分析入门》第二版
P116-P118
。
二．所编程序及编译：
clear
% ***********
读入地震记录
***********
fid = fopen(
''
);
[Accelerate,count] =
fscanf(fid,
'%g'
);
%count
读入的记录的量
time=0:0.02:(count-1)*0.02;
% ***********
线性加速度法计算各反应
***********
%
初始化各储存向量
Displace=zeros(1,count);
%
相对位移
Velocity=zeros(1,count);
%
相对速度
AbsAcce=zeros(1,count);
%
绝对加速
度
Damp=0.05;
%
结构阻
尼比取为
0.05
Tc=0.0:0.05:10;
%
结构自振周期
Dt=0.02;
%
地震记录的步长
< br>%
记录计算得到的反应，
MDis
为最大相对位移，
MVel
为最大相对速度
< br>%MAcc
某阻尼时最大
绝对加速度，用于画图
MDis=zeros(1,length(Tc));
MVel=zeros(1,length(Tc));
MAcc=zeros(1,length(Tc));
t=1;
%
在下一个循环中控制不同
的结构自振周期
for
T=0.0:0.05:10
Frcy=2*pi/T
%
结构自振频率
DamFrcy=Frcy*sqrt(1-Damp*Damp);
%
计算公式化简
e_t=exp(-Damp*Frcy*Dt);
s=sin(DamFrcy*Dt);
c=cos(DamFrcy*Dt);
A=zeros(2,2);
A(1,1)=e_t*(s*Damp/sqrt(1-D
amp*Damp)+c);
A(1,2)=e_t*s/DamFrcy;
A(2,1)=-Frcy*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)/(Frcy^2*Dt);
%
计算公式化简
d_3t=Damp/(Frcy^3*Dt);
B=zeros(2,2);
B(1,1)=e_t*((d
_f+Damp/Frcy)*s/DamFrcy+(2*d_3t+1/Frcy^2)*c)-2*d_3
t;
B(1,2)=-e_t*(d_f*s/DamFrcy+2*d_3t*c)
-1/Frcy^2+2*d_3t;
B(2,1)=-e_t*(((Damp/(
Frcy*Dt)+1)*s/DamFrcy)+(1/(Frcy^2*Dt))*c)+1/(Frcy^
2*Dt);
B(2,2)=e_t*((Damp/(Frcy*Dt)*s/Da
mFrcy)+(1/(Frcy^2*Dt))*c)-1/(Frcy^2*Dt);
for
i=1:(count-1)
%
根据地震记录
,
计算不同的反应
Dis
place(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)*Displa
ce(i)+A(2,2)*Velocity(i)+B(2,1)*Accelerate(i)+B(2,
2)*Accelerate(i+
1);
Abs
Acce(i+1)=-2*Damp*Frcy*Velocity(i+1)-Frcy^2*Displa
ce(i+1);
end
MDis(1,t)=max(abs(Displace));
MVel(1,t)=max(abs(Velocity));
if
T==0.0
MAcc(1,t)=max(abs(Accelerate));
%
当结构的自振周期为
0
时，
其绝对加速度应等于地面加速
度
else
MAcc(1,t)=max(abs(AbsAcce));
end
Displace=zero
s(1,count);
%
初始化各储存向量避免下次计算时引
用到前面结果
Velocity=zeros(1,count);
AbsAcce=zeros(1,count);
t=t+1;
%t=length
(Tc)
，即所求结构自振周期有多少个，对应就运行多少次。
end
%
***********
绘图
***********
figure
%
绘制位移反应谱
plot(Tc,MDis)
ti
tle(
'DisplacementResponseSpectrum'
)
xlabel(
'Tn(s)'
) <
/p>
ylabel(
'Displacement(m)'
)
legend(
'ζ=0
.05'
)
grid
figure
%
绘制速度反应谱
plot(Tc,MVel)
ti
tle(
'VelocityResponseSpectrum'
)
xlabel(
'Tn(s)'
) <
/p>
ylabel(
'Velocity(m)'
)
legend(
'ζ=0.05'
)
grid
figure
%
绘制速度反应谱
plot(Tc,MAcc)
ti
tle(
'absAccelerateResponseSpectrum'
)
xlabel(
'Tn(s)'
) <
/p>
ylabel(
'absAccelerate(m)'
)
legend(
'ζ=0.05'
)
grid
三．地震反应谱：
1.
相对位移反应谱（
DisplacementResponseSpectr
um
）