《5平面分析整体刚度矩阵》由会员分享,可在线阅读,更多相关《5平面分析整体刚度矩阵(46页珍藏版)》请在金锄头文库上搜索。
1、2024/9/9平面问题有限元分析-总刚15 平面问题有限元分析平面问题有限元分析整体刚度矩阵整体刚度矩阵曹国华曹国华5.1 整体刚度矩阵整体刚度矩阵5.2 整体刚度矩阵的特点整体刚度矩阵的特点5.3 边界条件边界条件5.4 计算结果整理计算结果整理5.5 收敛准则收敛准则 前 文 对 单 元 体 进 行 了 分 折 , 得 到 了 单 元 刚 度 方 程 ,但要解决问题,还必须进一步建立整个计算模型的整体刚度方程。完成这一步的关键,在于怎样将单元的刚度矩阵和节点荷载列阵,分别“组装”成整体刚度矩阵和整体节点荷载列阵。 这里通过研究任意节点的平衡来建立整体刚度矩阵,该方法不但比较直观、易懂,而
2、且对怎样编写计算机程序是很有帮助的。2024/9/92平面问题有限元分析-总刚5.1 整体刚度矩阵整体刚度矩阵 整体节点载荷列阵:由各节点所受等效节点力按节点号码以从小到大的顺序排列组成的列阵。等效节点力是由集中力、表面力和体积力共同移置构成的,其中集中力包括直接作用在弹性体上的外力和边界约束力,如支座反力。为了研究整体刚度矩阵的组装过程,先引入两个概念。 整体节点位移列阵:由各节点位移按节点号码以从小到大的顺序排列组成的列阵。2024/9/93平面问题有限元分析-总刚5.1 整体刚度矩阵整体刚度矩阵式中: , 不失一般性,仅考虑计算模型中有4个单元,如图所示。四个单元的整体节点位移列阵为20
3、24/9/94平面问题有限元分析-总刚5.1 整体刚度矩阵整体刚度矩阵 对每个单元都可以写出相应的单元刚度方程 ,即单元节点平衡方程。例如,对号单元,有 式中:号单元中第i(i=1,2,3)节点所受力。为了便于组装整体刚度矩阵,将上式以整体节点位移 表示,即 号单元的扩大刚度矩阵或称为单元贡献矩阵。2024/9/95平面问题有限元分析-总刚5.1 整体刚度矩阵整体刚度矩阵 同理,对于单元,有式中:号单元中第i(i=1,3,4)节点所受力; 号单元的扩大刚度矩阵。2024/9/96平面问题有限元分析-总刚5.1 整体刚度矩阵整体刚度矩阵对于3单元,有式中:3号单元中第i(i=3,4,5)节点所受
4、力; 3号单元的扩大刚度矩阵。2024/9/97平面问题有限元分析-总刚5.1 整体刚度矩阵整体刚度矩阵对于单元,有式中:号单元中第i(i=3,4,5)节点所受力; 号单元的扩大刚度矩阵。2024/9/98平面问题有限元分析-总刚5.1 整体刚度矩阵整体刚度矩阵 对于任意一个节点,可能承受两种力的作用,一种是其它单元给予该节点的反作用力;另一种是作用在节点上的等效节点力。对整体而言,前者属于内力,后者属于外力,每个节点在两种力的作用下处于平衡。 将各单元刚度方程左边相加,即将各节点所受力相加,由于对于整体而言,单元给予节点的反作用力属于内力,在相加过程中相互抵消,所以各节点所受力相加的结果只有
5、外力,即等效节点力,从而得到整体节点荷载列阵,如下2024/9/99平面问题有限元分析-总刚5.1 整体刚度矩阵整体刚度矩阵 将各单元刚度方程右边相加,从而得到整体刚度矩阵,如下2024/9/910平面问题有限元分析-总刚5.1 整体刚度矩阵整体刚度矩阵 通过以上分析得,整体节点载荷与整体节点位移之间的关系式,即结构整体有限元方程,如下式中:整体刚度矩阵。2024/9/911平面问题有限元分析-总刚5.1 整体刚度矩阵整体刚度矩阵 整体刚度矩阵组装的基本步骤:1)将单元刚度矩阵中的每个子块放到在整体刚度矩阵中的对应位置上,得到单元的扩大刚度矩阵。注意对于单元刚度矩阵是按照局部编码排列的,即对应
6、单元刚度矩阵中的i、j、m;对于整体刚度矩阵是按照整体编码排列的,即按节点号码以从小到大的顺序排列。在组装过程中,必须知道单元节点的局部编码与该节点在整体结构中的整体编码之间的关系,才能得到单元刚度矩阵中的每个子块在整体刚度矩阵中的位置。将单元刚度矩阵中的每个子块按总体编码顺序重新排列后,可以得到单元的扩大矩阵。例如在图中,单元的局部编码为i、j、m,对应整体编码为1、3、4,然后将单元刚度矩阵中的每个子块按总体编码顺序重新排列后,可以得到单元的扩大矩阵。注意有些书籍中将局部编码表示为1、2、3或1,2,3等;2)将全部单元的扩大矩阵相加得到整体刚度矩阵。2024/9/912平面问题有限元分析
7、-总刚5.1 整体刚度矩阵整体刚度矩阵通过以上组装过程可以得到组装整体刚度矩阵的一般规则: 1)结构中的等效节点力是相关单元结点力的叠加,整体刚度矩阵的子矩阵是相关单元的单元刚度矩阵子矩阵的集成; 2)当整体刚度矩阵中的子矩阵 中r=s时,该节点(节点r或s)被哪几个单元所共有,则 就是这几个单元的刚度矩阵中的子矩阵 的相加。如 应该是单元中对应子矩阵的集成,即 2024/9/913平面问题有限元分析-总刚5.1 整体刚度矩阵整体刚度矩阵 3)当 中 时,若rs边是组合体的内边,则 就是共用该边的两相邻单元刚度矩阵中的子矩阵 的相加。如13边为单元和的共用边,则 4)当 中r和s不同属于任何单
8、元时,则 =0。如节点r1和 s5不同属于任何单元,此时 =0。2024/9/914平面问题有限元分析-总刚5.1 整体刚度矩阵整体刚度矩阵 上述组装基本步骤和规则具有普遍意义 对于不同类型、不同形式的单元,只是相应节点的子矩阵的阶数(节点自由度节点自由度)可能不同,至于组装整体刚度矩阵的规律仍是相同的。正是因为有了这种组装规律,使得有限元法能够很方便地应用电子计算机进行计算。2024/9/915平面问题有限元分析-总刚5.1 整体刚度矩阵整体刚度矩阵 例:例: 如图所示有限元模型,弹性模量为 ,厚度为 ,为简化计算取 ,求整体刚度矩阵。2024/9/916平面问题有限元分析-总刚5.1 整体
9、刚度矩阵整体刚度矩阵单元元编号号整体整体编码1、2、32、4、55、3、23、5、6局部局部编码i、j、mi、j、mi、j、mi、j、m以整体以整体编码表表示的示的单元元刚度度矩矩阵子子块解:解:该模型中共有6个节点,4个单元,各单元的信息如表所示。 各单元信息2024/9/917平面问题有限元分析-总刚5.1 整体刚度矩阵整体刚度矩阵同上例类似的分析,得 根据单元刚度矩阵的性质,可知 ,若3单元5,3,2,则 整体刚度矩阵中的各子块是对所有单元相应的子块求和得到的(实际只是对相关单元求和),其中各子块矩阵均为2行2列,整体刚度矩阵用子块矩阵可以表示为2024/9/918平面问题有限元分析-总
10、刚5.1 整体刚度矩阵整体刚度矩阵上式中任意一子块矩阵均为2行2列,在计算过程中,无需将每个单元刚度矩阵进行扩大,只需判断整体刚度矩阵子块的下标,然后利用组装整体刚度矩阵的一般规则进行计算,如 ,由图形可知,节点2由单元、和所共有,则 2024/9/919平面问题有限元分析-总刚5.1 整体刚度矩阵整体刚度矩阵,由图形可知,25边为单元和的共用边,则,由图形可知,节点1、5不同属于任何单元,则采用同样的方法进行计算,得到整体刚度矩阵为2024/9/920平面问题有限元分析-总刚5.1 整体刚度矩阵整体刚度矩阵2024/9/921平面问题有限元分析-总刚5.1 整体刚度矩阵整体刚度矩阵MATHC
11、AD例子例子公共邮箱:公共邮箱:密码:密码: cncumtcncumt 作业:作业: 采用采用MATHCAD求解求解如图所示有限元模型,弹性模量为 ,厚度为 ,为简化计算取 ,E=1,t=1 ,求整体刚度矩阵。2024/9/924平面问题有限元分析-总刚5.1 整体刚度矩阵整体刚度矩阵每个班由负责人压缩成一每个班由负责人压缩成一个压缩包后发到个压缩包后发到每个人提交的文件名为:每个人提交的文件名为:班级班级-姓名姓名-序号(序号(1)比如:比如:机自机自10(1)班)班-某某某某-12(1)4. 4. 是奇异矩阵,在排除刚体位移后,它是正定阵是奇异矩阵,在排除刚体位移后,它是正定阵 3. 3.
12、 是稀疏矩阵,非零元素呈带状分布是稀疏矩阵,非零元素呈带状分布 用有限元方法分析复杂工程问题时,节点的数目比较多,整体刚度矩阵的阶数通常也是很高的。那么在进行计算时,如果存储整体刚度矩阵的全部元素,将会浪费较大的资源、降低计算效率。如果根据整体刚度矩阵的特点进行编写程序,可以大大节省资源、并提高计算效率。因此有必要了解和掌握整体刚度矩阵的特点,整体刚度矩阵具有以下几个显著的特点:1. 1. 是对称矩阵是对称矩阵2. 2. 中主对角元素总是正的中主对角元素总是正的2024/9/925平面问题有限元分析-总刚5.2 整体刚度矩阵的特点整体刚度矩阵的特点1. 是对称矩阵 由单元刚度矩阵的对称性和整体
13、刚度矩阵的组装过程,可知整体刚度矩阵必为对称矩阵。利用对称性,在计算机编写程序时,只存储整体刚度矩阵上三角或下三角部分即可。2. 中主对角元素总是正的 例如,刚度矩阵 中的元素 表示节点2在x方向产生单位位移,而其它位移均为零时,在节点2的x方向上必须施加的力; 表示节点2在y方向产生单位位移,而其它位移均为零时,在节点2的y方向上必须施加的力。很显然在此情况下力的方向应该与位移方向一致,故应为正号。2024/9/926平面问题有限元分析-总刚5.2 整体刚度矩阵的特点整体刚度矩阵的特点3. 是稀疏矩阵,非零元素呈带状分布 如果遵守一定的节点编号规则,就可使矩阵的非零元素都集中在主对角线附近呈
14、带状。整体刚度矩阵中的子矩阵 只有当下标s等于r或者s与r同属于一个单元时才不为零,这就说明,在第r双行中非零子矩阵的块数,应该等于节点r周围直接相邻的节点数目加1。可见, 中元素一般都不是填满的,是稀疏矩阵,且非零元素呈带状分布。 以下图所示的单元网格为例,其整体刚度矩阵中的非零子块(每个子块为2行2列)的分布情况如下图所示,图中阴影部分表示该子块不为零,其它子块部位均为零。2024/9/927平面问题有限元分析-总刚5.2 整体刚度矩阵的特点整体刚度矩阵的特点2024/9/928平面问题有限元分析-总刚5.2 整体刚度矩阵的特点整体刚度矩阵的特点2024/9/929平面问题有限元分析-总刚
15、5.1 整体刚度矩阵整体刚度矩阵 显然,带状刚度矩阵的带宽取决于单元网格中相邻节点号码的最大差值D。把半个斜带形区域中各行所具有的非零元素的最大个数叫做刚度矩阵的半带宽(包括主对角元),用B表示,如下 B=2(D+1) 通常的有限元程序,一般都利用刚度矩阵的对称和稀疏带状的特点,在计算求解中,只存储上半带的元素,即所谓的半带存储。因此,在划分完有限元网格进行节点编号时,要采用合理的编码方式,使同一单元中相邻两节点的号码差尽可能小,以便节省存储空间、提高计算效率。2024/9/930平面问题有限元分析-总刚5.2 整体刚度矩阵的特点整体刚度矩阵的特点 对于同样的有限元单元网格,按照图(a)的结点
16、编码,最大的半带宽为 按照图(b)的结点编码,最大的半带宽为B 按照图(c)的结点编码,最大的半带宽为B (a) (b) (c)2024/9/931平面问题有限元分析-总刚2(106)105.2 整体刚度矩阵的特点整体刚度矩阵的特点B = 2(104)=142(102)184. 是奇异矩阵,在排除刚体位移后,它是正定阵 无约束的弹性体(或结构物)的整体刚度矩阵是奇异的,不存在逆矩阵,即关于位移的解不唯一。这是因为弹性体在外力的作用下处于平衡,外力的分量应该满足三个静力平衡方程。这反映在整体刚度矩阵中就意味着存在三个线性相关的行或列,所以是个奇异阵,不存在逆矩阵。例如:设弹性体在外力的作用下处于
17、平衡,这时相应的解为 ,然后在给予弹性体以刚体位移而相应的节点位移 ,这时, 仍是问题的解,因为刚体位移不会破坏平衡。注:当排除刚体位移后,整体刚度矩阵是正定矩阵。2024/9/932平面问题有限元分析-总刚5.2 整体刚度矩阵的特点整体刚度矩阵的特点 前文已经提到在排除刚体位移后,整体刚度矩阵是正定的,方程才可求得唯一解。排除刚体位移可以通过引入边界约束条件来实现,这里介绍两种比较简单的引入已知位移的方法。1代入法2乘大数法2024/9/933平面问题有限元分析-总刚5.3 边界条件边界条件1代入法 该方法保持方程组仍为2n2n系统,仅对整体刚度矩阵 和整体载荷列阵进行修正。下面以一个只有四
18、个方程的简单例子加以说明,方程如下 假定系统中节点位移 、 ,则当引入这些节点的已知位移之后,方程就变成2024/9/934平面问题有限元分析-总刚5.3 边界条件边界条件若 ,则 然后,用这组维数不变的方程来求解所有的节点位移。显然,其解答仍为原方程的解答。 在手算时,可直接将零位移约束所对应的整体刚度矩阵中的行和列直接划去,使得整体刚阵的维数变小,更便于手算。2024/9/935平面问题有限元分析-总刚5.3 边界条件边界条件2乘大数法 将 中与指定的节点位移相对应的主对角元素乘上一个大数,同时将 中的对应元素换成指定的节点位移值、该大数与节点位移相对应的主对角元素三者的乘积。若把此方法用
19、于上面的例子,则方程就变成该方程组的第一个方程为解得 ,这种方法就是使 中相应行的修正项远大于非修正项。2024/9/936平面问题有限元分析-总刚5.3 边界条件边界条件1代入法代入法2乘大数法乘大数法 在以上的两种方法中,代入法接近人工解法,虽然该方法比较直观,但该方法对刚度矩阵改变较多,程序效率不高。乘大数法对刚度矩阵改变较少,工作量较小,但相乘的“大大数数”若取得过大,求解时会发生“溢溢出出”、若取得太小则会引起较大的误差误差。2024/9/937平面问题有限元分析-总刚5.3 边界条件边界条件 对于三节点三角形单元,单元内各点的应力值相等,算出的应力一般作为单元形心处的应力。由于单元
20、应力为常数,整个结构的应力场呈阶梯状,在单元之间不连续。而工程上往往更加关心边界和节点上的受力情况,因此,必须对所得到的应力再次进行处理,得到更加合理的应力场,并得到所需点上的应力值。 这里介绍两种简单的方法,一种方法称为节点平均法,即把环绕某一节点的各单元的应力加以平均作为该节点的应力值。例如图中节点3的应力为2024/9/938平面问题有限元分析-总刚5.4 计算结果整理计算结果整理 为了使通过这样平均得来的应力比较接近实际情况,要求环绕节点的各单元尺寸不应相差太大。这种做法,对内节点比较好,对边界点则可能很差。 因此,边界节点处的应力不宜直接由单元应力平均来获得,而要根据内节点的应力构造
21、插值函数推算出来。例如图中边界点1的应力,可以先用节点平均法求得节点2、3、4处的应力,在构造相应的插值函数推算边界点1的应力,如常用的抛物线插值公式如下5.4 计算结果整理计算结果整理2024/9/939平面问题有限元分析-总刚 另一种方法称为单元平均法,即把两相邻单元的应力加以平均,用以表示公共边界中点的应力。为了由这样平均所得到的应力具有较好的精度,两相邻单元的面积不应相差太大。如图中单元和边界的中点处的应力为 在不同的有限元软件中均具有各自的后处理方法,但无论后期怎样处理,应力场来源于应力的计算,应力的精度主要依赖于单元的尺寸、单元的类型(位移模式)。2024/9/940平面问题有限元
22、分析-总刚5.4 计算结果整理计算结果整理 对于有限元这种数值计算方法,一般总是希望随着网格的逐步细分所得到的解能够收敛于问题的精确解。根据前面的分析,可知在有限元分析中,一旦确定了单元的形状之后,位移模式的选择将是非常关键的。由于载荷的移置、应力矩阵和刚度矩阵的建立等等,都依赖于单元的位移模式,所以,如果所选择的位移模式与真实的位移分布有很大的差别,那么就很难获得良好的数值解。为了保证解答的收敛性,要求位移模式必须满足以下三个条件,即 1)位移模式必须包含单元的刚体位移。 2)位移模式必须包含单元的常应变。 3)位移模式在单元内要连续、且在相邻单元之间的位移必须协调。 5.5 收敛准则收敛准
23、则2024/9/941平面问题有限元分析-总刚 1)位移模式必须包含单元的刚体位移。也就是说,当节点位移是由某个刚体位移所引起时,弹性体内将不会产生应变。所以,位移模式不但要具有描述单元本身形变的能力,而且还要具有描述由于其它单元形变而通过节点位移引起单元刚体位移的能力。 2)位移模式必须包含单元的常应变。每个单元的应变一般都是包含着两个部分:一部分是与该单元中各点的坐标位置有关的应变(即所谓各点的变应变);另一部分是与位置坐标无关的应变(即所谓的常应变)。从物理意义上看,当单元尺寸无限缩小时,每个单元中的应变应该趋于常量。因此,在位移模式中必须包含有这些常应变,否则就不可能使数值解收敛于正确
24、解。2024/9/942平面问题有限元分析-总刚5.5 收敛准则收敛准则 3)位移模式在单元内要连续、且在相邻单元之间的位移必须协调。当选择多项式来构成位移模式时,单元内的连续性要求总是得到满足的,单元间的位移协调性,就是要求单元之间既不会出现开裂也不会出现重叠的现象。通常,当单元交界面上的位移取决于该交界面上节点的位移时,就可以保证位移的协调性。 在有限单元法中,把能够满足条件1和2的单元,称为完备单元,完备单元是收敛的必要条件。满足条件3的单元,叫做协调单元或保续单元。前面讨论过的三节点三角形单元,均能同时满足上述三个条件,因此属于完备的协调单元,完备的协调单元是收敛的充分条件。在某些梁、
25、板及壳体分析中,要使单元满足条件3比较困难,所以实践中有时也出现一些只满足条件1和2的单元,其收敛性往往也能够令人满意。目前,完备而不协调的单元,已获得了很多成功应用。2024/9/943平面问题有限元分析-总刚5.5 收敛准则收敛准则 在有限元分析中,将实际连续体分成许多单元体的组合后,根据线性或非线性位移的假定,人为地选择一个位移场,通过这些措施所得到的模型比实际连续体的刚性要高。 因而,近似模型的刚度是实际连续体刚度的上界。若选择不协调单元,那么这种模型可能由于单元分离、叠加或单元之间形成铰而降低刚性,这种影响就有可能使得不协调元比应用协调元所得的结果要好。不过,应用不协调单元事先不能肯
26、定所得的刚度是真实刚度的上界。换句话说,不协调元不一定象协调单元那样刚硬,可能比较柔软,因此有可能会比协调单元收敛得快。2024/9/944平面问题有限元分析-总刚5.5 收敛准则收敛准则 经验证明,根据巴斯卡(Pascal)三角形来选择二维多项式的各项。在二维多项式中,如果包含有对称轴一边的某一项,那么就必须同时包含有另一边的对称项。选择多项式位移模式时,还应该要考虑的一个因素是,多项式中的项数必须等于或稍大于单元边界上的外节点的自由度数。通常是取项数与单元的外节点的自由度数相等,取过多的项数是不恰当的。2024/9/945平面问题有限元分析-总刚收敛准则收敛准则2024/9/946平面问题有限元分析-总刚5.5 收敛准则收敛准则