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

第第第第4 4 4 4章章章章 杆系结构的有限元分析原理杆系结构的有限元分析原理杆梁单元概述杆梁单元概述讨论杆梁单元和由它们组成的平面和空间杆梁结构系统讨论杆梁单元和由它们组成的平面和空间杆梁结构系统.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 1）用标准化的分段小单元来逼近原结构）用标准化的分段小单元来逼近原结构2 2）寻找能够满足位移边界条件的许可位移场）寻找能够满足位移边界条件的许可位移场3 3）基于位移场的最小势能原理来求解）基于位移场的最小势能原理来求解 基本变量为：基本变量为：节点位移内部各点位移应变应力(1)(3)(2)完整的求解过程1）离散化）离散化 该构件由两根杆件做成，因此可以自然离散成构件由两根杆件做成，因此可以自然离散成2个杆个杆单元。元。假定以假定以这类单元位移的特征元位移的特征为两个端点位移，就两个端点位移，就这两个离散两个离散单元元给出出节点点编号和号和单元元编号。号。单元元1：i=1，j=2 单元元2：i=2，j=32 2）单元分析）单元分析 单元位移模式：单元位移模式：u u(x x)=)=a a0 0+a a1 1x x 单元节点条件：单元节点条件：u u(0)=(0)=u u1 1，u u(1)=(1)=u u2 2 从而得：从而得：回代得回代得 写成矩阵形式为写成矩阵形式为其中其中Ni，Nj是形函数。是形函数。形函数矩阵形函数矩阵根据几何方程可得根据几何方程可得应变的表达的表达写成矩写成矩阵形式形式为简记为几何函数矩阵或者是应变转换矩阵几何函数矩阵或者是应变转换矩阵根据物理方程可得应力的表达根据物理方程可得应力的表达写成矩阵形式为写成矩阵形式为简记为简记为应力矩阵或者是应力转换矩阵应力矩阵或者是应力转换矩阵势能的表达势能的表达 写成矩阵形式为写成矩阵形式为刚度矩阵刚度矩阵节点力列阵节点力列阵3 3）离散单元的装配）离散单元的装配 在得到各个单元的势能表达式后，需要进行离散单元的装配，以求在得到各个单元的势能表达式后，需要进行离散单元的装配，以求出整个系统的总势能，对于该系统，总势能包括两个单元部分出整个系统的总势能，对于该系统，总势能包括两个单元部分4 4）边界条件的处理）边界条件的处理 处理边界条件是获取可能位移场，将左端的约束条件，处理边界条件是获取可能位移场，将左端的约束条件，即即u u1 1=0=0代入上式可以得到简化的势能表达式代入上式可以得到简化的势能表达式5）建立刚度方程）建立刚度方程 由于上式是基于许可位移场的表达的系统势能，这是由全由于上式是基于许可位移场的表达的系统势能，这是由全部节点位移分段所插值出的位移场为全场许位移场，且基本部节点位移分段所插值出的位移场为全场许位移场，且基本未知量为节点位移，根据最小势能原理（即针对未知位移求未知量为节点位移，根据最小势能原理（即针对未知位移求一阶导数）有一阶导数）有6）求解节点位移求解节点位移 将结构参数和外载荷代入上式有将结构参数和外载荷代入上式有求解得（单位求解得（单位m m）7 7）计算单元应变）计算单元应变8 8）计算单元应力）计算单元应力9 9）计算支反力）计算支反力 对于单元势能的表达，对其取极值有对于单元势能的表达，对其取极值有具体地对于单元具体地对于单元1 1，有，有其中其中R R1 1是节点是节点1 1的支反力，的支反力，P P2 2是单元是单元1 1的节点的节点2 2所受的力，即所受的力，即单元单元2 2对该节点的作用力，将前面求得的节点位移代入上式对该节点的作用力，将前面求得的节点位移代入上式可得支反力大小。可得支反力大小。以上是一个简单结构有限元方法求解得完整过程，对于以上是一个简单结构有限元方法求解得完整过程，对于复杂结构，其求解过程完全相同，由于每一个步骤都具备复杂结构，其求解过程完全相同，由于每一个步骤都具备标准化和规范性的特征，所以可以在计算机上编程而自动标准化和规范性的特征，所以可以在计算机上编程而自动实现。实现。讨论讨论1 1：对于一个单元的势能取极值，所得到的方程为：对于一个单元的势能取极值，所得到的方程为节点的位移和节点力之间的关系，也称为单元的平衡关系，节点的位移和节点力之间的关系，也称为单元的平衡关系，由此可以求出每一个单元所受的节点力。由此可以求出每一个单元所受的节点力。讨论讨论2 2：由前面的步骤，我们也可以直接将各个单元的刚：由前面的步骤，我们也可以直接将各个单元的刚度矩阵按照节点编号的对应位置来进行装配，即在未处理边度矩阵按照节点编号的对应位置来进行装配，即在未处理边界条件之前，先形成整体刚度矩阵。界条件之前，先形成整体刚度矩阵。其物理意义是，表示在未处理边界条件前的基于节点描述其物理意义是，表示在未处理边界条件前的基于节点描述的总体平衡关系。在对该方程进行位移边界条件的处理后就的总体平衡关系。在对该方程进行位移边界条件的处理后就可以求解，这样与先处理边界条件再求系统势能的最小值所可以求解，这样与先处理边界条件再求系统势能的最小值所获得的方程完全相同。获得的方程完全相同。4.2有限元分析的基本步骤及表达式有限元分析的基本步骤及表达式1、物体几何区域的离散化物体几何区域的离散化2 2、单元的研究（所有力学信息都用节点位移）来表达、单元的研究（所有力学信息都用节点位移）来表达3 3、装配集成、装配集成4 4、边界条件的处理并求解节点位移、边界条件的处理并求解节点位移5 5、支反力的求取以及其它力学量（应力、应变及位移三大、支反力的求取以及其它力学量（应力、应变及位移三大物理量）的计算物理量）的计算4.2有限元分析的基本步骤及表达式有限元分析的基本步骤及表达式4.34.3杆单元及其坐标变换杆单元及其坐标变换局部坐标系中的单元描述局部坐标系中的单元描述5.25m3.75m24mF6m3mF24mE=3E7pa=0.2836kg/m3F=100N变截面杆单元的推导变截面杆单元的推导单元的位移模式形状函数矩阵单元的几何矩阵变截面杆单元的推导变截面杆单元的推导单元刚度矩阵为4.3杆单元及其坐标变换杆单元及其坐标变换局部坐标系中的单元描述局部坐标系中的单元描述E=2E10paF=60kNA=250mm2150mm150mmF1.2mm4.3杆单元及其坐标变换-局部坐标o由于杆单元只有两个节点位移，故可以设杆单元的位移模式为之包含两个待定常数的形式 u(x)=a1+a2x根据有限元法的基本思路，将弹性体离散成有限个单元体的组合，以结点的位移作为未知量。弹性体内实际的位移分布可以用单元内的位移分布函数来分块近似地表示。在单元内的位移变化可以假定一个函数来表示，这个函数称为单元位移函数单元位移函数、或单元位移模式单元位移模式。o回代得 写成矩阵形式为其中Ni，Nj是形函数。o根据位移条件有u(0)=u0，u(l)=ul，从而得o根据几何方程得p根据物理方程得p从而，根据单元分析结果，进行整体分析，求解整体方程组，进行结果分析4.3.2杆单元的坐标变换 规定：规定：杆端位移和杆端力取在截面形心上，符号以与单元系坐标正向相同为正，相反为负。下面讨论整体坐标系下与局部坐标系下的转换关系式。整体坐标系单元杆端位移和杆端力仍定义在截面形心上，符号以与坐标正向同向为正反之为负。局部坐局部坐标标系系 整体坐整体坐标标系系4.3.2杆单元的坐标变换-平面问题其中是一个单位正交矩阵，单位正交矩阵的逆即等于其转置。从上图可以得出，整体坐标系逆针旋转角后与单元系相重合。写成矩阵形式为 由于单元的势能是一个标量(能量)，不会因坐标系的不同而改变，因此，可将节点位移的坐标变换关系代入原来基于局部坐标系的势能表达式中，整体坐标系下的刚度方程整体坐标系下的刚度方程根据根据得得其中其中单刚的性质：性质：是对称矩阵。是奇异矩阵。坐标变换并不改变矩阵的奇异性质。1 结构的离散化与编号 2各个单元的矩阵描述各个单元的矩阵描述 结构包括有斜杆，所以必须在总体坐标下对节点位移进行表达，所推导的单结构包括有斜杆，所以必须在总体坐标下对节点位移进行表达，所推导的单元刚度矩阵也要进行变换元刚度矩阵也要进行变换3 建立整体刚度方程 1.将所得到的各个单元刚度矩阵按节点编号进行组装，可以形成整体刚度矩阵;2.同时将所有节点载荷也进行组装。4 边界条件的处理及刚度方程求解 5 各单元应力的计算 6 支反力的计算 将节点位移的结果代入整体刚度方程中基于基于MATLAB平台求解该平台求解该（1）结构的离散化与编号结构的离散化与编号（2）计算各单元的刚度矩阵）计算各单元的刚度矩阵1.建立一个工作目录，将所编制的用于平面桁架单元分析的四个建立一个工作目录，将所编制的用于平面桁架单元分析的四个MATLAB函数函数（1.单元刚度；单元刚度；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/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(n1),DOF(n2)+k(n1,n2);endendz=KK;（3）建立整体刚度方程（4）边界条件的处理及刚度方程求解（高斯消去法）（5）支反力的计算（6）各单元的应力计算 基于基于MATLAB平台求解该平台求解该基于ANSYS求解该1.前处理2.求解器的设定3.后处理对对于于单单元元2 2：取：取i i=1=1，j j=2=2，则则 ，故，故 对于单元对于单元1：取：取i=3，j=1，则，则c=1，s=0，故，故 对于单元对于单元3：取：取i=2，j=3，则，则c=0，s=1，故，故 平面杆单元应用平面杆单元应用2 2整体编号，对号入座得总刚整体编号，对号入座得总刚杆单元的坐标变换-空间整体和局部的坐标转换关系与平面问题一致。ANSYS应用实例24.4梁单元及其坐标变换一般平面梁单元的描述 纯弯梁单元由于单元有四个位移分量，可设梁单元的位移模式v(x)为包含4个待定常数的三次多项式：由材料力学知,各截面的转角:根据边界条件可以确定待定系数，将其进一步回代，可以得到用节点位移表示的梁单元位移。式中单元的应力应变单元的应力应变 在弹性范围内，并且不考虑剪力的影响时，平面刚架单元内任一点的轴向线应变由两部分组成，即轴向应变与弯曲应变之和，其轴向应变与平面桁架轴向应变相同。轴向应变为 弯曲应变为 y为梁单元任意截面上任意点至中性轴(x轴)的距离。得出平面刚架单元应变 图3-5 弯曲应变计算示意图 则平面刚架梁单元的应变转换矩阵。根据梁的平面假定可知梁单元的轴向应变为：这里利用平面假设（这里利用平面假设（变形后横截面仍保持平面，与纵线正交）如图：从而可以由单向虎克定律得出单元的轴向应力：由虚功原理可以推得 组装总刚仍用后处理法，“对号入座，子块搬家”的方法。如：对于单元1，我们取i=1，j=2。故 对于单元2，取i=2，j=3。故由于I1=2I2=2I，按照“整体编号，对号入座”的原则，得总刚为 对于此，列出总刚度方程为考虑到边界条件，修正后的刚度方程为解之得 4.5平面刚架的有限元法小变形情况下，可以把平面刚架单元看成是发生轴向位移的杆单元和发生挠度和转角的梁单元的组合。将杆单元刚度矩阵与纯弯梁单元刚度矩阵进行组合可得到平面梁单元刚度矩阵。将杆单元刚度矩阵与纯弯梁单元刚度矩阵进行组合可得到平面梁单元刚度矩阵。4.5平面刚架的有限元法单元位移模式（1）平面桁架的单元位移模式（2）平面梁的单元位移模式其中：综合平面桁架和平面梁单元，得到平面刚架单元的单元位移模式。以下简记为 单元的应力和应变杆单元的轴向应变：梁单元的轴向应变：综合平面桁架和平面梁单元，得到平面刚架单元的应力和应变。简记为：局部坐标系下的单元刚度矩阵局部坐标系服从右手法则，考虑如图所示的典型单元。利用虚功原理得局部坐标系下的单刚，其中每个元素都有明确的物理意义 从上图可以得出，整体坐标系逆针旋转角后与单元系相重合。写成矩阵形式为 坐标变换（平面）其中T是一个单位正交矩阵，单位正交矩阵的逆即等于其转置。简记为从而得 对称矩阵。奇异矩阵。分块性质。坐标变换（空间）o纯轴向拉压o纯扭转oxoy面内弯曲oxoz面内弯曲将对应各部分刚度矩阵进行组合以完成完整的单元刚度矩阵o其中分别表示局部坐标轴对整体坐标轴的方向余弦。将所有的物理量写在一起，就得到梁单元的常用等效节点载荷梁单元的常用等效节点载荷 梁单元在承受非节点载荷下的节点载荷等效值，该等效值一般是根据外力功的计算公式得到的，因此，它与梁单元的边界条件没有关系。