东南大学 计算力学 中期作业
2)-2，3*Elem.Def(i，2)-1，3*Elem.Def(i，2)，3*Elem.Def(i，3)-2，3*Elem.Def(i，3)-1，3*Elem.Def(i，3)]); k = Beam_Elem_Stiff(Elem.Def(i，:)，flag); InFL=k*Displ; if ~isempty(Elem.Load) L=Elem.Def(i，5);p=Elem.Load(i，2); if flag==1; %经典梁 Pe=[0;-p*L/2 ;-p*L^2/12 ;0; -p*L/2 ;p*L^2/12]; else Pe=[0;-p*L/2 ;0 ;0; -p*L/2 ;0]; %timoshenko end InFL=InFL-Pe;%一定要注意求单元杆端力时，减去的是等效结点荷载列阵 end end function [ TotalLoad ] = PF_Load( Joint，Elem ，flag) %计算与集成结构综合荷载列阵 %包括集成：直接结点荷载，等效结点荷载 TotalLoad=zeros(3*Joint.NJoint，1); if flag==1 %经典梁的等效结点荷载 %把直接结点荷载集成到综合荷载列阵 if ~isempty(Joint.Load) for i=1:length(Joint.Load(:，1)) TotalLoad(3*Joint.Load(i，1)-3+Joint.Load(i，2))=Joint.Load(i，3); end end %计算单元等效结点荷载，并集成到综合荷载列阵中去 if ~isempty(Elem.Load) for i=1:length(Elem.Load(:，1)) L=Elem.Def(i，5);p=Elem.Load(i，2); Pe=[0;-p*L/2 ;-p*L^2/12 ;0; -p*L/2 ;p*L^2/12]; TotalLoad=AssemElem( Elem.Def(Elem.Load(i，1)，:)，Pe，TotalLoad ); end end else %Timosenko梁的等效结点荷载 %把直接结点荷载集成到综合荷载列阵 if ~isempty(Joint.Load) for i=1:length(Joint.Load(:，1)) TotalLoad(3*Joint.Load(i，1)-3+Joint.Load(i，2))=Joint.Load(i，
11
东南大学 计算力学 中期作业
3); end end %计算单元等效结点荷载，并集成到综合荷载列阵中去 if ~isempty(Elem.Load) for i=1:length(Elem.Load(:，1)) L=Elem.Def(i，5);p=Elem.Load(i，2); Pe=[0;-p*L/2 ;0 ;0; -p*L/2 ;0]; TotalLoad=AssemElem( Elem.Def(Elem.Load(i，1)，:)，Pe，TotalLoad ); end end end 在matlab命令窗口中输入process，以上程序运行结果：
分析其跨中挠度:?1、?2分别为经典梁理论和Timoshenko梁理论求得的中点挠度。用matlab计算出的结果与用ansys采用beam3和beam188的计算结果一致（这里忽略剪切系数的影响）。 计算结果对比： L/m?1/m ?3/m ?1/?3 ?2/m ?2/?3 0.2 0.3 0.5 0.7
-3.3028×10-7 -8.4107×10-7 -3.0884×10-6 -8.1379×10-6 -3.2719×10-7 -8.3380×10-7 -3.0664×10-6 -8.0912×10-6 -3.3010×10-7 -8.4058×10-7 -3.0862×10-6 -8.1312×10-6 1.0005 1.0006 1.0007 1.0008 0.9912 0.9919 0.9936 0.9951 12
东南大学 计算力学 中期作业
1 1.2 1.5 1.8 2 2.2 5 10 15 -2.5498×10-5 -4.7713×10-5 -1.0613×10-4 -2.0838×10-4 -3.0989×10-4 -4.4535×10-4 -1.1024×10-2 -1.7391×10-1 -8.7808×10-1 -2.5390×10-5 -4.7542×10-5 -1.0582×10-5 -2.0786×10-4 -3.0918×10-4 -4.4444×10-4 -1.1007×10-2 -1.7367×10-1 -8.769×10-1 -2.5474×10-5 -4.7664×10-5 -1.0602×10-4 -2.0815×10-4 -3.0954×10-4 -4.4484×10-4 1.11×10-2 1.737×10-1 8.77×10-1 1.0009 1.0010 1.0010 1.0011 1.0011 1.0011 0.9932 1.0012 1.0012 0.9967 0.9974 0.9981 0.9986 0.9988 0.9991 0.9996 0.9998 0.9999 （3）根据能量法求中点位移： ?k???kFQFQPFNFNPMMPds???ds???ds（*） EAGAEI实际荷载作用下的内力虚设单位荷载
px2Mp??Fx?M??x2FN?0 FNP?0FQP?px?FFQ?1利用（*）公式可求得滑动铰支座的反力F，再用同种方法可求得该超静定结构
中点挠度w。
取剪切系数k=1.2（这里k不随截面大小、跨度变化，为了手算简便），
(3.113L4?0.2728L2)?10?4w??(1.22L?0.1281)??L2(1.037L2?0.1921)?10?4 24.395L?0.28952代入L的长度，得到中点挠度?3，列于上表。
结论：
注：以下曲线图的横轴代表跨高比，此处定义为0.52/2/L。 ①
13
东南大学 计算力学 中期作业
?1/?3比值
从?1/?3比值都很接近1，用经典梁计算得到的结果很接近真实解，这是因为在经典梁基础上引入了剪切变形的影响。求出来的结果普遍比真实解略大，原因可能是用能量法手算时，只能得到一定位数的有效数，从而导致这种结果，所以上图中的曲线应该整体下移至红色曲线，最后一点应该越接近于1，这样结果就与已有结论（经典梁理论对长细比越大的梁越适用）相同。
高跨比越小，有限元解与真实解越接近，但是所有跨度中有限元解与真实解的差距都是可以忽略的。
②
?2/?3比值
?1/?3也都比较接近于1，Timoshenko由于一阶剪切变形理论的限制，只有
14
东南大学 计算力学 中期作业
适度的“粗”梁可以分析（这里的粗梁是指网格化后的梁单元的长度要小，故在用Timoshenko梁元是一定要划充分数量的单元）。梁的长细比GAL2/EI可以用来判定单元的适用性，其中，G为剪切模量，A 为横截面积，L为长度，EI为抗弯刚度。长细比的计算需要用全局尺寸，不是基于独立的单元尺度，这在上图中可以看出，L越大，曲线越接近于1。上述结论与ansys的beam188帮助文件中的说明一致。
③所有解都比1小，说明有限元解偏刚。
（二）通过实体单元分析梁截面的剪力滞效应——纵向弯曲 考虑闭口截面和开口截面两种情况。
单元采用shell181，采用映射法划分网格，综合考虑计算精度和计算机运行能力，单元边长不超过50mm，由于是采用完全积分（KEYOPTION=2），所以每个薄壁的厚度方向只需有一个单元。
梁两端简支（一定是将最下端的节点施加边界条件），面荷载P = 0.02N /mm2，弹性模量E = 2. 1×1011 Pa，泊松比u=0.3，跨度L=2000mm。 ①开口截面梁（以T型梁为例）：截面尺寸见下图
图.T型梁截面
15