第4章-杆梁结构的有限元分析原理.课件.ppt

1、杆梁单元概述讨论杆梁单元和由它们组成的平面和空间杆梁结构系统.o从构造上来说其长度远大于其截面尺寸的一维构件o承受轴力或扭矩的杆件成为杆o杆梁问题都有精确解o承受横向力和弯矩的杆件称为梁o平面桁架 平面刚架 连续梁 空间刚架 空间桁架等o变截面杆和弯曲杆件本章主要内容o4.1有限元分析的完整过程o4.2有限元分析的基本步骤及表达式o4.3杆单元及其坐标变换o4.4梁单元及其坐标变换4.1有限元分析的完整过程E1=E2=2E7PaA1=A2=2cm2l1=l2=10cmP3为10N作用下二杆结构的变形。o问题的解题思路：1）用标准化的分段小单元来逼近原结构2）寻找能够满足位移边界条件的许可位移场

2、3）基于位移场的最小势能原理来求解 基本变量为：节点位移内部各点位移应变应力(1)(3)(2)完整的求解过程1）离散化 该构件由两根杆件做成，因此可以自然离散成2个杆单元。假定以这类单元位移的特征为两个端点位移，就这两个离散单元给出节点编号和单元编号。 单元1：i=1，j=2 单元2：i=2，j=32）单元分析 单元位移模式：u(x)=a0+a1x 单元节点条件：u(0)=u1， u(l)=u2 从而得01,1,2jiieuuauaijl回代得 写成矩阵形式为01( )1jiieijeeiijju xaa xuuuxlxxuullN uN uiiujujeuqNNuNq其中Ni，Nj是形函数。

3、形函数矩阵形函数矩阵根据几何方程可得应变的表达写成矩阵形式为简记为1d1dxjieuauuxl11 1iiiujuejjuuNNuuleBq几何函数矩阵或者是应变转换矩阵几何函数矩阵或者是应变转换矩阵根据物理方程可得应力的表达写成矩阵形式为简记为ddxjieuEEuuxleSq应力矩阵或者是应力转换矩阵应力矩阵或者是应力转换矩阵1 1iiiujuejjuuEE NNuul节点位移列阵势能的表达 1 122T1 1220TT1 1220TT1 1221d21d21d212eeeeeeijijleeeleeeeee eUWPuPuBqSqAxPuPuqB EBq AxPuPuqB EBq A lP

4、uPu 写成矩阵形式为TT1 1221112122211121222TT1211111 1121212eeee ee eeeeeeeeeeeeeeeeqB EBq A lPuPuuuuuEA lPPuullEAEAuulluuPPuuEAEAllqK qP q 刚度矩阵刚度矩阵节点力列阵节点力列阵3）离散单元的装配 在得到各个单元的势能表达式后，需要进行离散单元的装配，以求出整个系统的总势能，对于该系统，总势能包括两个单元部分 121T112T221T12T21122112222111212331122332211221111111231112110022012eq K qqK qP qPqE

5、AEAEAEAuuuulllluuRuuFuuuuEAEAEAEAllllEAEAllEAEAuuull 112221322222332200uuEAEAuRFulluuEAEAll4）边界条件的处理 处理边界条件是获取可能位移场，将左端的约束条件，即u1=0代入上式可以得到简化的势能表达式 121T112T221T12T21221222223322332212102eq K qqK qP qPqEAEAEAuullluuFuuEAEAll 5）建立刚度方程 由于上式是基于许可位移场的表达的系统势能，这是由全部节点位移分段所插值出的位移场为全场许位移场，且基本未知量为节点位移，根据最小势能原理

6、（即针对未知位移求一阶导数）有12212222233220EAEAEAullluFEAEAll6）求解节点位移 将结构参数和外载荷代入上式有求解得（单位m）2222222332230EAEAulluFEAEAll233102E41110uu232.5E47.5E4uu7）计算单元应变1111211 12.5E3iiujujuNNuuul22222311 15E3iiujujuNNuuul8）计算单元应力111121 10.05MpaiiujujuE NNuuEul2222231 10.1MpaiiujujuE NNuuEul9）计算支反力 对于单元势能的表达，对其取极值有具体地对于单元1，有其

