第四章平面问题的有限单元法课件.ppt

《第四章平面问题的有限单元法课件.ppt》由会员分享，可在线阅读，更多相关《第四章平面问题的有限单元法课件.ppt（55页珍藏版）》请在淘文阁 - 分享文档赚钱的网站上搜索。

1、第四章第四章 平面问题的有限单元法平面问题的有限单元法弹性力学的平面应力问题基本条件基本条件（1 1）等厚度的）等厚度的薄板薄板；（2 2）体力体力作用于体内，平行于板的中面，沿板厚不变；作用于体内，平行于板的中面，沿板厚不变；（3 3）面力面力作用于板边，平行于板的中面，沿板厚不变；作用于板边，平行于板的中面，沿板厚不变；（4 4）约束约束作用于板边，平行于板的中面，沿板厚不变。作用于板边，平行于板的中面，沿板厚不变。弹性力学的平面应力问题坐标系：坐标系：由于两板面上无面力和由于两板面上无面力和约束作用：约束作用：由于薄板很薄，应力是连续变化的，又无由于薄板很薄，应力是连续变化的，又无z z

2、向外力，可认为：向外力，可认为：简化为平面应力问题，仅剩：简化为平面应力问题，仅剩：其值与其值与z z无关无关弹性力学的平面应变问题基本条件基本条件（1 1）很长的）很长的常截面柱常截面柱；（2 2）体力体力作用于体内，平行于横截面，沿柱体长度方向不变；作用于体内，平行于横截面，沿柱体长度方向不变；（3 3）面力面力作用于柱面，平行于横截面，沿柱体长度方向不变；作用于柱面，平行于横截面，沿柱体长度方向不变；（4 4）约束约束作用于柱面，平行于横截面，沿柱体长度方向不变。作用于柱面，平行于横截面，沿柱体长度方向不变。弹性力学的平面应变问题坐标系：坐标系：由于由于截面、外力、约截面、外力、约束沿束

3、沿z z 向不变，外力、向不变，外力、约束平行约束平行xyxy面，柱体面，柱体非常长非常长：故任何故任何z z 面面（截面）均为对称面。（截面）均为对称面。简化为平面应变问题：简化为平面应变问题：其值与其值与z z无关无关以平面应力问题为例介绍以平面应力问题为例介绍平面问题的有限单元法平面问题的有限单元法平面应力单元类型平面应力单元类型简介平面应力单元类型简介3 3节点三角形单元节点三角形单元4 4节点节点4 4边形单元边形单元8 8节点节点4 4边形曲边单元边形曲边单元节点位移分量节点位移分量每节点每节点2 2个位移分量（自由度）个位移分量（自由度）x x方向的位移方向的位移u u，y y方

4、向的位移方向的位移v v单元位移分量（单元位移分量（4 4节点）节点）jik三角形单元三角形单元单元单元ekijl单元单元e四边形单元四边形单元123456788节点单元节点单元单元单元e平面应力单元网格划分应力梯度变化比较大的地方，网格应密一些应力梯度变化比较大的地方，网格应密一些有应力集中的地方，网格应密一些有应力集中的地方，网格应密一些单元边界长度不要相差过大单元边界长度不要相差过大单元各边夹角不要太大单元各边夹角不要太大集中载荷处要设置节点集中载荷处要设置节点结构不同材料交界面处要设置节点并作为单元边界结构不同材料交界面处要设置节点并作为单元边界结构厚度突变处要设置节点并作为单元边界结

5、构厚度突变处要设置节点并作为单元边界分布载荷突变处要设置节点分布载荷突变处要设置节点施加位移约束处要设置节点施加位移约束处要设置节点注意单元间的连接注意单元间的连接平面应力单元网格划分设置节点设置节点设置节点设置节点材料材料A材料材料B界面界面这样不行这样不行病态单元病态单元a-边长差别太大边长差别太大b-边长差别太大边长差别太大c-边夹角太大边夹角太大 abc单元节点信息节点信息节点信息节点号xyz100021003200421051106010700.50810.50920.50单元拓扑信息单元拓扑信息单元号节点i节点j节点k节点l材料编号其它常数112871278561358941439

