《结构动力学大作业(共10页).docx》由会员分享,可在线阅读,更多相关《结构动力学大作业(共10页).docx(10页珍藏版)》请在taowenge.com淘文阁网|工程机械CAD图纸|机械工程制图|CAD装配图下载|SolidWorks_CaTia_CAD_UG_PROE_设计图分享下载上搜索。
1、精选优质文档-倾情为你奉上结构动力学大作业班级: 学号:姓名: 目 录1. Wilson-法原理简介图1-1 Wilson-法示意图Wilson-法是基于对加速度a的插值近似得到的,图1-1为Wilson-法的原理示意图。推导由t时刻的状态求t+t时刻的状态的递推公式: (1-1)对积分可得速度与位移的表达式如下: (1-2) (1-3)其中=t,由式(1-2)、(1-3)可以解出: (1-4) (1-5)将式(1-4)、(1-5)带入运动方程: (1-6) (1-7)注意到此时的式子为和上一个时刻、以及t+t时刻的荷载相关,可以运用迭代的思想来求解,下图给出线弹性条件下Wilson-法的流程
2、图:图1-2 Wilson-法流程图2. Wilson-程序验算对线弹性条件下的Wilson-法进行MATLAB编程,源代码见附录。选取如下算例进行验证。对于一个单自由度的无阻尼结构,当其受到一个周期荷载时,其结构响应分为稳态解和瞬态解,由于没有阻尼的影响,其瞬态解并不会衰减,其理论表达式为: (2-1)式中,为位移响应,为激励,为刚度,为荷载频率与固有振动频率之比,为荷载频率,为结构固有频率。现令为1,为1,则为1,取为2/3。程序求得的解与解析解对比如图2-1所示(由于理论解与程序基本重合,所以将理论解乘以-1,方便比较): a)位移 b)速度c)加速度图2-1 Wilson-法结果验证2
3、.1 t的影响上述算例验证时选择的t非常小,因此看不出理论解与Wilson-法的求解区别,以下改变t的取值,探讨t对迭代的影响。图2-2 t对位移曲线的影响可以看出并不是t太大时计算结果很不准确,偏小,反映不出周期特征;当t合适时正好基本和理论解重合,也不是t越小越好越小时越能反映出一些细部特征,但这也不是很准确。2.2 的影响当1.37时,该算法是无条件稳定的算法,以下探讨对算法的影响。图2-3 对位移曲线的影响由上图可知随着值越大,位移的周期变大。3. 非线性问题求解由于实际结构并不一定为线性,其刚度会随着位移的的变化而改变,下图为求解非线性问题时的Wilson-法流程。此处要说明的是,刚
4、度矩阵Ky(t)是与位移相关的量,判断那时候速度的大小是为了确定其是否处于卸载段。具体可能得根据实际情况求解。图3-1 Wilson-法解非线性问题修改MATLAB程序,并用该程序来计算如下例题:对该问题采用Wilson-法非线性方式计算,采用t=0.1s和t=0.05s两种方式,计算位移、速度和加速度曲线如下图所示: a)位移 b)速度 c)加速度图3-2 非线性分析结果由上图可知,结构在0.6s时达到位移极值,在t=0.1s和t=0.05s算得的值分别为0.096m和0.108m,速度极值在0.9s取到分别为-0.468和-0.580,加速度极值在t=0.1s时为0.7s时取到,为-2.1
5、27,在t=0.05s极值在0.75s时取到,为-2.9。4. 附录Wilson-法源程序function y_1,y_2,y_3=wilson_theta(p,m,c,k,dt,v0,y0,a_0,theta)%p代表输入的荷载,c为阻尼矩阵,dt为时间间隔,m为质量矩阵,k为刚度矩阵%v0为初始的速度,y0为初始的位移.a_0为初始加速度%输出的矩阵y_1代表位移,y_2代表速度,y_3代表加速度if nargin9 theta=1.4;endL,r=size(p);y_1=NaN(L,r);y_2=NaN(L,r);y_3=NaN(L,r);y_1(:,1)=y0;y_2(:,1)=v0
6、;y_3(:,1)=a_0;%计算积分常数a0=6/(theta*dt)2);a1=3/theta/dt;a2=2*a1;a3=theta*dt/2;a4=a0/theta;a5=-a2/theta;a6=1-3/theta;a7=dt/2;a8=dt2/6;%计算拟刚度矩阵k0=k+a0*m+a1*c;%计算拟荷载for i=1:r-1 R=p(:,i)+theta*(p(:,i+1)-p(:,i)+m*(a0*y_1(:,i)+a2*y_2(:,i)+2*y_3(:,i)+c*(a1*y_1(:,i)+2*y_2(:,i)+a3*y_3(:,i); y_theta=k0R; y_3(:,i
7、+1)=a4*(y_theta-y_1(:,i)+a5*y_2(:,i)+a6*y_3(:,i); y_2(:,i+1)=y_2(:,i)+a7*(y_3(:,i)+y_3(:,i+1); y_1(:,i+1)=y_1(:,i)+y_2(:,i)*dt+a8*(y_3(:,i+1)+2*y_3(:,i);end上述代码只适合分析线弹性结构,对于非线性结构,编写起来比较繁琐,针对不同的情况可能需要具体处理,所以本文只给出了针对本文例题的代码。与上述代码不同的是以下代码增加了一个判断的语句。function y_1,y_2,y_3=wilson_theta2(p,m,c,dt,v0,y0,a_0,
8、theta)%p代表输入的荷载,c为阻尼矩阵,dt为时间间隔,m为质量矩阵%v0为初始的速度,y0为初始的位移.a_0为初始加速度%输出的矩阵y_1代表位移,y_2代表速度,y_3代表加速度if nargin0 k=60*(y_1(:,i)0.05); else k=60; end %计算拟刚度矩阵 k0=k+a0*m+a1*c; %计算拟荷载 R=p(:,i)+theta*(p(:,i+1)-p(:,i)+m*(a0*y_1(:,i)+a2*y_2(:,i)+2*y_3(:,i)+c*(a1*y_1(:,i)+2*y_2(:,i)+a3*y_3(:,i); y_theta=k0R; y_3(:,i+1)=a4*(y_theta-y_1(:,i)+a5*y_2(:,i)+a6*y_3(:,i); y_2(:,i+1)=y_2(:,i)+a7*(y_3(:,i)+y_3(:,i+1); y_1(:,i+1)=y_1(:,i)+y_2(:,i)*dt+a8*(y_3(:,i+1)+2*y_3(:,i);end专心-专注-专业