7、中R1是节点1的支反力，P2是单元1的节点2所受的力，即单元2对该节点的作用力，将前面求得的节点位移代入上式可得支反力大小。eeeK qP1111221111uREAuPl以上是一个简单结构有限元方法求解得完整过程，对于复杂结构，其求解过程完全相同，由于每一个步骤都具备标准化和规范性的特征，所以可以在计算机上编程而自动实现。讨论1：对于一个单元的势能取极值，所得到的方程为节点的位移和节点力之间的关系，也称为单元的平衡关系，由此可以求出每一个单元所受的节点力。讨论2：由前面的步骤，我们也可以直接将各个单元的刚度矩阵按照节点编号的对应位置来进行装配，即在未处理边界条件之前，先形成整体刚度矩阵。其物

8、理意义是，表示在未处理边界条件前的基于节点描述的总体平衡关系。在对该方程进行位移边界条件的处理后就可以求解，这样与先处理边界条件再求系统势能的最小值所获得的方程完全相同。KqP4.2有限元分析的基本步骤及表达式1、物体几何区域的离散化2、单元的研究（所有力学信息都用节点位移）来表达3、装配集成4、边界条件的处理并求解节点位移5、支反力的求取以及其它力学量（应力、应变及位移三大物理量）的计算eBqeqNqeSqeeeK qP4.2有限元分析的基本步骤及表达式4.3杆单元及其坐标变换4.3.1局部坐标系中的单元描述5.25m3.75m24mF6m3mF24mE=3E7pa=0.2836kg/m3F

9、=100N变截面杆单元的推导变截面杆单元的推导1d1dxjieuauuxl 1212121( )TeTexxA xAAllquuPPP节点位移列阵节点力列阵单元的位移模式形状函数矩阵eBq单元的几何矩阵变截面杆单元的推导变截面杆单元的推导单元刚度矩阵为4.3杆单元及其坐标变换4.3.1局部坐标系中的单元描述E=2E10paF=60kNA=250mm2150mm150mmF1.2mm4.3杆单元及其坐标变换-局部坐标 i j xyo由于杆单元只有两个节点位移，故可以设杆单元的位移模式为之包含两个待定常数的形式 u(x)=a1+a2x根据有限元法的基本思路，将弹性体离散成有限个单元体的组合，以结点

10、的位移作为未知量。弹性体内实际的位移分布可以用单元内的位移分布函数来分块近似地表示。在单元内的位移变化可以假定一个函数来表示，这个函数称为单元位移函数单元位移函数、或单元位移模式单元位移模式。 o回代得 写成矩阵形式为12( )1jiiijiuijuju xaa xuuuxlxxuullN uN uiiujujuuNNu其中Ni，Nj是形函数。o根据位移条件有u(0)=u0， u(l)=ul，从而得12,1,2jiiuuauaijlo根据几何方程得p根据物理方程得p从而，根据单元分析结果，进行整体分析，求解整体方程组，进行结果分析2d1dxjiuauuxl11 1iiiujujjuuNNuul

11、ddxjiuEEuuxl1 1iiiujujjuuEE NNuul4.3.2杆单元的坐标变换 i j xy i j xy 规定：规定：杆端位移和杆端力取在截面形心上，符号以与单元系坐标正向相同为正，相反为负。下面讨论整体坐标系下与局部坐标系下的转换关系式。整体坐标系单元杆端位移和杆端力仍定义在截面形心上，符号以与坐标正向同向为正反之为负。局部坐标系局部坐标系 整体坐标系整体坐标系4.3.2杆单元的坐标变换-平面问题其中是一个单位正交矩阵，单位正交矩阵的逆即等于其转置。 o R A x y xyB uuv从上图可以得出，整体坐标系逆针旋转角后与单元系相重合。cossincossincossinc