6、881582331以最经典的三角形单元为例以最经典的三角形单元为例单元位移函数（位移模式）单元位移模式概念单元位移模式概念单元内任一点的位移要用节点上的位移值近似表达出来，这单元内任一点的位移要用节点上的位移值近似表达出来，这就需要假定一个近似函数来表示单元内的位移分布，所选择就需要假定一个近似函数来表示单元内的位移分布，所选择的近似函数就称为单元位移函数或单元位移模式。的近似函数就称为单元位移函数或单元位移模式。对于弹性力学平面问题，一般选择多项式对于弹性力学平面问题，一般选择多项式 (polynomial)(polynomial)来作为单元内的位移解或插值函数或位移模式。来作为单元内的位移

7、解或插值函数或位移模式。一维单元一维单元一维单元一维单元二维单元二维单元二维单元二维单元多项式的项数越多，结果就越精确，但取多少项由单元形式决定。多项式的项数越多，结果就越精确，但取多少项由单元形式决定。多项式的项数越多，结果就越精确，但取多少项由单元形式决定。多项式的项数越多，结果就越精确，但取多少项由单元形式决定。三角单元的位移函数节点上只有六个位移分量，所以节点上只有六个位移分量，所以单元内部位移函数的待定参数不单元内部位移函数的待定参数不能超过这个数目。可假设单元内能超过这个数目。可假设单元内部位移为部位移为x x、y y的线性函数：的线性函数：参数参数ai由位移边界条件确定。由位移边

8、界条件确定。三角单元的位移函数节点节点i i节点节点j j节点节点k k于是：于是：三角单元的位移函数如果令如果令则：则：根据线性代数的知识，可知：根据线性代数的知识，可知：三角单元的位移函数T T*为为T T的伴随矩阵的伴随矩阵其中：其中：三角单元的位移函数把求得的系数把求得的系数代入位移函数公式：代入位移函数公式：得到：得到：三角单元的位移函数表达为矩阵形式：表达为矩阵形式：这里：这里：N Ni i,N,Nj j,N,Nk k是坐标的函数，它们反映了单元的位移形态，故是坐标的函数，它们反映了单元的位移形态，故称为三角单元的称为三角单元的形态函数形态函数（或（或形函数形函数）三角单元的位移函

9、数形函数具有明确的几何意义：形函数具有明确的几何意义：如图所示三角单元如图所示三角单元IJKIJK，P P为三角单为三角单元内任意一点，其坐标为元内任意一点，其坐标为(x x,y y)P P点在三角单元各角点上产生的形点在三角单元各角点上产生的形函数分别是函数分别是N Ni i,N,Nj j,N,Nk k同理：同理：三角单元的位移函数位移函数运用示例：位移函数运用示例：已知各节点位移为：已知各节点位移为：求求P P点位移点位移P P P P点的位移可由节点位移近似表示为点的位移可由节点位移近似表示为点的位移可由节点位移近似表示为点的位移可由节点位移近似表示为三角单元的位移函数于是：于是：三角单

10、元的位移函数形函数的本质形函数的本质计算点计算点计算点计算点(x x x x,y y y y)的位移的位移的位移的位移u u u u(x x x x,y y y y)、v v v v(x x x x,y y y y)可用单元内各节点的可用单元内各节点的可用单元内各节点的可用单元内各节点的位移值位移值位移值位移值u u u ui i i i,v v v vi i i i的加权之和来近似表示，其中，各节点位移的加权之和来近似表示，其中，各节点位移的加权之和来近似表示，其中，各节点位移的加权之和来近似表示，其中，各节点位移加权系数为关于计算点加权系数为关于计算点加权系数为关于计算点加权系数为关于计算

11、点(x,yx,yx,yx,y)的函数，即为形函数的函数，即为形函数的函数，即为形函数的函数，即为形函数三角形单元形函数的性质三角形单元形函数的性质1 1 1 1、单元节点产生的形函数值为、单元节点产生的形函数值为、单元节点产生的形函数值为、单元节点产生的形函数值为1 1 1 1或或或或0 0 0 02 2 2 2、形函数之和等于、形函数之和等于、形函数之和等于、形函数之和等于1 1 1 1位移函数应满足的条件应满足：应满足：单元内位移模式必须是连续的，公共边上位移必单元内位移模式必须是连续的，公共边上位移必须协调须协调位移模式必须反映单元的刚体位移位移模式必须反映单元的刚体位移位移模式必须反映

