用matlab编程实现法计算多自由度体系的动力响应样本.doc

用matlab编程实现法计算多自由度体系的动力响应

姓名:
学号:
班级:
专业:

用matlab编程实现法
计算多自由度体系的动力响应
一、 法的基本原理
Newmark-b法是一种逐步积分的方法, 避免了任何叠加的应用, 能很好的适应非线性的反应分析。
Newmark-b法假定:
(1-1)
(1-2)
式中, b和g是按积分的精度和稳定性要求进行调整的参数。当b=0.5, g=0.25时, 为常平均加速度法, 即假定从t到t+Dt时刻的速度不变, 取为常数。研究表明, 当b≥0.5, g≥0.25(0.5+b)2时, Newmark-b法是一种无条件稳定的格式。
由式(2-141)和式(2-142)可得到用及, , 表示的, 表示式, 即有
(1-3)
(1-4)
考虑t+Dt时刻的振动微分方程为:
(1-5)
将式(2-143)、 式(2-144) 代入(2-145), 得到关于ut+Dt的方程
(1-6)
式中

求解式(2-146)可得, 然后由式(2-143)和式(2-144)可解出和。
由此, Newmark-b法的计算步骤如下:
1.初始计算:
( 1) 形成刚度矩阵[K]、 质量矩阵[M]和阻尼矩阵[C];
( 2) 给定初始值, 和;
( 3) 选择积分步长Dt、 参数b、 g, 并计算积分常数
, , , ,
, , , ;
( 4) 形成有效刚度矩阵;
2.对每个时间步的计算:
( 1) 计算t+Dt时刻的有效荷载:

( 2) 求解t+Dt时刻的位移:

( 3) 计算t+Dt时刻的速度和加速度:

Newmark-b方法是一种无条件稳定的隐式积分格式, 时间步长Dt的大小不影响解的稳定性, Dt的选择主要根据解的精度确定。

二、 本文用法计算的基本问题
四层框架结构在顶部受一个简谐荷载的作用, 力的作用时间=5s, 计算响应的时间为100s, 分 步完成。阻尼矩阵由Rayleigh阻尼构造。
具体数据如下图:

图一: 结构基本计算简图
三、 计算法的源程序

m=[1,2,3,4];
m=diag(m);
k= [800 -800 0 0;
-800 2400 -1600 0;
0 -1600 4800 -3200;
0 0 -3200 8000];
c=0.05*m+0.02*k;
f0=100;
t1=5;
nt= ;
dt=0.01;

alfa=0.25;
beta=0.5;
a0=1/alfa/dt/dt;
a1=beta/alfa/dt;
a2=1/alfa/dt;
a3=1/2/alfa-1;
a4=beta/alfa-1;
a5=dt/2*(beta/alfa-2);
a6=dt*(1-beta);
a7=dt*beta;

d=zeros(4,nt);
v=zeros(4,nt);
a=zeros(4,nt);

for i=2:nt
t=(i-1)*dt;
if (t<t1), f=[f0*sin(4*pi*t/t1);0;0;0]; else f=[0;0;0;0]; end
ke=k+a0*m+a1*c;
fe=f+m*(a0*d(:,i-1)+a2*v(:,i-1)+a3*a(:,i-1))+c*(a1*d(:,i-1)+a4*v(:,i-1)+a5*a(:,i-1));

d(:,i)=inv(ke)*fe;
a(:,i)=a0*(d(:,i)-d(:,i-1))-a2*v(:,i-1)-a3*a(:,i-1);
v(:,i)=v(:,i-1)+a6*a(:,i-1)+a7*a(:,i);
end

四、 计算结果截图
最后程序分别计算出四个质点的位移、 速度、 加速度响应。
现将部分截图如下:

1、 位移响应:

图二: 1质点的位移响应

图三: 4质点的位移响应

2、 速度响应

图四: 1质点的速度响应

图五: 4质点的速度响应

3、 加速度响应

图六: 1质点的加速度响应

图七: 4质点的加速度响应