12、ossiniiijjjiiijjjuuvuuvXXYXXY写成矩阵形式为 eeeeqTqPTP0000csTcs 由于单元的势能是一个标量(能量)，不会因坐标系的不同而改变，因此，可将节点位移的坐标变换关系代入原来基于局部坐标系的势能表达式中， 整体坐标系下的刚度方程整体坐标系下的刚度方程根据根据eeeeqTqPTP得得其中其中22222222-eeeeTeeeccsccscsscssE AkT k Tlccsccscsscss单刚的性质：性质：是对称矩阵。是奇异矩阵。坐标变换并不改变矩阵的奇异性质。 eeTeeTeeeTeeeePTPTk qTk Tqk q1 结构的离散化与编号 2各个单元

13、的矩阵描述 结构包括有斜杆，所以必须在总体坐标下对节点位移进行表达，所推导的单元刚度矩阵也要进行变换22222222-eeeeTeeeccsccscsscssE AkT k Tlccsccscsscss3 建立整体刚度方程 1.将所得到的各个单元刚度矩阵按节点编号进行组装，可以形成整体刚度矩阵;2.同时将所有节点载荷也进行组装。 4 边界条件的处理及刚度方程求解 5 各单元应力的计算 6 支反力的计算 将节点位移的结果代入整体刚度方程中基于MATLAB平台求解该（1） 结构的离散化与编号 （2）计算各单元的刚度矩阵1.建立一个工作目录，将所编制的用于平面桁架单元分析的四个MATLAB函数（1.

14、单元刚度；2.总刚矩阵的组装；3.单元应力的求解；4.支反力的求解）2.在MATLAB环境中，输入弹性模量E、横截面积A，各点坐标、角度3.调用四次单元刚度矩阵计算函数，得到各个单元的刚度矩阵单元的刚度矩阵的计算function k=Bar2D2Node_Stiffness(E,A,x1,y1,x2,y2,alpha)%该函数计算单元的刚度矩阵%输入弹性模量E，横截面积A%输入第一个节点坐标（x1,y1），第二个节点坐标（x2,y2），角度alpha（单位是度）%输出单元刚度矩阵k(4X4)。%-L=sqrt(x2-x1)*(x2-x1)+(y2-y1)*(y2-y1);x=alpha*pi/