12、单元的常应变位移模式必须反映单元的常应变可以证明三节点三角形单元是收敛的可以证明三节点三角形单元是收敛的完备单元和协调单元三条准则：三条准则：1 1、位移模式必须包含单元的刚体位移、位移模式必须包含单元的刚体位移2 2、位移模式必须能包含单元的常应变、位移模式必须能包含单元的常应变3 3、位移模式在单元内要连续、并使相邻单元间的位移必、位移模式在单元内要连续、并使相邻单元间的位移必须协调须协调满足条件满足条件1 1、2 2的的单元为单元为完备单元完备单元满足条件满足条件3 3的的单元为单元为协调单元协调单元应变的离散过程应变的离散过程应变的离散过程根据弹性力学中的几何关系，单元内任一点根据弹性

13、力学中的几何关系，单元内任一点(x,yx,y)的应变表达式的应变表达式为为矩阵形式矩阵形式应变的离散过程应变的离散过程应变的离散过程单元内任一点单元内任一点(x,yx,y)的位移的位移(u,vu,v)可以采用节点位移近似表示：可以采用节点位移近似表示：将其代入应变表达式，则将其代入应变表达式，则应变的离散过程应变的离散过程应变的离散过程 为书写方便，应变分量矩阵可用分块矩阵表示为书写方便，应变分量矩阵可用分块矩阵表示为书写方便，应变分量矩阵可用分块矩阵表示为书写方便，应变分量矩阵可用分块矩阵表示简写为：简写为：B B也称为也称为“应变矩阵应变矩阵”应变的离散过程应变的离散过程应变的离散过程 由

14、于形函数由于形函数由于形函数由于形函数 所以刚才的矩阵事实上可表示为所以刚才的矩阵事实上可表示为所以刚才的矩阵事实上可表示为所以刚才的矩阵事实上可表示为B B矩阵中的所有元素已经由三角形单元的节点坐标确定。矩阵中的所有元素已经由三角形单元的节点坐标确定。应变在单元内为常数，所以又称为常应变单元。应变在单元内为常数，所以又称为常应变单元。应力的离散过程应力的离散过程应力的离散过程 根据广义虎克定律，对于平面应力问题：根据广义虎克定律，对于平面应力问题：根据广义虎克定律，对于平面应力问题：根据广义虎克定律，对于平面应力问题：矩阵形式矩阵形式矩阵形式矩阵形式物理方程物理方程物理方程物理方程应力的离散

15、过程应力的离散过程应力的离散过程 如果令：如果令：如果令：如果令：物理方程简写为物理方程简写为物理方程简写为物理方程简写为 D D又称为又称为又称为又称为“弹性矩阵弹性矩阵弹性矩阵弹性矩阵”将前面应变的表达式代入，则将前面应变的表达式代入，则将前面应变的表达式代入，则将前面应变的表达式代入，则虚位移与虚应变我们已经知道了应变与位移的关系我们已经知道了应变与位移的关系 那么很自然的那么很自然的那么很自然的那么很自然的 如果发生了虚位移如果发生了虚位移如果发生了虚位移如果发生了虚位移则会发生虚应变则会发生虚应变虚功原理建立控制方程外力虚功等于内力虚功。外力虚功等于内力虚功。外力虚功外力虚功外力虚功

16、外力虚功 内力虚功内力虚功内力虚功内力虚功 其具体可计算为：其具体可计算为：其具体可计算为：其具体可计算为：虚功原理建立控制方程外力虚功等于内力虚功。外力虚功等于内力虚功。结果：结果：结果：结果：考虑到节点虚位移的任意性：考虑到节点虚位移的任意性：考虑到节点虚位移的任意性：考虑到节点虚位移的任意性：上式即为上式即为上式即为上式即为有限元控制方程有限元控制方程有限元控制方程有限元控制方程。此处此处此处此处KK称为称为称为称为“刚度矩阵刚度矩阵刚度矩阵刚度矩阵”刚度矩阵如果将求解域划分为多个单元，则如果将求解域划分为多个单元，则 即即即即总体刚度矩总体刚度矩总体刚度矩总体刚度矩阵（总刚）阵（总刚）

