平面问题有限元分析.ppt

有限元分析,平面问题的有限单元法,第三章 平面问题的有限单元法,3-2、有限单元法的计算步骤3-3、单元位移函数3-4、单元载荷移置3-5、单元应力矩阵3-6、单元刚度矩阵3-7、单元刚度矩阵的物理意义及其性质3-8、整体分析3-9、整体刚度矩阵的形成3-10、支承条件的处理3-11、整体刚度矩阵的特点,2-2 有限单元法的计算步骤,弹性力学平面问题的有限单元法包括三个主要步骤：1、离散化 2、单元分析 3、单元综合,2-2 有限单元法的计算步骤,1、离散化 有限单元法的基础是用所谓有限个单元的集合体来代替原来的连续体，因而必须将连续体简化为由有限个单元组成的离散体。对于平面问题，最简单，因而最常用的单元是三角形单元。这些单元在结点处用铰相连，荷载也移置到结点上，成为结点荷载。在结点位移或其某一分量可以不计之处，就在结点上安置一个铰支座或相应的连杆支座。,5.相邻单元的尺寸尽可能接近。,6.结点所连接的单元个数尽可能一致。,有限元分析应注意的问题,7、充分利用结构的对称性,2-2 有限单元法的计算步骤,2、单元分析 对三角形单元，建立结点位移与结点力之间的转换关系。,结点位移,结点力,2-2 有限单元法的计算步骤,2、单元分析-单元刚度矩阵 取结点位移作基本未知量。由结点位移求结点力：其中，转换矩阵称为单元刚度矩阵。单元分析的主要目的就是要求出单元刚度矩阵。单元分析的步骤可表示如下：,2-2 有限单元法的计算步骤,3、单元综合 将离散化了的各个单元合成整体结构，利用结点平衡方程求出结点位移。在位移法中，主要的任务是求出基本未知量-结点位移。为此需要建立结点的平衡方程。,2-2 有限单元法的计算步骤,3、单元综合 i点总的结点力应为：根据结点的平衡条件，得 单元e的结点力，可按式(2-2)用结点位移表示，代入得到用结点位移表示的平衡方程。每个可动结点有两个未知位移，有两个平衡方程，所以方程总数与未知位移总数相等，可以求出所有的结点位移。单元综合的目的就是要求出结点位移。结点位移求出后，可进一步求出各单元的应力。,2-3 单元位移函数,如果弹性体的位移分量是座标的已知函数，则可用几何方程求应变分量，再从物理方程求应力分量。但对一个连续体，内部各点的位移变化情况很难用一个简单函数来描绘。有限单元法的基本原理是分块近似，即将弹性体划分成若干细小网格，在每一个单元范围内，内部各点的位移变化情况可近似地用简单函数来描绘。对每个单元，可以假定一个简单函数，用它近似表示该单元的位移。这个函数称为位移函数，或称为位移模式、位移模型、位移场。对于平面问题，单元位移函数可以用多项式表示，多项式中包含的项数越多，就越接近实际的位移分布，越精确。但选取多少项数，要受单元型式的限制。,2-2 单元位移函数,三结点三角形单元,六个节点位移只能确定六个多项式的系数，所以平面问题的3结点三角形单元的位移函数如下，所选用的这个位移函数，将单元内部任一点的位移定为座标的线性函数，位移模式很简单。,位移函数写成矩阵形式为：,将水平位移分量和结点坐标代入位移函数第一式，,写成矩阵形式，,令,则有,A为三角形单元的面积。,T的伴随矩阵为，,令,则,同样，将垂直位移分量与结点坐标代入位移插值公式,2-3 单元位移函数,最终确定六个待定系数,2-3 单元位移函数,令（下标i，j，m轮换）简写为,I是单位矩阵，N称为形态矩阵，Ni称为位移的形态函数,2-3 单元位移函数,选择单元位移函数时，应当保证有限元法解答的收敛性，即当网格逐渐加密时，有限元法的解答应当收敛于问题的正确解答。因此，选用的位移模式应当满足下列两方面的条件：(1)必须能反映单元的刚体位移和常量应变。6个参数 到 反映了三个刚体位移和三个常量应变。(2)必须保证相邻单元在公共边界处的位移连续性。(线性函数的特性),形态函数Ni具有以下性质：1)在单元结点上形态函数的值为1或为0。2）在单元中的任意一点上，三个形态函数之和等于1。,用 来计算三角形面积时，要注意单元结点的排列顺序，当三个结点i，j，m取逆时针顺序时，,当三个结点i，j，m取顺时针顺序时，,2-3 单元位移函数,作业：图示等腰三角形单元，求其形态矩阵N,应变矩阵,应力矩阵,单元刚度矩阵。,2-3 单元位移函数,由三角形的面积,2-4 单元载荷移置,连续弹性体离散为单元组合体时，为简化受力情况，需把弹性体承受的任意分布的载荷都向结点移置(分解)，而成为结点载荷。如果弹性体受承受的载荷全都是集中力，则将所有集中力的作用点取为结点，就不存在移置的问题，集中力就是结点载荷。但实际问题往往受有分布的面力和体力，都不可能只作用在结点上。因此，必须进行载荷移置。如果集中力的作用点未被取为结点，该集中力也要向结点移置。将载荷移置到结点上，必须遵循静力等效的原则。静力等效是指原载荷与结点载荷在任意虚位移上做的虚功相等。在一定的位移模式下，移置结果是唯一的，且总能符合静力等效原则。,单元的虚位移可以用结点的虚位移 表示为,表示为，,（3-15）令结点载荷为,令结点载荷为,集中力的移置如图所示，在单元内任意一点作用集中力,由虚功相等可得，,由于虚位移是任意的，则,例题1：在均质、等厚的三角形单元ijm的任意一点p（xp，yp）上作用有集中载荷。,体力的移置令单元所受的均匀分布体力为,由虚功相等可得，,分布面力的移置设在单元的边上分布有面力,，同样可以得到结点载荷，,例题：设有均质、等厚的三角形单元ijm，受到沿y方向的重力载荷qy的作用。求均布体力移置到各结点的载荷。,同理，,例题：在均质、等厚的三角形单元ijm的ij边上作用有沿x方向按三角形分布的载荷，求移置后的结点载荷。,取局部坐标s，在i点s=0，在j点s=l，L为ij边的长度。在ij边上，以局部坐标表示的插值函数为，,，,，,载荷为,例：,总载荷的2/3移置到结点i，1/3移置到结点j，与原载荷同向,2-5 单元应力矩阵,本节利用几何方程、物理方程，实现用结点位移表示单元的应变和单元的应力。用结点位移表示单元的应变的表达式为，B矩阵称为几何矩阵。,2-5 单元应力矩阵,由物理方程，可以得到单元的应力表达式 为应力矩阵,2-6 单元刚度矩阵,讨论单元内部的应力与单元的结点力的关系，导出用结点位移表示结点力的表达式。由应力推算结点力，需要利用平衡方程。第一章中已经用虚功方程表示出平衡方程。,2-6 单元刚度矩阵,考虑上图三角形单元的实际受力，结点力和内部应力为：任意虚设位移，结点位移与内部应变为,2-6 单元刚度矩阵,令实际受力状态在虚设位移上作虚功，外力虚功为,2-6 单元刚度矩阵,计算内力虚功时，从弹性体中截取微小矩形，边长为dx和dy，厚度为t，图示微小矩形的实际应力和虚设变形。,2-6 单元刚度矩阵,微小矩形的内力虚功为 整个弹性体的内力虚功为,2-6 单元刚度矩阵,根据虚功原理，得 这就是弹性平面问题的虚功方程，实质是外力与应力之间的平衡方程。虚应变可以由结点虚位移求出：代入虚功方程,2-6 单元刚度矩阵,接上式，将应力用结点位移表示出 有 令 则 建立了单元的结点力与结点位移之间的关系，称为单元刚度矩阵。它是6*6矩阵，其元素表示该单元的各结点沿坐标方向发生单位位移时引起的结点力，它决定于该单元的形状、大小、方位和弹性常数，而与单元的位置无关，即不随单元或坐标轴的平行移动而改变。,2-6 单元刚度矩阵,由于D中元素是常量，而在线性位移模式下，B中的元素也是常量，且 因此 可以进一步得出平面应力问题和平面应变问题中的单元刚度矩阵。,2-7 单元刚度矩阵的物理意义及其性质,已经求出了下列关系,2-7 单元刚度矩阵的物理意义及其性质,结点力和结点位移的关系：(以简单平面桁架为例)平面问题中，离散化的单元组合体极为相似，单元组合体在结点载荷的作用下，结点对单元、单元对结点都有作用力与反作用力存在，大小相等方向相反，统称为结点力。结点力和结点位移的关系前面已经求出：,2-7 单元刚度矩阵的物理意义及其性质,单元刚度矩阵的物理意义：将 写成分块矩阵 写成普通方程 其中 表示结点S(S=i,j,m)产生单位位移时，在结点r(r=i,j,m)上所需要施加的结点力的大小。,2-7 单元刚度矩阵的物理意义及其性质,单元刚度矩阵的物理意义：将结点力列矩阵 与结点位移列矩阵 均展开成(6*1)阶列矩阵，单元刚度矩阵相应地展开成(6*6)阶方阵：元素K的脚码，标有“-”的表示水平方向，没有标“-”的表示垂直方向。,2-7 单元刚度矩阵的物理意义及其性质,单元刚度矩阵的物理意义：单元刚度矩阵的每一个元素都有明显的物理意义。表示结点S(S=i,j,m)在水平方向、垂直方向产生单位位移时，在结点r(r=i,j,m)上分别所要施加的水平结点力和垂直结点力的大小。例如 表示结点j在垂直方向产生单位位移时，在结点i所需要施加的水平结点力的大小。,2-7 单元刚度矩阵的物理意义及其性质,单元刚度矩阵的性质：1)对称性：是对称矩阵 2)奇异性：是奇异矩阵，单元刚度矩阵所有奇数行的对应元素之和为零，所有偶数行的对应元素之和也为零。由此可见，单元刚度矩阵各列元素的总和为零。由对称性可知，各行元素的总和也为零。,2-7 单元刚度矩阵的物理意义及其性质,单元刚度矩阵的性质：例题：求下图所示单元的刚度矩阵，设,1、求B2、求D3、求S4、求,2-8 整体分析,将各单元组合成结构，进行整体分析。,整体分析分4个步骤1、建立整体刚度矩阵；2、根据支承条件修改整体刚度矩阵；3、解方程组，求出结点位移；(消去法与叠加法)4、根据结点位移求出应力。,2-8 整体分析,图示结构的网格共有四个单元和六个结点。在结点1、4、6共有四个支杆支承。结构的载荷已经转移为结点载荷。整体分析的四个步骤：1、建立整体刚度矩阵；2、根据支承条件修改整体刚度矩阵；3、解方程组，求结点位移；4、根据结点位移求出应力。,对单元的分析得出单元刚度矩阵，下面，将各单元组合成结构，进行整体分析。,2-8 整体分析,1、建立整体刚度矩阵(也叫作结构刚度矩阵)上图中的结构有六个结点，共有12个结点位移分量和12个结点力分量。由结构的结点位移向量求结构的结点力向量时，转换关系为：分块形式为：其中子向量 和 都是二阶向量，子矩阵 是二行二列矩阵。整体刚度矩阵K是12*12阶矩阵。,2-8 整体分析,2、根据支承条件修改整体刚度矩阵。建立整体刚度矩阵时，每个结点的位移当作未知量看待，没有考虑具体的支承情况，因此进行整体分析时还要针对支承条件加以处理。在上图的结构中，支承条件共有四个，即在结点1、4、6的四个支杆处相应位移已知为零：建立结点平衡方程时，应根据上述边界条件进行处理。3、解方程组，求出结点位移。通常采用消元法和迭代法两种方法。4、根据结点位移求出应力。,2-9 整体刚度矩阵的形式,整体刚度矩阵 是单元刚度矩阵 的集成。1、刚度集成法的物理概念：刚度矩阵中的元素是刚度系数，即由单位结点位移引起的结点力。由2-8节的例题可见，与结点2和3相关的单元有单元和，当结点3发生单位位移时，相关单元和同时在结点2引起结点力，将相关单元在结点2的结点力相加，就得出结构在结点2的结点力。由此看出，结构的刚度系数是相关单元的刚度系数的集成，结构刚度矩阵中的子块是相关单元的对应子块的集成。,2-9 整体刚度矩阵的形式,2、刚度矩阵的集成规则：先对每个单元求出单元刚度矩阵，然后将其中的每个子块 送到结构刚度矩阵中的对应位置上去，进行迭加之后即得出结构刚度矩阵K的子块，从而得出结构刚度矩阵K。关键是如何找出 中的子块在K中的对应位置。这需要了解单元中的结点编码与结构中的结点编码之间的对应关系。,2-9 整体刚度矩阵的形式,2、刚度矩阵的集成规则：,结构中的结点编码称为结点的总码，各个单元的三个结点又按逆时针方向编为i,j,m,称为结点的局部码。单元刚度矩阵中的子块是按结点的局部码排列的，而结构刚度矩阵中的子块是按结点的总码排列的。因此，在单元刚度矩阵中，把结点的局部码换成总码，并把其中的子块按照总码次序重新排列。,2-9 整体刚度矩阵的形式,以单元为例，局部码i,j,m对应于总码5,2,4，因此 中的子块按照总码重新排列后，得出扩大矩阵 为：,2-9 整体刚度矩阵的形式,用同样的方法可得出其他单元的扩大矩阵 将各单元的扩大矩阵迭加，即得出结构刚度矩阵K：集成规则包含搬家和迭加两个环节：1、将单元刚度矩阵 中的子块搬家，得出单元的扩大刚度矩阵。2、将各单元的扩大刚度矩阵 迭加，得出结构刚度矩阵K。(例题略),2-10 支承条件的处理,整体刚度矩阵K求出后，结构的结点力F可表示为 在无支杆的结点处，结点力就等于已知的结点载荷。在有支杆的结点处，则求结点力时，还应把未知的支杆反力考虑在内。如果用P表示结点载荷和支杆反力组成的向量，则结点的平衡方程为 根据支承条件对平衡方程加以处理。先考虑结点n有水平支杆的情况。与结点n水平方向对应的平衡方程是第2n-1个方程，根据支承情况，上式应换成，即在K中，第2n-1行的对角线元素 应改为1，该行全部非对角线元素应改为0。在P中，第2n-1个元素 应改为0。此外，为了保持矩阵K的对称性，则第2n-1列全部非对角线元素也改为0。,2-10 支承条件的处理,同理，如果结点n有竖向支杆，则平衡方程的第2n个方程应改为，为此，在矩阵K中，第2n行的对角线元素改为1，该行全部非对角线元素改为0，同时，第2n列全部非对角线元素也改为0。在P中，第2n个元素改为0。,2-10 支承条件的处理,2-8节中的结构，结点1有水平支杆，结点2有两个支杆，结点3有竖向支杆。对支承条件处理后，矩阵修改为：,2-10 支承条件的处理,最后考虑支点n的水平位移 为已知非零值 的情况，这时的支承条件为 对平衡方程的第2n-1个方程作如下修改：对角线系数 乘以一个大数A(例如A=1010)；右边自由项换成；其余各项保持不变，即有 实际上，上式中除包含大数A的两项外，其他各项相对地都比较小，可以忽略不计，因此与支承条件等价。对于支点的竖向位移 的支承条件，也可用同样的方法进行修改。,2-11 整体刚度矩阵的特点,在有限元法中，整体刚度矩阵的阶数通常是很高的，在解算时常遇到矩阵阶数高和存贮容量有限的矛盾。找到整体刚度矩阵的特性达到节省存贮容量的途径。1、对称性。只存贮矩阵的上三角部分，节省近一半的存贮容量。2、稀疏性。矩阵的绝大多数元素都是零，非零元素只占一小部分。,2-11 整体刚度矩阵的特点,2、稀疏性。矩阵的绝大多数元素都是零，非零元素只占一小部分。,结点5只与周围的六个结点(2、3、4、6、8、9)用三角形单元相连，它们是5的相关结点。只有当这七个相关结点产生位移时，才使该结点产生结点力，其余结点发生位移时并不在该结点处引起结点力。因此，在矩阵K中，第5行的非零子块只有七个(即与相关结点对应的七个子块)。,2-11 整体刚度矩阵的特点,2、稀疏性。,一般，一个结点的相关结点不会超过九个，如果网格中有200个结点，则一行中非零子块的个数与该行的子块总数相比不大于9/200，即在5%以下，如果网格的结点个数越多，则刚度矩阵的稀疏性就越突出。利用矩阵K的稀疏性，可设法只存贮非零元素，从而可大量地节省存贮容量。,2-11 整体刚度矩阵的特点,3、带形分布规律。上图中，矩阵K的非零元素分布在以对角线为中心的带形区域内，称为带形矩阵。在半个带形区域中(包括对角线元素在内)，每行具有的元素个数叫做半带宽，用d表示。半带宽的一般计算公式是：半带宽 d=(相邻结点码的最大差值+1)*2 上图中相邻结点码的最大差值为4，故d=(4+1)*2=10 利用带形矩阵的特点并利用对称性，可只存贮上半带的元素，叫半带存贮。,2-11 整体刚度矩阵的特点,图(a)中的矩阵K为n行n列矩阵，半带宽为d。半带存贮时从K中取出上半带元素，按图(b)中的矩阵 的排列方式进行存贮，即将上半部斜带换成竖带。存贮量n*d，存贮量与K中元素总数之比为d/n，d值越小，则存贮量约省。,2-11 整体刚度矩阵的特点,同一网格中，如果采用不同的结点编码，则相应的半带宽d也可能不同。如图，是同一网格的三种结点编码，相邻结点码的最大差值分别为4、6、8，半带宽分别为10、14、18。因此，应当采用合理的结点编码方式，以便得到最小的半带宽，从而节省存贮容量。,