newmark法程序法计算多自由度体系的动力响应_精品文档Word格式文档下载.doc
《newmark法程序法计算多自由度体系的动力响应_精品文档Word格式文档下载.doc》由会员分享，可在线阅读，更多相关《newmark法程序法计算多自由度体系的动力响应_精品文档Word格式文档下载.doc（9页珍藏版）》请在冰豆网上搜索。
由此，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，分2000步完成。
阻尼矩阵由Rayleigh阻尼构造。
具体数据如下图：
图一：
结构基本计算简图
三、计算法的源程序
clc
clear
m=[1,2,3,4];
m=diag（m）;
%计算质量矩阵
k=[800-80000;
-8002400-16000;
0-16004800-3200;
00-32008000];
%计算刚度矩阵
c=2*m*0.05*sqrt（k/m）;
t1=5;
%力的作用时间
nt=2000;
%分2000步完成
dt=0.01;
%时间步长
alfa=0.25;
%g=0.25
beta=0.5;
%b=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）;
%初位移为0
v=zeros（4,nt）;
%初速度为0
a=zeros（4,nt）;
%初加速度为0
fori=2:
nt
t=（i-1）*dt;
if（t<
t1）,
f=[300*sin（6*pi*t）-50*cos（3*pi*t）;
0;
0];
%力作用时间内对结构进行加载
else
f=[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））;
%t+Dt时刻的有效荷载
d（:
i）=inv（ke）*fe;
%求解t+Dt时刻的位移
a（:
i）=a0*（d（:
i）-d（:
i-1））-a2*v（:
i-1）-a3*a（:
i-1）;
%计算t+Dt时刻的加速度
v（:
i）=v（:
i-1）+a6*a（:
i-1）+a7*a（:
i）;
%计算t+Dt时刻的速度
end
T=[0:
dt:
19.99];
%离散系统dt为采样周期19.99为终端时间
closeall
figure%控制窗口数量
plot（T,[d]）%绘制位移函数图像
title（'
¸
÷
Ö
Ê
µ
ã
Î
»
Ò
Æ
×
Ü
Í
¼
'
）%添加标题为各质点位移总图
legend（'
'
¶
þ
È
ý
Ë
Ä
）%添加图例的标注
xlabel（'
±
ä
（s）'
）%对x轴进行标注为时间（s）
ylabel（'
（m）'
）%对y轴进行标注为位移（m）
grid%显示画图中的个网线
figure
plot（T,[a]）%绘制加速度函数图像
Ó
Ù
）%添加标题为各质点加速度总图
）
）%对x轴进行标注为时间（s）
（m/s^2）'
）%对y轴进行标注为加速度（m/s2）
grid
%
figure
plot（T,[v]）%绘制速度函数图像
）%添加标题为各质点速度总图
）%对x轴进行标注为时间（s）
）%对y轴进行标注为加速度（m/s2）
四、计算结果截图
最后程序分别计算出四个质点的位移、速度、加速度响应。
现将部分截图如下：
1、位移响应：
图二：
1质点的位移响应
图三：
4质点的位移响应
2、速度响应
图四：
1质点的速度响应
图五：
4质点的速度响应
3、加速度响应
图六：
1质点的加速度响应
图七：
4质点的加速度响应
9
升级会员 