17、阵（总刚）阵（总刚）单元刚度单元刚度单元刚度单元刚度矩阵（单矩阵（单矩阵（单矩阵（单刚）刚）刚）刚）单元刚度矩阵三节点等厚三角形单元中三节点等厚三角形单元中B B和和D D的分量均为常量，的分量均为常量，则单元刚度矩阵可以表示为则单元刚度矩阵可以表示为其具体形式为其具体形式为其具体形式为其具体形式为单元刚度矩阵对于平面应力问题，其具体可计算如下：对于平面应力问题，其具体可计算如下：对于平面应力问题，其具体可计算如下：对于平面应力问题，其具体可计算如下：单元刚度矩阵的物理意义把前面获得的有限元控制方程展开：把前面获得的有限元控制方程展开：把前面获得的有限元控制方程展开：把前面获得的有限元控制方程

18、展开：那么那么那么那么事实上就是当节点事实上就是当节点事实上就是当节点事实上就是当节点j j j j产生单位位移时，在节点产生单位位移时，在节点产生单位位移时，在节点产生单位位移时，在节点i i i i上需要施加的上需要施加的上需要施加的上需要施加的节点力。节点力。节点力。节点力。单元刚度矩阵的物理意义更具体一点：更具体一点：更具体一点：更具体一点：当节点当节点当节点当节点i i i i在垂直方向产生单位位移时，在节点在垂直方向产生单位位移时，在节点在垂直方向产生单位位移时，在节点在垂直方向产生单位位移时，在节点i i i i上需要施加的垂直节点力上需要施加的垂直节点力上需要施加的垂直节点力上

19、需要施加的垂直节点力 当节点当节点当节点当节点j j j j在水平方向产生单位位移时，在节点在水平方向产生单位位移时，在节点在水平方向产生单位位移时，在节点在水平方向产生单位位移时，在节点j j j j上需要施加的水平节点力上需要施加的水平节点力上需要施加的水平节点力上需要施加的水平节点力 当节点当节点当节点当节点j j j j在垂直方向产生单位位移时，在节点在垂直方向产生单位位移时，在节点在垂直方向产生单位位移时，在节点在垂直方向产生单位位移时，在节点i i i i上需要施加的水平节点力上需要施加的水平节点力上需要施加的水平节点力上需要施加的水平节点力 当节点当节点当节点当节点i i i i

20、在水平方向产生单位位移时，在节点在水平方向产生单位位移时，在节点在水平方向产生单位位移时，在节点在水平方向产生单位位移时，在节点j j j j上需要施加的垂直节点力上需要施加的垂直节点力上需要施加的垂直节点力上需要施加的垂直节点力单元刚度矩阵的性质性质性质性质性质1 1 1 1：对称性：对称性：对称性：对称性单元刚度矩阵的性质性质性质性质性质2 2 2 2：对角线上元素恒为正：对角线上元素恒为正：对角线上元素恒为正：对角线上元素恒为正单元刚度矩阵的性质性质性质性质性质3 3 3 3：此矩阵为奇异矩阵此矩阵为奇异矩阵意义：没有对节点施加位移约束，所以单元产生任何的刚性位移意义：没有对节点施加位移

21、约束，所以单元产生任何的刚性位移都是可以的，由力得不到位移的唯一解。都是可以的，由力得不到位移的唯一解。性质性质4 4：此矩阵的各行元素之和为零，由于对称性，各列元素之和：此矩阵的各行元素之和为零，由于对称性，各列元素之和也为零。也为零。整体刚度矩阵的形成单元刚度矩阵形成后，要将单元组成一个整体结构，即整体分析，单元刚度矩阵形成后，要将单元组成一个整体结构，即整体分析，基本方法是基本方法是刚度集成法刚度集成法，即整体刚度矩阵是单元刚度矩阵的集成。，即整体刚度矩阵是单元刚度矩阵的集成。整体刚度矩阵的集成是按整体刚度矩阵的集成是按对号入座对号入座的方式叠加的。的方式叠加的。用下面的三角形薄板作为示

