角形单元有限元

单击此处编辑母版标题样式,单击此处编辑母版文本样式,第二级,第三级,第四级,第五级,*,结构有限元法,第1章 三角形常应变单元的有限元法,第2章 有限元程序设计与分析软件,第3章 平面问题高阶单元的有限元法,第4章 空间实体的有限元法,第5章 杆系结构的有限元法,第6章 板壳问题的有限元法,第7章 结构动力问题的有限元法？,第8章 弹塑性问题的有限元法,结 构 有 限 元 分 析,第1章 三角形单元的有限元法,1.1 有限元法的基本思想,有限元法在20世纪50年代起源于飞机结构的矩阵分析，其,基本思想,是用有限个离散单元的集合体代替原连续体，采用能量原理研究单元及其离散集合体的平衡，以计算机为工具进行结构数值分析。它避免了经典弹性力学获得连续解的困难（建立和求解偏微分方程），使大型、复杂结构的计算容易地在计算机上完成，应用十分广泛。ANSYS, SAP2K,把整体结构离散为有限个单元，研究单元的平衡和变形协调；再把这有限个离散单元集合还原成结构，研究离散结构的平衡和变形协调。,划分的单元大小和数目根据计算精度和计算机能力来确定。,1,2,3,4,5,6,7,8,9,10,P,5,7,6,4,5,6,3,4,5,6,7,8, , ,弹性悬臂板,剖分与集合,单元、节点需编号,有限元法主要优点：,（1） 概念浅显，容易掌握。（离散、插值、能量原理、数学分析）,（2）适用性强，应用范围广，几乎适用于所有连续体和场问题的分析。（结构、热、流体、电磁场和声学等问题）,（3）计算规格化（采用矩阵表示），便于计算机编程。,1.1.1 有限元法的分析步骤,（1）结构,离散化：,用点、线或面把结构剖分为有限个离散单元体，并在单元指定点设置,节点,。研究单元的平衡和变形协调，形成单元平衡方程。,l,/2,l,/2,P,1,2,3,1,、F,1,2,、F,2,3,、F,3,4,、F,4,l,/2,1,2,l,/2,2,3,1,、F,1,2,、F,2,3,、F,3,4,、F,4,单元的节点上有位移,和力F,（2）单元集合：把所有离散的有限个单元集合起来代替原结构，形成离散结构节点平衡方程。,（3）由平衡方程求解得节点位移和计算单元应力。,1,、F,1,2,、F,2,3,、F,3,4,、F,4,l,/2,1,2,l,/2,2,3,1,、F,1,2,、F,2,3,、F,3,4,、F,4,l,/2,l,/2,P,1,2,3,1.1.2 有限元法分析思路流程,解综合方程K= P,求结构节点位移,计算结构内力和应力,系统分析,(把单元刚度矩阵集合成结构刚度矩阵K,形成等价节点荷载P ),离散（剖分）结构,为若干单元,单元分析,(建立单元刚度矩阵k,e,形成单元等价节点力),（1-1）,2、单元内任意点的,体积力,列阵,q,V,（1-2）,1、单元表面或边界上任意点的,表面力,列阵,q,s,i,j,m,x,y,i,j,m,x,y,q,V,q,s,1.2 基本力学量矩阵表示,图1-1,i,j,m,x,y,u,v,3、单元内任意点的位移列阵,f,（1-3）,4、单元内任意点的应变列阵,（1-4）,i,j,m,x,y,5、单元内任意点的应力列阵,（1-5）,6、,几何方程,（1-6）,将上式代入式(1-4），,i,j,m,x,y,（1-4）,7、物理方程矩阵式,（1-7）,式中 E、,弹性模量、泊松比。,上式可简写为,（1-8）,其中,对于弹性力学的平面应力问题，,物理方程,的矩阵形式可表示为：,(1-9）,矩阵D称为弹性矩阵。,对于平面应变问题，将式(1-9）中的E换为 ，,换为 。,(1-8）,各种类型结构的弹性物理方程都可用式(1-8）描述。但结构类型不同，力学性态 (应力分量、应变分量)有区别， 弹性矩阵D的体积和元素是不同的。,1.3 位移函数和形函数,1,、位移函数概念,由于有限元法采用能量原理进行单元分析，因而必须事先设定位移函数。 “位移函数”也称 “位移模式”，是,单元内部位移变化的数学表达式，设为坐标的函数,。,一般而论，位移函数选取会影响甚至严重影响计算结果的精度。在弹性力学中，恰当选取位移函数不是一件容易的事情；但,在有限元中，当单元划分得足够小时，把位移函数设定为简单的多项式就可以获得相当好的精确度。,这正是有限单元法具有的重要优势之一。,不同类型结构会有不同的位移函数。这里，仍以平面问题三角形单元（图1-2）为例，说明设定位移函数的有关问题。,图1-2是一个三节点三角形单元，其节点,i,、,j,、,m,按,逆时针,方向排列。每个节点位移在单元平面内有两个分量：,（1-10）,一个三角形单元有3个节点（以,i,、,j,、,m为,序），共有6个节点位移分量。其,单元位移或单元节点位移列阵,为：,图1-2,i,j,m,u,i,u,j,u,m,v,i,v,j,v,m,x,y,2、位移函数设定,本问题选位移函数（单元中任意一点的位移与节点位移的关系）为简单多项式：,(1-12）,式中：,a,1,、,a,2,、,a,6,待定常数，由单元位移的,6个分量确定。,a,1,、,a,4,代表刚体位移，,a,2,、,a,3,、,a,5,、,a,6,代表单元中的常应变，而且，位移函数是连续函数。,(1-11）,i,j,m,u,i,u,j,u,m,v,i,v,j,v,m,x,y,u,v,选取位移函数应考虑的问题,（1）,位移函数的个数,等于,单元中任意一点的位移分量个数。本单元中有,u,和,v,，与此相应，有2个位移函数；,（3）,位移函数中待定常数个数,待定常数个数应等于,单元节点自由度总数,，以便用单元节点位移确定位移函数中的待定常数。本单元有6个节点自由度，两个位移函数中共包含6个待定常数。,（2）,位移函数是坐标的函数,本单元的坐标系为：,x,、,y,；,（4）,位移函数中必须包含单元的刚体位移。,（5）,位移函数中必须包含单元的常应变。,（6）,位移函数在单元内要连续。相邻单元间要尽,量协调。,条件（4）、（5）构成单元的,完备性,准则。,条件（6）是单元的位移,协调性,条件。,理论和实践都已证明，完备性准则是有限元解收敛于真实解的必要条件。单元的位移协调条件构成有限元解收敛于真实解的充分条件。,容易证明，三角形三节点常应变单元满足以上必要与充分条件。,（7）,位移函数的形式,一般选为完全多项式。为实现（4）（6）的要求，根据Pascal三角形由低阶到高阶按顺序、对称地选取；多项式的项数等于（或稍大于）单元节点自由度数。,例：平面应力矩形板被划分为若干三角形单元。,位移函数中包含了单元的常应变。,(a,2, a,6, a,3,+a,5,),位移函数中包含了单元的刚体位移。,(a,1, a,4,),2,5,4,1,3,6,对任一单元，如单元，取位移函数：,、单元的位移函数都是,可以看出：,位移函数在单元内是连续的；,以、的边界26为例,2,5,6,2,6,3,5,6,2,3,x,y,u,u,6,u,2,u,u,6,u,2,两条直线上有两个点重合，此两条直线必全重合。,位移函数在单元之间的边界上也连续吗？是。,3、形函数,形函数是用单元节点位移分量来描述位移函数的插值函数,。,(1-13）,（1）形函数确定,现在，通过单元节点位移确定位移函数中的待定常数a,1,、a,2,、a,6,。设节点,i,、,j,、,m,的坐标分别为（,x,i,、,y,i,）、（,x,j,、,y,j,）、（,x,m,、,y,m,），节点位移分别为（,u,i,、,v,i,）、 （,u,j,、,v,j,） 、 （,u,m,、,v,m,）。将它们代入式（1-12），有,从式(1-13）左边3个方程中解出待定系数a,1,、a,2,、a,3,为,(1-14）,式中， A为三角形单元的面积，有,(1-15）,特别指出：,为使求得面积的值为正值，本单元节点号的次序必须是,逆时针,转向，如图所示。至于将哪个节点作为起始节点,i,，则没有关系。,将式(1-14）代入式(1-12）的第一式，整理后得,同理,i,j,m,x,y,（2）,（1）,（7）,(1-16）,式中,(1-17）,i,j,m,式(1-17）中（,i,、,j,、,m,）意指：按,i,、,j,、,m,依次轮换下标，可得到,a,j,、,b,j,、,c,j,a,m,、,b,m,、,c,m,。后面出现类似情况时，照此推理。式(1-17）表明：,a,j,、,b,j,、,c,j,a,m,、,b,m,、,c,m,是单元三个节点坐标的函数。,(1-16）,令,(1-18）,位移模式(1-16）可以简写为,(1-19）,式(1-19）中的,N,i,、,N,j,、,N,m,是坐标的函数，反应了单元的位移形态，称为单元位移函数的形函数。数学上它反应了节点位移对单元内任一点位移的插值，又称插值函数。,(1-16）,用形函数把式(1-16）写成矩阵，有,缩写为,(1-20）,形函数是有限单元法中的一个重要函数，它具有以下性质：,N,为形函数矩阵，写成分块形式：,(1-21）,其中子矩阵,(1-22）,I,是22的单位矩阵。,（2）形函数性质,性质1,形函数,N,i,在节点,i,上的值等于1，在其它节点,上的值等于0。对于本单元，有,（,i,、,j,、,m,）,性质2,在单元中任一点，所有形函数之和等于1。对,于本单元，有,x,y,N,（i，j，m）,N,i,=1,i,j,m,图1-3,？,x,y,N,（I，j，m）,N,i,=1,i,j,m,N,j,=1,i,j,m,N,m,=1,i,j,m,N,i,=1,i,j,m,N,j,=1,N,m,=1,图1-4,也可利用行列式代数余子式与某行或列元素乘积的性质（等于行列式值或0）证明。,性质3,在三角形单元的边界,ij,上任一点（,x,，,y,），有,x,x,i,x,j,x,y,N,i,(,x,i，,y,i,),j,(,x,j，,y,j,),m,(,x,m，,y,m,),N,i,(,x、y,),1,证,图1-5,（1）,性质4 形函数在单元上的面积分和在边界上的线积分公式为,(1-23）,式中 为 边的长度。,1.4 单元应变和应力,根据几何方程(1-6）和位移函数(1-16）可以求得单元应变。,1、单元应变,(1-6）,对位移函数（式(1-16）,(1-24）,(1-16）,求导后代入式(1-6），得到应变和节点位移的关系式。,上式简写一般式：,(1-25）,式中， ,B,单元应变矩阵。,对本问题，维数为36。它的分块形式为：,子矩阵,(1-26）,由于 与,x,、,y,无关，都是常量，因此B矩阵也是常量。单元中任一点的应变分量是,B,矩阵与单元位移的乘积，因而也都是常量。因此，这种单元被称为常应变单元。,2、单元应力,将式(1-25）代入物理方程式(1-8），得 单元应力,(1-27）,也可写为,(1-28）,其中：,S,称为,单元应力矩阵,，并有,(1-29）,这里，,D,是33矩阵，,B,是36矩阵，因此,S,也是36矩阵。它可写为分块形式,(1-30）,将弹性矩阵（式(1-9） 和应变矩阵（式(1-26）代入，得子矩阵,S,i,由式(1-29）,(1-31）,式(1-31）是平面应力的结果。对于平面应变问题，只要将上式中的E换成 ，,换成 即得。,(1-32）,由于同一单元中的,D,、,B,矩阵都是常数矩阵，所以,S,矩阵也是常数矩阵。也就是说，三角形三节点单元内的应力分量也是常量。,当然，相邻单元的b,i,、c,i,(,i,，,j,，,m,)一般不完全相同，因而具有不同的应力，这就造成在相邻单元的公共边上存在着,应力突变现象,。但是随着网格的细分，这种突变将会迅速减小，收敛于平衡被满足。,1.5 单元平衡方程,1、 单元应变能,对于平面应力问题中的三角形单元，设单元厚度为h 。,将式(1-25）和(1-8）代入上式进行矩阵运算，并注意到弹性矩阵,D,的对称性，有,应变能,U,为,i,j,m,x,y,h,(1-25),(1-8）,由于,和,T,是常量，提到积分号外，上式可写成,引入矩阵符号,k,，且有,(1-33a）,式(1-33a）是针对平面问题三角形单元推出的。注意到其中,h,d,x,d,y,的实质是任意的微体积d,v,，于是得计算k的,一般式,。,(1-33）,式(1-33）不仅适合于平面问题三角形单元，也是计算各种类型单元,k,的一般式。,d,v,1.6节中将明确,k,的力学意义是单元刚度矩阵。式(1-33）便是计算单元刚度矩阵的基本矩阵式。它适合于各种类型的单元。,单元应变能写成,(1-34）,2、 单元外力势能,单元受到的外力一般包括,体积力,、,表面力,和,集中力,。自重属于体积力范畴。表面力指作用在单元表面的分布载荷，如风力、压力，以及相邻单元互相作用的内力等。,(1-33）,（1） 体积力势能,单位体积中的体积力,如式(1-35）所示。,单元上体积力具有的势能,V,v,为,(1-35）,i,j,m,x,y,q,V,x,q,V,y,i,j,m,x,y,u,v,注意到式(1-20）,有,(1-20）,（2） 表面力势能,面积力虽然包括单元之间公共边上互相作用的分布力，但它们属于结构内力，成对出现，集合时互相抵消，在结构整体分析时可以不加考虑，因此单元分析时也就不予考虑。,现在，只考虑弹性体边界上的表面力，它只在部分单元上形成表面力（右下图）。设,边界面上单位面积受到的表面力,如下式：,l,单元边界长度,h,单元厚度,A,表面力作用面积,q,s,q,s,沿厚度均匀分布，,则单元表面力的势能V,s,为,（3） 集中力势能,当结构受到集中力时，通常在划分单元网格时就把集中力的作用点设置为节点。于是单元集中力,P,c,的势能,V,c,为,p,p/2,C,（4）,总势能,把(1-35）式中原括符内的部分用列阵,F,d,代替，,综合以上诸式，单元外力的总势能,V,为,(1-35）,F,d,具有和,相同的行、列数。则,(1-36）,由单元的应变能,U,(1-34）和外力势能,V,(1-36），可得单元的总势能,(1-37）,将式(1-37）代入，,根据弹性力学最小势能原理：结构处于,稳定平衡,的必要和充分条件是,总势能有极小值,。,3、单元平衡方程,于是有，,(1-34）,(1-36）,式(1-38）是从能量原理导出的单元平衡方程。这个方程表达了单元力与单元位移之间的关系。其中，,F,d,和单元节点力,F,具有相同的意义。,(1-38）,即得单元平衡方程,1.6 单元刚度矩阵,平衡方程(1-38）中的矩阵,k,是单元力和单元位移关系间的系数矩阵，代表了单元的刚度特性，称为单元刚度矩阵。单元刚度矩阵的体积为,n,j,n,j,，,n,j,是单元位移总数。其一般计算公式为：,1、一般计算公式,它与单元应变矩阵,B,和弹性矩阵,D,有关。,(1-33）,对于平面应力三角形单元，应变矩阵,B,是常数矩阵，同时弹性矩阵,D,也是常数矩阵，于是式(1-33）可以化简为,式中,A,表示三角形单元的面积。,h,是单元厚度。,2、平面问题三角形单元刚度矩阵,（1）平面应力三角形单元,(1-39）,将式(1-9）和(1-26）代入上式，,即得平面应力三角形单元刚度矩阵。写成分块形式，有,(1-40）,(1-9）,(1-26）,式(1-40)中子矩阵,k,rs,为22矩阵，有,(1-41）,（2）平面应变三角形单元,对于平面应变问题，须将上式中的E换为 ， 换为,，于是有,，组合见式(1-40）,b,i(j,m),、c,i(j,m),是形函数式(1-16）中的系数(式2-17)。,(1-42）,平面问题的单元刚度矩阵,k,不随单元（或坐标轴）的平行移动或作,n,角度（,n,为整数）的转动而改变。,由公式(1-41）、(1-42）知，,k,rs,矩阵和其中的,b,r,、,c,r,、,b,s,、,c,s,（,r,、,s,=,i,、,j,、,m,）有关。, 单元平移时，,b,i,、,c,i,不变。,，组合见式(1-40）,(3) 三角形单元刚度矩阵与坐标系无关,i,j,m,x,y,o, 单元转动时，,b,i,、,c,i,不变。,当单元旋转时，各节点的编号保持不变。如图1-7所示，图,a,所示的单元旋转,时，到达图,b,所示位置。,(1-,17）,i,j,m,y,j,y,m,i,j,m,图1-7,x,y,o,(b),x,y,o,(a),j,i,m,可以证明，这两种情形的,k,是相同的。,其实，推演公式(1-40）、(1-41）、(1-42）时并没有规定坐标系的方位，当坐标系旋转任意角度时，也不影响刚度矩阵的结果。因此，,平面问题的单元刚度矩阵可以认为是结构坐标系中的单元刚度矩阵，没有坐标变换问题。,(1-38）,（1）单元刚度矩阵中每个元素有明确的物理意义,例如，,k,ij,表示单元第,j,个自由度产生单位位移(,j,=1)，,其他自由度固定（=0）时，在第,i,个自由度产生的节点力F,i,。,主对角线上元素,k,ii,(,i,=1,n,j,)恒为正值。,3、单元刚度矩阵性质,（2）,k,的每一行或每一列元素之和为零,F,1,=0,F,2,=0,F,3,=0,F,i,=0,F,j,=0,F,nj,=0,r,s,t,1,1,以上式中第,i,行为例，,当所有节点沿,x,向或,y,向,都产生单位位移时，,单元作平动运动，无应变，也无应力。,则有：,即：,k,的每一行元素之和为零。根据对称性，每一列元素之和也为零。,r,s,t,x,y,图1-6,（3）,k,是对称矩阵,由,k,各元素,的表达式，可知k具有对称性。,n,j,n,j,对于主对角线元素对称。对称表达式：,k,ij,=,k,ji,证明,k,ij,表示当单元位移中第,j,个元素为1(,j,=1),其余元素为零时，引起的单元力中的第,i,个节点力F,i,k,ji,表示当单元位移中第,i,个元素为1(,i,=1),其余元素为零时，引起的单元力中的第,j,个节点力F,j,第,i,自由度,第,j,自由度,位移,i,=1,j,=1,力,F,i,=,k,ij,F,j,=,k,ji,虚功,F,i,i,=,k,ij,F,j,j,=,k,ji,由虚功原理，得,k,ij,=,k,ji,（4）单元刚度矩阵是奇异矩阵,即,k,的行列式为零（由行列式性质） 。,单元刚度矩阵是在单元处于平衡状态的前提下得出的。单元作为分离体看待，作用在它上面的外力（单元力）必定是平衡力系。然而，研究单元平衡时没有引入约束。,承受平衡力系作用的无约束单元，其变形是确定的，但位移不是确定的。,所以出现性质（3）中的“平动问题”，即单元可以发生任意的刚体运动。从数学上讲，方程(1-28）的解不是唯一的或不能确定的。由此，单元刚度矩阵一定是奇异的。,（5）单元刚度矩阵是常量矩阵,单元力和单元位移成线性关系是基于弹性理论的结果。,4、例：平面应力直角三角形单元刚度矩阵,图1-8示出一平面应力直角三角形单元，直角边长分别为,a,、,b,，厚度为,h,，弹性模量为,E,，泊松比为,，计算单元刚度矩阵。,图1-8,i,j,m,a,b,x,y,第一步：计算,b,i,、,c,i,和单元,面积A。,图1-8,（1-17）,i,j,m,a,b,x,y,x,i(j,m),y,i(j,m),b,i(j,m),c,i(j,m),i,a,0,b,0,j,0,b,0,a,m,0,0,-,b,-a,表2-1 单元节点坐标和,b,i,、,c,i,值（,i,、,j,、,m,）,参数,节点,单元面积:,A=ab/2, 计算步骤,第二步：求子矩阵,由式(1-41），算得,其他从略。,第三步：形成,k,将,k,ii,等按式(1-40）组集成,k, 。,(1-43a）,2i-1 2i 2j-1 2j 2m-1 2m,2i-1,2i,2j-1,2j,2m-1,2m,红色号码,是单元位移（,1,、,2,、）在结构中对应的节点位移的序号。,i,j,m,i,j,m,i,、,j,、,m,表示单元中3个节点在结构系统中的编号。,当a=b时，即等腰直角三角形单元，有,(1-43b）,i j m,i,j,m,1.7 等价节点力,从前面单元分析可以看出：单元平衡所用到的的量均要属于节点的量，如单元位移、单元力。载荷亦应如此，必须将体积力、表面力转化到节点上去，成为等价节点力（载荷）。在第2.5节中已经得到了公式(1-35）和(1-36） 。,这里，,F,d,就是体积力、表面力和集中力之和的总等价节点力。,(1-35）,(1-36）,(1-44）,把总等价节点力,F,d,分解成体积力、表面力和集中力的等价节点力之和，有,F,V,单元上体积力的等价节点力,F,S,单元上表面力的等价节点力,p,C,单元上节点上的集中力,注意到式(1-35），得体积力等价节点力计算公式：,表面力的等价节点力计算公式：,(1-45）,(1-46）,1、体积力的等价节点力,2、表面力的等价节点力,3、等价节点力计算举例,（1）单元自重,图1-9所示平面应力三角形单元，单元厚度为h。单元单位体积自重为,，自重指向y轴的负方向。,P,vix,P,viy,P,vjx,P,vjy,P,vmx,P,vmy,(1-45）, 计算式,(1-21）,图1-9,x,y,i,j,m,-,注意到形函数的性质4：,(1-23）,得自重荷载的等价节点力,(1-22）,（,i,j,m,）,根据体积力和式(1-45）、(1-21）、(1-22），得,(1-47）,上式表明：自重载荷的等价节点力为单元重量的1/3。,（2）均布面力,i,j,m,图1-10,x,y,q,s,单元边界上作用了均匀的分布力，如图1-10所示，其集度为,q,s,。,(1-46）,(1-21）,根据式(1-46）、(1-21）和(1-22）, 计算式,注意到形函数性质4 ：,(1-23）,得,(1-48）,(1-22）,均匀分布力的等价节点力为,式(1-48）表明：在,ij,边上受均布面力的平面问题三角形单元，其等价节点力等于将均布面力合力之半简单地简化到,i,、,j,节点上，方向与分布力方向相同。,m,节点上为零。,(1-48）,i,j,m,x,y,q,sx,F,s,1,F,s,3,i,j,m,x,y,q,sy,F,s,2,F,s4,（3）线性分布面力,i,j,m,图1-11,x,y,s,表面力集度在,i,点为,q,sx,q,sy,T,，而在,j,点为0。设坐标轴s的原点取在,j,点，沿,ji,为正向， 。,ij,边上任一点的面力集度,q,s,s,q,si,q,s,i,j,m,图1-12,x,y,s,l,在,ij,边上有：,将,q,s,和上式代入式(1-46），有,由形函数的性质3：,(1-49）,式(1-49）表明：,ij,边受线性分布面力：,i,点为q,sx, q,sy,T,，,j,点为0,时，其等价节点力可将总载荷的2/3分配给,i,点，1/3分配给,j,点，,m,点为零得出。,x,y,i,j,m,q,si,q,s,体积力和表面力向节点的移置符合静力等效原理的前提条件是：线性位移模式。,1.7 系统分析,1.7.1 坐标系,研究各离散单元集合成整体结构，集合整体结构的平衡和变形协调，建立整体结构平衡方程。,单元分析时采用的坐标系成为局部坐标或单元坐标（单元刚度矩阵的通用性）。而结构系统分析时，必须在统一的坐标系内进行（各力学量才能叠加），,称为“结构坐标”或“整体坐标”,，如图1-13所示。,单元坐标系下，单元位移、单元力、单元刚度矩阵表示为：,整体坐标系下，单元位移、单元力、单元刚度矩阵表示为：,X,Y,X,Y,P,P,图（1-13）,(a) 平面桁架(杆件单元),悬臂深梁,(,平面三角形单元,),x,y,x,y,x,y,x,y,如何从单元坐标转化为结构坐标将在第4章中讨论。,1.7.2 整体刚度矩阵,假设整体结构被划分为,n,e,个单元和,n,个节点，在整体坐标系下，对于每个单元均有：,将上述这些方程集合起来（整体坐标下叠加），便可得到整个结构的平衡方程。为此，需要将,k,、,、,F,体积膨胀，分别扩大为,n1,n1,、,n1,1,和,n1,1,的矩阵才能相加。膨胀后，原有节点号对应位置的元素不变，而其它元素均为零。,组装方法：,建立一个体积为n1n1的方阵，按单元序号依次把结构坐标单元刚度矩阵的元素放入该方阵中。,放入方法：,（1）按单元节点编码,对号入座,； （2）,同位置元素累加,。,式中：,K,为整体刚度矩阵，,为整体节点位移列阵；,P,为整体等价节点荷载列阵。如下：,（1-50）,i,j,m,i,j,m,例：平面三角单元,双行双列,1.7.3 结构刚度矩阵特性,1、结构刚度矩阵元素的力学意义,把方程(1-50）写开，,=1,=0,=0,=0,=0,=0,(1-51）,2、结构刚度矩阵是对称矩阵,已知单元刚度矩阵是对称矩阵（1.7节），用单元刚度矩阵组集,结构刚度矩阵的,过程又没有破坏,其对称性，结构,刚度矩阵必然也,是对称的。当然，,对称性也可以通,过虚功原理得到证明。,结构刚度矩阵中的任一元素,k,ij,是,j,为单位位移（,j,=1）,，其它位移为零时的P,i,。,3、结构刚度矩阵主对角线上的元素恒为正值,由性质(1)可知，任一主对角线上元素,k,ii,是使节点位移,i,为一单位位移，其它节点位移为零时必须在第,i,号位移方向施加的力P,i,。它的方向自然应与位移方向相同，因而是正值。,4、结构刚度矩阵是一个稀疏矩阵,稀疏矩阵指：存在大量零元素。非零元素稀疏排列。,矩阵的每一列都有很多零元素。考察矩阵中第,j,列。,再分析图（1-14）。设节点b发生单位位移,j,=1，其它位移为零时，,j,只能在与点节b有直接联系的 q 、 r节点引起节点力，不能,在其它节点引起节点,力。所以式（1-52）,中，只有和q、p、r、,b节点位移的相关元素,才不为零，其余的元素,都是零元素。,任一元素,k,ij,是,j,=1（其它=0）引起的P,i,（,i,=1、2）,（1-52）,j,=1,t,图（1-14）,p,q,r,s,c,b,b,其它各列的情况也是类似的。,结构的节点总数通常都比直接环绕于任何一个节点的节点数大得多，因而，结构刚度矩阵中很大一部分元素是零，即所谓的稀疏矩阵。,5、结构刚度矩阵是一个奇异矩阵,从单元刚度矩阵的奇异性讨论中知，,处于静力平衡状态的无约束单元可以发生任意的刚体位移,。与单元刚度矩阵是奇异矩阵的理由一样，无约束结构的结构刚度矩阵K也是奇异矩阵，即K的行列式为零。,1.7.6 引入支承约束的结构节点平衡方程,6、结构刚度矩阵是常量矩阵,结构刚度矩阵是常量矩阵。结构的节点力和节点位移成线性关系都是基于弹性理论的结果。,（1-53）,用平衡方程（1-53）是解不出结构的节点位移,的，因为结构刚度矩阵是奇异矩阵。因此，必须,引入约束，排除任何刚体位移，使结构为几何不变体系。,方程（1-53）中的刚度矩阵K和节点荷载向量列阵,P可分割为约束和自由两部分：,（1-54）,式中，,P,r,是支承反力，约束位移,自由,约束,(1-55）,(1-56）,展开（154），有：,K,ff,引入约束后的结构刚度矩阵。它通对K,引入约束后获得，,具体方法: 从无约束的结构刚度矩阵K中删去与受约束位移号对应的行和列，再将矩阵压缩排列成nn阶方阵，即为约化后的结构刚度矩阵K,ff, 。,K,ff,这是一个非奇异矩阵，它存在逆矩阵。,方程(1-55）是引入约束后的结构节点平衡方程，用于计算结构所有非刚性约束节点的节点位移。而,方程(1-60）可以用来计算结构所有受刚性约束节点的反力。,(1-61）,由式(1-55）即可解出全部未知的节点位移：,1.7.7 节点位移和单元力的解答,1、结构节点位移,2、支座反力,把解出的,f,代入(1-56），即得支座反力P,r,：,关于方程 (1-61）的解算方法，当K,ff,采用本章中上述方法组集时，可直接采用结构力学中的高斯（Gauss）法求解。,(1-56）,至此，我们可以看出：系统分析的主要任务是：,（1）组集引入约束后的结构刚度矩阵K,ff,；,（2）求解式(1-55）给出的线性代数方程组。算出全部未知的节点位移。,至于支座反力的计算，实际计算时，根本不去组集式(1-56）中的矩阵K,rf,，不用式(1-62）而是直接对,支承点使用节点平衡,方程计算。,3、单元力,(1-62）,（1）全解公式,全解是指结构在未经简化的实际荷载作用下（全解系统）的解答。,余解,是结构在,节点荷载,作用下（余解系统）的解答。节点荷载包括原来给定的节点荷载和由分布荷载简化而来的等价节点荷载。,特解,是单元从实际分布荷载到等价节点力过程中（特解系统）导致的单元节点力。,全解 = 余解 + 特解,(1-63）,通过式(1-55）解出获得全部未知节点位移,，据此找出任何一个单元的单元位移（单元节点位移），然后可以求得单元力。,图(1-15）,i,j,M,i,M,j,R,i,R,j,=,+,(e),i,j,（a）全解系统,（b）特解系统,（c）余解系统,（2）余解,-M,i,-M,j,-R,i,-R,j,i,j,（3）特解,把单元等价节点力反向作用在单元的相应节点，即为单元力的特解。,4、结论,（1）全解的节点位移等于余解的节点位移,因为，特解给出的节点位移总是为零。,（2）对不承受分布荷载的单元，内力由余解直接给出。,（3）对承受分布荷载作用的单元，其内力由余解和特解的迭加给出。,F,全,e,= F,余,e, F,q,e,(1-64）,等价节点力,3.9 举例,如图1-17a所示两端固支的矩形深梁，跨度为2a，梁高为a，截面宽度为h，已知，,=0、,E，承受均布压力q。试用有限元法解此平面应力问题。,a,2a,q,图1-17,(a),利用对称性，可取梁的一半分析，共分2个单元，4个结点，支座约束如(a)图。,解,1、划分单元。,1,2,4,3,1,2,3,4,x,y,a,a,qa/2,(b),图1-18,a,a,1,2,3,4,7,8,5,6,（C)自由度,等腰直角三角形单元dde ,k,(1-43b）,i j m,i,j,m,2、单元刚度矩阵,由书公式（1-46b）,3 4 5 6 1 2,3,4,5,6,1,2,1,3,5 6 3 4 7 8,5,6,3,4,7,8,3,3,3、结构刚度矩阵,引入约束前，结构刚度矩阵是88阶的。受约束的自由度号为1，3，4，5，7，8。删除相应行和列后，结构刚度矩阵剩下2行2列，相应的自由度号是：2，6。,2 6,2,6,3/2,-1,1,-1,1/2,+,3/2,4、结构节点荷载向量,一般情况节点荷载向量有8个元素。引入约束后剩下2个，相应序号为：2，6。,2,6,5、结构无约束节点平衡方程,解得,根据结点位移,找单元位移,，代入（1-38）式可求各单元的单元力,F,。,（1-38）,利用单元力,F计算支点反力。,1,2,3,4,a,a,1,2,3,4,7,8,5,6,将单元位移,代入式（1-28），其中S可由式（1-30）、（1-31）求得，从而可求出单元应力。,（1-28）,习题：,算出式（1-28）的结果.,剖分结构为有限个单元，对节点、单元编号,组集单元刚度矩阵和单元等价节点力向量,组集结构刚度矩阵并引入约束,组集结构等价节点荷载向量,从节点平衡方程解未知节点位移,计算结构内力、应力,