newmark法程序法计算多自由度体系的动力响应(论文资料).docVIP

用matlab编程实现法计算多自由度体系的动力响应 用matlab编程实现法 计算多自由度体系的动力响应 一、法的基本原理 Newmark-?法是一种逐步积分的方法，避免了任何叠加的应用，能很好的适应非线性的反应分析。 Newmark-?法假定： (1-1) (1-2) 式中，?和?是按积分的精度和稳定性要求进行调整的参数。当?=0.5，?=0.25时，为常平均加速度法，即假定从t到t+?t时刻的速度不变，取为常数。研究表明，当?≥0.5，??≥0.25(0.5+?)2时，Newmark-?法是一种无条件稳定的格式。 由式(2-141)和式(2-142)可得到用及，，表示的，表达式，即有 (1-3) (1-4) 考虑t+?t时刻的振动微分方程为： (1-5) 将式(2-143)、式(2-144) 代入(2-145)，得到关于ut+?t的方程 (1-6) 式中 求解式(2-146)可得，然后由式(2-143)和式(2-144)可解出和。 由此，Newmark-?法的计算步骤如下： 1.初始计算： （1）形成刚度矩阵[K]、质量矩阵[M]和阻尼矩阵[C]； （2）给定初始值， 和； （3）选择积分步长?t、参数?、?，并计算积分常数，，，， ，，，； （4）形成有效刚度矩阵； 2.对每个时间步的计算： （1）计算t+?t时刻的有效荷载： （2）求解t+?t时刻的位移： （3）计算t+?t时刻的速度和加速度： Newmark-?方?t的大小不影响解的稳定性，?t的选择主要根据解的精度确定。 本文用法计算的基本问题 四层框架结构在顶部受一个简谐荷载的作用，力的作用时间=5s，计算响应的时间为100s，分2000步完成。阻尼矩阵由Rayleigh阻尼构造。 具体数据如下图： 图一：结构基本计算简图 计算法的源程序 clc clear 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=2*m*0.05*sqrt(k/m) ; t1=5; %力的作用时间 nt=2000; %分2000步完成 dt=0.01; %时间步长 alfa=0.25; %?=0.25 beta=0.5; %?=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 for i=2:nt t=(i-1)*dt; if (tt1), f=[300*sin(6*pi*t)-50*cos(3*pi*t);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