22、例：用下面的三角形薄板作为示例：共计共计4 4个单元，单元节点编号为：个单元，单元节点编号为：整体刚度矩阵的形成各个单元的刚度矩阵为：各个单元的刚度矩阵为：整体刚度矩阵的形成设单元节点总数为设单元节点总数为N N，每个节点的自由度数为，每个节点的自由度数为NDOFNDOF。(对于一维情况，对于一维情况，对于一维情况，对于一维情况，NDOF=1NDOF=1NDOF=1NDOF=1；对于二维情况，；对于二维情况，；对于二维情况，；对于二维情况，NDOF=2NDOF=2NDOF=2NDOF=2；三维，；三维，；三维，；三维，NDOF=3NDOF=3NDOF=3NDOF=3)。如果是二维问题，则总自由

23、度数为如果是二维问题，则总自由度数为2N2N个，个，相应的整体刚度矩阵大小为相应的整体刚度矩阵大小为2N2N2N2N阶方阵。阶方阵。整体刚度矩阵的意义与性质Kijij 表示表示j自由度发生单位位移，其他位移为零时，自由度发生单位位移，其他位移为零时，第第i个自由度上必须施加的节点力。个自由度上必须施加的节点力。总体刚度矩阵中的元素具有如下性质：总体刚度矩阵中的元素具有如下性质：总体刚度矩阵中的元素具有如下性质：总体刚度矩阵中的元素具有如下性质：（1 1 1 1）主对角元素主对角元素主对角元素主对角元素K Kiiii 00（2 2）总体刚度矩阵总体刚度矩阵总体刚度矩阵总体刚度矩阵K K K K是

24、对称的奇异矩阵是对称的奇异矩阵是对称的奇异矩阵是对称的奇异矩阵（3 3 3 3）总体刚度矩阵总体刚度矩阵K K是带状稀疏矩阵。是带状稀疏矩阵。整体刚度矩阵的存储由于整体刚度矩阵具有对称性、稀疏性和非零元素带状分布的由于整体刚度矩阵具有对称性、稀疏性和非零元素带状分布的特点，所以没有必要将全部的整体刚度矩阵进行存储。特点，所以没有必要将全部的整体刚度矩阵进行存储。（1 1 1 1）利用对称性：利用对称性：利用对称性：利用对称性：只保存整体刚度矩阵上只保存整体刚度矩阵上只保存整体刚度矩阵上只保存整体刚度矩阵上(下下下下)三角的带宽内元素即可；三角的带宽内元素即可；三角的带宽内元素即可；三角的带宽内

25、元素即可；（2 2）利用稀疏性：利用稀疏性：利用稀疏性：利用稀疏性：在用分块表示的整体刚度矩阵中，与相关节点对应的分块才在用分块表示的整体刚度矩阵中，与相关节点对应的分块才在用分块表示的整体刚度矩阵中，与相关节点对应的分块才在用分块表示的整体刚度矩阵中，与相关节点对应的分块才能具有非零的元素，其他位置上的分块矩阵的元素为零能具有非零的元素，其他位置上的分块矩阵的元素为零能具有非零的元素，其他位置上的分块矩阵的元素为零能具有非零的元素，其他位置上的分块矩阵的元素为零（3 3 3 3）利用带状分布）利用带状分布）利用带状分布）利用带状分布：整体刚度矩阵的非零元素分布在以对角线为中心的带状区域整体刚

26、度矩阵的非零元素分布在以对角线为中心的带状区域整体刚度矩阵的非零元素分布在以对角线为中心的带状区域整体刚度矩阵的非零元素分布在以对角线为中心的带状区域中，每行具有元素的元素的个数叫做中，每行具有元素的元素的个数叫做中，每行具有元素的元素的个数叫做中，每行具有元素的元素的个数叫做“半带宽半带宽半带宽半带宽”，用，用，用，用D D D D表示。表示。表示。表示。整体刚度矩阵的存储各行的半带宽各行的半带宽D D怎么计算：怎么计算：整体刚度矩阵的存储可用一维数组可用一维数组A A来存储半带宽内的元素，而不必储存所有元素。来存储半带宽内的元素，而不必储存所有元素。本例中：总带宽本例中：总带宽则可以采用如