15、180;C=cos(x);S=sin(x);k=E*A/L*C*C C*S -C*C -C*S; C*S S*S -C*S -S*S; -C*C -C*S C*C C*S; -C*S -S*S C*S S*S;总刚矩阵的组装function z = Bar2D2Node_Assembly(KK,k,i,j)%该函数进行单元刚度矩阵的组装%输入单元刚度矩阵k，单元的节点编号i、j%输出整体刚度矩阵KK%-DOF(1)=2*i-1;DOF(2)=2*i;DOF(3)=2*j-1;DOF(4)=2*j;for n1=1:4 for n2=1:4 KK(DOF(n1),DOF(n2)= KK(DOF(

16、n1),DOF(n2)+k(n1,n2); endendz=KK;（3） 建立整体刚度方程 （4） 边界条件的处理及刚度方程求解（高斯消去法）（5）支反力的计算 （6）各单元的应力计算 基于MATLAB平台求解该基于ANSYS求解该1.前处理2.求解器的设定3.后处理对于单元2：取i=1，j=2，则 ，故 x 1 2 3 1 1 p 45o y 对于单元1：取i=3，j=1，则c=1，s=0，故 111333111111131110-1 00000-1 0100000kkE Aklkk22,22sc22211122222212221-1 -11-111-1-111-12 21-1 -11kkE

17、 Akkkl对于单元3：取i=2，j=3，则c=0，s=1，故 333333233333232200-0-001-0-10-000-0-101kkE Aklkk平面杆单元应用2整体编号，对号入座得总刚1 22111121322 33212223131 331323311111102 22 22 22 21111002 22 22 22 21111002 22 22 22 211111012 22 22 22 2100010000101kkkKkkkkkk11112222333311111102 22 22 22 21111002 22 22 22 21111002 22 22 22 21111

18、1012 22 22 22 2100010000101xyxyxyFuFvFuFvFuFv11112233111102 22 22 211102 22 22 2111112 22 22 20001xyyyFuFvFvFv 杆单元的坐标变换-空间12Tequu111222Tequvwuvw212121cos,cos,cos,xxx xlyyx ylzzx zl11122212cos,cos,cos,000000cos,cos,cos,eeeux xx yx zqux xx yxuvTuzwwqv整体和局部的坐标转换关系与平面问题一致。eeeeeqT qPTPANSYS应用实例24.4梁单元及其坐

19、标变换一般平面梁单元的描述 纯弯梁单元由于单元有四个位移分量，可设梁单元的位移模式v(x)为包含4个待定常数的三次多项式：2301232123( )( )23v xaa xa xa xxaa xa x由材料力学知,各截面的转角:xv 根据边界条件可以确定待定系数，将其进一步回代，可以得到用节点位移表示的梁单元位移。 ( )v xN q 1234( )iijjvv xNNNNNqv式中232323231234232232132,2,32,xxxxxxxxNNxNNllllllll 单元的应力应变单元的应力应变 在弹性范围内，并且不考虑剪力的影响时，平面刚架单元内任一点的轴向线应变由两部分组成，即

20、轴向应变与弯曲应变之和，其轴向应变与平面桁架轴向应变相同。 轴向应变为 弯曲应变为 y为梁单元任意截面上任意点至中性轴(x轴)的距离。 得出平面刚架单元应变 xulx22xvybx图3-5 弯曲应变计算示意图 22xvyxubxlxx exB则平面刚架梁单元的应变转换矩阵。 B exxBEE根据梁的平面假定可知梁单元的轴向应变为： ( )yv xy NqBq这里利用平面假设（这里利用平面假设（变形后横截面仍保持平面，与纵线正交）如图： 从而可以由单向虎克定律得出单元的轴向应力：EEy NSq 123436124661226ySSSSSElxlx llxx ll 由虚功原理可以推得 zzzz32

21、32zzzz22zzzz3232zzzz221261266462d1261266264eTVEIEIEIEIllllEIEIEIEIllllKB EB VEIEIEIEIllllEIEIEIEIlllleeePK q组装总刚仍用后处理法，“对号入座，子块搬家”的方法。如：对于单元1，我们取i=1，j=2。故 1111323211111221112121221111323211112212612664621261266264EIEIEIEIllllEIEIEIEIkkllllkkkEIEIEIEIllllEIEIEIEIllll对于单元2，取i=2，j=3。故222232322222222222

22、2323222222212612664621261266264EIEIEIEIllllEIEIEIEIllllkEIEIEIEIllllEIEIEIEIllll由于I1=2I2=2I，按照“整体编号，对号入座”的原则，得总刚为 32322232323222232322224122412001281240024123661261246126212612600626400EIEIEIEIllllEIEIEIEIllllEIEIEIEIEIEIllllllKEIEIEIEIEIEIllllllEIEIEIEIllllEIEIEIEIllll对于此，列出总刚度方程为323222323232222323

23、22224122412001281240024123661261246126212612600626400EIEIEIEIllllEIEIEIEIllllEIEIEIEIEIEIllllllEIEIEIEIEIEIllllllEIEIEIEIllllEIEIEIEIllll111122223333yyyvFMvFMvFM考虑到边界条件，修正后的刚度方程为32322232323222232322224122412001281240024123661261246126212612600626400EIEIEIEIllllEIEIEIEIllllEIEIEIEIEIEIllllllEIEIEIEIE

24、IEIllllllEIEIEIEIllllEIEIEIEIllll111122223333yyyvFMuFMuFM32222366-6120EIEIvpllEIEIll 解之得 32223366plvEIplEI4.5平面刚架的有限元法小变形情况下，可以把平面刚架单元看成是发生轴向位移的杆单元和发生挠度和转角的梁单元的组合。将杆单元刚度矩阵与纯弯梁单元刚度矩阵进行组合可得到平面梁单元刚度矩阵。将杆单元刚度矩阵与纯弯梁单元刚度矩阵进行组合可得到平面梁单元刚度矩阵。4.5平面刚架的有限元法单元位移模式（1）平面桁架的单元位移模式（2）平面梁的单元位移模式iiujujuuNNuiiivijvjivi

25、jvjjjvNNNNvNNNNv 其中： 1,iujuxxNNll 232323232322321 32,2,32,ivijvjxxxxxxxxNNxNNllllllll 综合平面桁架和平面梁单元，得到平面刚架单元的单元位移模式。00000000iiiujuiivijvjjjvijvjjjuvuNNvNNNNuNNNNv 以下简记为 eqNq单元的应力和应变1( )xu x2( )xyv x杆单元的轴向应变：梁单元的轴向应变：综合平面桁架和平面梁单元，得到平面刚架单元的应力和应变。12232232161226161226Txxxiiijjjxxxxyyyyuvuvllllllllll 简记为：

26、eBqeeEEBqSq局部坐标系下的单元刚度矩阵局部坐标系服从右手法则，考虑如图所示的典型单元。 利用虚功原理得局部坐标系下的单刚， 其中每个元素都有明确的物理意义 xyzijlE A I323222323222000012612600646200000012612600626400eeiieeiieeiiejejejEAEAllEIEIEIEIUullllVvEIEIEIEIMllllUEAEAllVEIEIEIEIMllllEIEIEIEIllllejejejuvo R A x y xyB uuvv从上图可以得出，整体坐标系逆针旋转角后与单元系相重合。cossinsincoscossins

27、incosiiiiiiiijjjjjjjjuuvvuvuuvuuv写成矩阵形式为 eeeeqTqPTP坐标变换（平面）其中T是一个单位正交矩阵，单位正交矩阵的逆即等于其转置。 0000000000100000000000000001csscTcssc简记为eeqTqeeePT P从而得 对称矩阵。 奇异矩阵。 分块性质。 TeeeekT k T坐标变换（空间）111111222222Texyzxyzquvwuvw111222111222TeuvwxyzuvwxyzPPPPMMMPPPMMMo纯轴向拉压12,1111u uEAklo纯扭转12,1111GlkI oxoy面内弯曲112222z,3

28、2212612664621261266264vzvzllllllEIkllllllloxoz面内弯曲112222,32212612664621261266264wywyyllEIllllklllllll将对应各部分刚度矩阵进行组合以完成完整的单元刚度矩阵111111222222Texyzxyzquvwuvw111222111222TeuvwxyzuvwxyzPPPPMMMPPPMMM111111111111111cos( , )cos( , )cos( , )cos( , )cos( , )cos( , )cos( , )cos( , )cos( , )uux xux yux zuvvy xv

29、y yvy zvwwz xwz ywz zw111111111111111cos( , )cos( , )cos( , )cos( , )cos( , )cos( , )cos( , )cos( , )cos( , )xxyzxyxyzyzxyzzx xx yx zy xy yy zz xz yz z o其中分别表示局部坐标轴对整体坐标轴的方向余弦。将所有的物理量写在一起，就得到cos( , ),cos( , ).cos( , )x xx yz zeeeeeqT qPTP12 12000000000000eT梁单元的常用等效节点载荷梁单元的常用等效节点载荷 梁单元在承受非节点载荷下的节点载荷等效值，该等效值一般是根据外力功的计算公式得到的，因此，它与梁单元的边界条件没有关系。