27、下方式存储：则可以采用如下方式存储：整体刚度矩阵的存储最大半带宽是多少？最大半带宽是多少？相邻节点的编码最大差值相邻节点的编码最大差值相邻节点的编码最大差值相邻节点的编码最大差值+1+1+1+1）NDOFNDOFNDOFNDOFDmax=(10-6+1)2=10=(10-6+1)2=10设整体刚度矩阵设整体刚度矩阵设整体刚度矩阵设整体刚度矩阵K K K K是一个是一个是一个是一个nnnnnnnn的矩阵的矩阵的矩阵的矩阵,其最大半带宽为其最大半带宽为其最大半带宽为其最大半带宽为D D D D，那么，那么，那么，那么利用带状矩利用带状矩阵的特点和对称性，只需要存储以阵的特点和对称性，只需要存储以D

28、 D为固定宽度的元素，为固定宽度的元素，这种存储方法称这种存储方法称这种存储方法称这种存储方法称为为为为二维等带宽存储二维等带宽存储二维等带宽存储二维等带宽存储。整体刚度矩阵的存储元素元素元素元素K Krsrs在整刚矩阵在整刚矩阵在整刚矩阵在整刚矩阵K K K K中的行列编码记为中的行列编码记为中的行列编码记为中的行列编码记为r r，s s，在二维等带宽矩在二维等带宽矩在二维等带宽矩在二维等带宽矩阵阵阵阵K K K K*中的行列编码为中的行列编码为中的行列编码为中的行列编码为r r*，s s*边界条件的处理非节点载荷的移植：非节点载荷的移植：非节点载荷的移植：非节点载荷的移植：集中力：集中力：

29、集中力：集中力：最好在集中力处设置节点最好在集中力处设置节点分布面力：分布面力：分布面力：分布面力：分布体力：分布体力：分布体力：分布体力：边界条件的处理事实上，由于整体刚度矩阵的奇异性，仍然是没有办法求解。事实上，由于整体刚度矩阵的奇异性，仍然是没有办法求解。事实上，由于整体刚度矩阵的奇异性，仍然是没有办法求解。事实上，由于整体刚度矩阵的奇异性，仍然是没有办法求解。原因在于位移边界条件没有引入。原因在于位移边界条件没有引入。原因在于位移边界条件没有引入。原因在于位移边界条件没有引入。方法一：划零置一法方法一：划零置一法方法一：划零置一法方法一：划零置一法若已知边界条件：若已知边界条件：若已知

30、边界条件：若已知边界条件：则则则则边界条件的处理事实上，由于整体刚度矩阵的奇异性，仍然是没有办法求解。事实上，由于整体刚度矩阵的奇异性，仍然是没有办法求解。事实上，由于整体刚度矩阵的奇异性，仍然是没有办法求解。事实上，由于整体刚度矩阵的奇异性，仍然是没有办法求解。原因在于位移边界条件没有引入。原因在于位移边界条件没有引入。原因在于位移边界条件没有引入。原因在于位移边界条件没有引入。方法一：乘大数法方法一：乘大数法方法一：乘大数法方法一：乘大数法若已知边界条件：若已知边界条件：若已知边界条件：若已知边界条件：则则则则M是很大的数，远大于是很大的数，远大于其它元素如其它元素如M=1.0E+30由于

31、由于由于由于计算结果的整理（1 1）求解整体结构平衡方程的结构节点位移）求解整体结构平衡方程的结构节点位移 总刚用一维变带宽存储技术总刚用一维变带宽存储技术 矩阵求解方法可以采用多种方法矩阵求解方法可以采用多种方法（2 2）根据节点位移求单元应力）根据节点位移求单元应力 3 3节点三角形单元为常应力单元节点三角形单元为常应力单元 可以理解应力为单元中心处的应力可以理解应力为单元中心处的应力 结果的整理结果的整理 采用绕节点平均法求节点应力采用绕节点平均法求节点应力 插值法求边界节点应力插值法求边界节点应力 2 2点的应力可由点的应力可由3 3点和点和4 4点插值而得点插值而得教学程序学习阶段教学程序学习阶段