有限元分析教案.doc
《有限元分析教案.doc》由会员分享，可在线阅读，更多相关《有限元分析教案.doc（89页珍藏版）》请在冰点文库上搜索。
有限元分析基础
第一章有限元法概述
在机械设计中，人们常常运用材料力学、结构力学等理论知识分析机械零构件的强度、刚度和稳定性问题。
但对一些复杂的零构件，这种分析常常就必须对其受力状态和边界条件进行简化。
否则力学分析将无法进行。
但这种简化的处理常常导致计算结果与实际相差甚远，有时甚至失去了分析的意义。
所以过去设计经验和类比占有较大比重。
因为这个原因，人们也常常在设计中选择较大的安全系数。
如此也就造成所设计的机械结构整体尺寸和重量偏大，而局部薄弱环节强度和刚度又不足的设计缺陷。
近年来，数值计算机在工程分析上的成功运用，产生了一门全新、高效的工程计算分析学科——有限元分析方法。
该方法彻底改变了传统工程分析中的做法。
使计算精度和计算领域大大改善。
§1.1有限元方法的发展历史、现状和将来
一，历史
有限元法的起源应追溯到上世纪40年代（20世纪40年代）。
1943年R.Courant从数学的角度提出了有限元法的基本观点。
50年代中期在对飞机结构的分析中，诞生了结构分析的矩阵方法。
1960年R.W.Clough在分析弹性力学平面问题时引入了“FiniteElementMethod”这一术语，从而标志着有限元法的思想在力学分析中的广泛推广。
60、70年代计算机技术的发展，极大地促进了有限元法的发展。
具体表现在：
1）由弹性力学的平面问题扩展到空间、板壳问题。
2）由静力平衡问题——稳定性和动力学分析问题。
3）由弹性问题——弹塑性、粘弹性等问题。
二，现状
现在有限元分析法的应用领域已经由开始时的固体力学，扩展到流体力学、传热学和电磁力学等多个传统的领域。
已经形成了一种非常成熟的数值分析计算方法。
大型的商业化有限元分析软件也是层出不穷，如：
SAP系列的代表SAP2000（StructureAnalysisProgram）
美国安世软件公司的ANSYS大型综合有限元分析软件
美国航天航空局的NASTRAN系列软件
除此以外，还有MASTER、ALGO、ABIQUES、ADINA、COSMOS等。
三，将来
有限元的发展方向最终将和CAD的发展相结合。
运用“四个化”可以概括其今后的发展趋势。
那就是：
可视化、集成化、自动化和网络化。
§1.2有限元法的特点
机械零构件的受力分析方法总体说来分为解析法和数值法两大类。
如大家学过的材料力学、结构力学等就是经典的解析力学分析方法。
在这些解析力学方法中，弹性力学的分析方法在数学理论上是最为严谨的一种分析方法。
其解题思路是：
从静力、几何和物理三个方面综合考虑，建立描述弹性体的平衡、应力、应变和位移三者之间的微分方程，然后考虑边界条件，从而求出微分方程的解析解。
其最大的有点就是，严密精确。
缺点就是微分方程的求解困难，很多情况下，无法求解。
数值方法是一种近似的计算方法。
具体又分为“有限差分法”和“有限元法”。
“有限差分法”是将得到的微分方程离散成近似的差分方程。
通过对一系列离散的差分方程求解，得到最终的力学问题近似解。
其优点就是：
计算简单收敛性好。
缺点是：
计算程序无法标准化，在不能获得整个问题的微分方程时，该方法不能运用。
由于其是将微分方程转为差分方程，所以它是一种数学近似。
“有限元法”的基本思想就是“先分后合”或者“化整为零，又积零为整”。
与有限差分不同，它是在力学模型上进行近似处理，也就是（分块近似）。
具体做法：
把连续体模型转为由有限个单元组成的离散体模型，离散体模型之间通过一些节点联系。
对于每一个离散体个体选择简单的函数近似表示其中的物理变化规律（如位移等），运用力学方程推导单元的平衡方程组，然后集合所有的方程组形成表征整体结构的方程组，引入边界条件，求取最后问题的解。
优点：
概念清晰、易于学习理解，适用性强，便于电算化。
缺点：
计算精度受单元划分的影响较大。
§1.3有限元分析的一般过程
为了能够了解有限元分析的全貌，我们就一个简单的例子，来分析一下有限元分析的三个过程：
结构离散化、单元分析、整体分析。
一，结构离散化
在该阶段中，要完成把连续结构的力学模型转变为离散的力学模型。
处理的好坏，直接影响到最后分析结果的正确与否、计算的精度和计算的效率。
根据模型的传力特性和分析的目标，正确选择单元类型。
通常单元分为：
一维单元、二维单元和三维单元。
所谓一维单元就是指所求物理量仅随一个坐标变量而变化的单元。
如桁架、平面刚架和空间刚架单元。
一维单元：
杆单元、梁单元。
二维单元：
三角形单元、四边形单元（平面类问题）
三维单元：
四面体单元、六面体单元等（空间问题）
计算精度和计算效率：
取决于单元划分的形状、大小和分布状况。
通常单元愈多、愈密集，计算精度愈高，但计算效率愈低。
有限元分析工作就是要在精度和效率两者之间做到有机的统一。
二，单元分析
进行单元分析的目的是为了到处表征单元力学特性的“单元刚度矩阵”。
一般说来该过程有三种方法：
1，直接法。
2，虚功原理法（变分法）。
3，加权余数法。
直接法概念浅显，易于理解物理含义。
变分法需要泛函的数学知识，其推导过程具有严谨的数学概念。
加权余数法适用于泛函不存在的应用范围。
本教材将运用虚功原理方程结合弹性力学和材料力学中的知识来推求几种常见单元的单刚计算公式。
现在先看一个简单的阶梯轴的轴向拉伸问题
例：
如图所示的变截面直杆，受拉力P，运用有限元方法分析其变形。
取任意单元，长度为l，面积为A，
单元内任意一点的轴向位移和其位置坐标成正比，即
u=a0+a1x
其中a0、a1为待定系数。
由于杆的两个端点节点1、2是单元上的点，所以它们应该满足上述方程。
节点1，x1=0，∴u1=a0+a1×0=a0
节点2，x2=l，∴u2=a0+a1×la1=（u2-u1）/l
将求出的结果带入方程并整理，就得：
式中：
N1、N2是形函数
[N]形函数矩阵
{δ}e节点位移向量
由位移与应变的关系知道：
将上面推出的位移表达式代入，可得：
上式中的[B]称为应变矩阵或几何矩阵。
运用材料力学中的虎克定律，可以将应变和应力联系起来。
单向应力状态的虎克定律为
××
其中[S]称为应力矩阵。
利用虚功方程可以建立力与位移之间的关系，也就是单刚方程。
在后面我们将会推导出它的一般形式如下：
式中：
{F}e为单元节点力向量，对我们这个例子应为[U1U2]T。
[K]e为单元刚度矩阵。
后面将推导出它的计算公式为
[D]矩阵是弹性矩阵。
对于一维单元来说，就是E。
所以我们这儿讨论的例题：
求得单刚矩阵，也就完成了单元分析。
总结单刚矩阵推导的步骤，应该分为四步：
1）假定单元内位移变化的近似规律，即选择位移模式。
2）运用几何关系，推求位移与应变的关系。
3）应用物理规律，把应变与应力联系起来。
4）运用虚功方程的力与位移关系，求出单刚矩阵。
单元分析是整个有限元分析的核心。
不同的单元因为其力学特性不同，而具有不同的单元刚度矩阵，我们这本教材就是要学习几种常用单元矩阵的推导和计算。
了解各种单元的力学特性，为以后选择单元类型打好基础。
三，整体分析
1，由各单元刚度矩阵组集成整个结构的总刚度矩阵。
整个结构有三个节点，首先将单元刚度矩阵扩充为3X3的矩阵，移动各元素使之与单刚矩阵中的元素位置相对应，如下：
然后直接相加。
2，把各单元的节点力向量组集成总的节点载荷向量。
3，根据边界条件，修改总刚度矩阵，获得总刚方程组。
边界条件修改之前的总刚方程：
修改以后（采用置“0，1”法）
4，求解方程组，得出总的节点位移向量。
解得的解是：
有了节点位移，再回代到前面单元推导过程中的公式×和××，就可以求得每个单元的应变和应力了。
从这个简单的例子，我们了解了有限元法求解力学问题三大步骤中的内容，想必很多同学会说，这样复杂，如果运用材料力学的知识，我还来得快些。
但是大家不要忘记，有限元的计算很多都是编程完成，而且现在很多的商业软件都已经完成了很多的工作。
我们学习有限元主要是了解它的原理，并对常见单元的力学特性有所了解，这样对于以后运用有限元起到帮助作用。
所以下面章节的内容，就是围绕这个主题展开。
要达到这个目的，我们还必须学习必要的弹性力学知识。
对弹性力学知识的学习，也对我们以后把握问题的本质有帮助。
第二章平面问题
平面问题在力学研究的课题中属弹性力学的范畴。
该类问题不仅本身具有典型性，而且在机械零构件的分析中，也是应用得非常广泛。
所以这类问题也称之为经典的力学问题。
我们知道，实际的机械零构件都是具有三维空间尺寸的物体，理应作为三维对象处理，但是当物体的几何形状和受力状态处于某些特定的情况下，近似地简化为平面问题，不仅可以大大简化计算的工作量，而且其精度也完全能够满足所要求。
如：
直齿圆柱齿轮可在垂直与孔轴线的截平面内作平面应力分析就足以了解整个齿轮的受力状态；大坝的横断面可作平面应变分析来了解整个大坝受力情况等。
本章是全书的重点，在这里不仅介绍弹性力学的基本知识。
还将系统地讲解有限元的基本概念、原理和方法。
是学习以后各章节的基础。
§2.1外力、应力、应变和位移
外力、应力、应变和位移的概念在材料力学中已经学习过，由于这些概念在弹性力学、有限元法中具有和在材料力学中不同的规定，弹性力学中的规定和有限元法是完全相同的，所以在这里我们将按照弹性力学的习惯表达方法把他们集中的加以阐述。
一，外力
外界作用在物体上的作用力，可以分为两大类：
1）体积力
分布在物体体积内的力。
如：
重力、惯性力和磁性力等。
单位体积的体力在坐标轴上的分量X、Y、Z，称为体力分量。
符号规定为沿坐标轴正向的为正，沿负向的为负。
2）面力
作用在物体表面的力。
如轮压、水压等。
它又可细分为集中力与分布力。
面力在坐标轴上的投影，表示为X、Y、Z。
符号沿正轴为正，负轴为负。
二，应力
弹性体受到外力作用后，内部产生的抵抗变形的内力。
以弹性体中P点为定点的微单元体来考察。
所谓微单元体，就是图中PA、PB、PC的边长分别为dx、dy和dz。
以下简称这样的微单元体为微元体。
微元体每个面上的应力都可以分解为三个应力分量。
以图中红面为例，分别是σx、τxy、τxz。
应力命名的规则：
以应力所在面垂直的坐标轴为第一个下标，应力指向为第二下标。
如果下标相同就用一个下标表示。
符号规定：
正面上的应力与坐标轴同向为正，反之为负。
负面上的应力与坐标轴反向为正。
反之为负。
所谓正面就是面的外法向与坐标轴同向为正。
反之为负面。
作用在两个相互垂直的面上，并且垂直与该两面交线的剪应力互等。
即：
τxy=τyx；τyz=τzy；τxz=τzx
如此以来，代表P点应力状态的应力分量应有6个，它们是：
三，位移
任一点的位置移动用u、v、w表示它在坐标轴上的三个投影分量。
符号规定：
沿坐标轴正向为正，反之为负。
四，应变
弹性体内各点的位移在受力后一般是不相同的。
各点之间距离的改变，从而使物体形状发生变化，即所谓的变形。
而物体的形状总可以用它各部分的长度和角度表示。
长度的改变称为正应变ε，角度的改变称为剪应变γ。
以微元体三个棱边的线伸长和角度的变化，就分别有和6个应力分量相对应的6个应变分量，即：
为与前面符号规定一致，这里对剪应变的符号规定如下：
正应变伸长为正，缩短为负；剪应变使直角变小为正，变大为负。
§2.2两类平面问题
前面我们讲过，实际受力物体都是三维的空间物体，作用在其上的外力，通常也是一个空间力系，其应力分量、应变分量和位移也都是x、y、z三个变量的函数。
但是当所考察的物体具有某种特殊的形状和特殊的受力状态时，就可以简化为平面问题处理。
弹性力学中的平面问题有两类。
一，平面应力问题
当物体的长度与宽度尺寸，远大于其厚度（高度）尺寸，并且仅受有沿厚度方向均匀分布的、在长度和宽度平面内的力作用时，该物体就可以简化为弹性力学中的平面应力问题。
我们分析以下其应力特征。
当z=±t/2时，有σz=0、τzx=0、τzy=0。
由于板较薄（相对于长度和宽度尺寸），外力沿板厚又是均匀分布的，根据应力应连续的假定（弹性力学中的基本假定），所以可以认为，整个板的各点均有σz=0、τzx=0、τzy=0。
如此以来，描述空间问题的6个应力分量也就变为了3个，即
而且这些应力分量仅是x、y两个变量的函数。
二，平面应变问题
当物体是一个很长、很长的柱形体，其横截面沿长度方向保持不变，物体承受平行于横截面且沿长度方向均匀分布的力时，该问题就可以简化为平面应变问题处理。
分析其应力特征。
假定其长度方向为无限长，那么任一横截面都可以看作是物体的对称面，如此则有该面上的点都有w=0，也就是横截面上的所有点都不会发生Z方向的位移。
由这一点可以推出也就有εz=0、τzx=0、τzy=0。
和平面应力相比较，平面应变是εz=0，那么是否也就有σz=0呢？
可能有同学想σ=Eε，当然也就有σz=0，这是错误的。
平面应变状态下σz≠0的。
虽然不等于零，但它也不是一个独立的变量了，它由σx、σy的大小而决定。
如此以来，独立的应力分量同平面应力问题一样也是3个：
三，两类问题的比较
1，几何特征
平面应力厚度＜＜长度、宽度
平面应变厚度＞＞长度、宽度为便于说明可讲上述长度看作为厚度
2，受力特征
外力都必须在其面内且不沿厚度方向变化
3，应力特征
平面应力σz=0、τzx=0、τzy=0。
εz≠0自由变形（无约束）
平面应变σz≠0但不是自变量、τzx=0、τzy=0。
εz=0
总以上比较可以看出，平面应力是真正的2维（平面）应力状态，而平面应变却不是，而是3维应力状态，只不过σz不是独立变量而是随横截面平面应力分量而定。
独立变化的应力分量只有3个，类似于平面应力状态。
§2.3平衡微分方程
弹性力学求解问题是从静力学、几何学和物理学三方面综合考虑的。
所以我们首先微元体应该满足的平衡条件——平衡微分方程。
我们以平面问题为例推导，看看它应该具有什么形式。
首先对平面问题的微元体进行受力分析图，如左所示。
物体静力平衡的条件是：
∑Fx（y）=0；∑M=0。
先看∑Fx=0
展开化简得
同理可求得∑Fy=0满足得条件，
由∑M=0，列出方程如下：
化简后得：
略去微量项，可得：
。
这就是前面所将的剪应力互等。
对于平面应变问题，微元体的前后面还有正应力σz，不过它们是互等的。
对于推导出来的结果，没有任何影响。
所以平面问题的平衡微分方程就是：
写成矩阵形式
§2.4几何方程
考察平衡微分方程，其中具有三个未知变量σx、σy、τxy，，而只有两个方程，方程具有无数个解。
表明仅从静力学关系无法求解该方程。
我们必须从其它方面寻求帮助。
弹性体在受到外力后，会发生位移和形变，从几何上描述弹性体各点位移于应变之间的关系，就是弹性力学中的又一个重要方程——几何方程。
仍然取截面的微元体ABCD，AB、CD边长为dx、dy，厚度为“1”。
位移u、v都是x、y的函数，即u（x,y）、v（x,y），偏导数、表示位移分量u、v沿坐标轴x的变化率，偏导数、表示位移分量u、v沿坐标轴y的变化率，设A点的位移为u、v，那么B’点的位移就是：
同理的D’点的位移分量
由于α角在位移和形变很微小的情况下非常小，所以
A’B’≈A’B”
线段AB位移后的总伸长量为
A’B’-AB=A’B”-AB=uB’-uA=-u=
∴，同理可得
剪应变由α、β两个角度组成
由于，所以，同理可得
∴
综合以上几何方程，并将它们写成矩阵形式：
由以上方程可以看出，当弹性体的位移分量确定以后，由几何方程可以完全确定应变，反过来，已知应变却不能完全确定弹性体的位移。
这是因为物体产生位移的原因有两点：
1）变形产生的位移。
2）因运动产生的位移。
因此弹性体有位移不一定有应变，有应变就一定有位移。
§2.5物理方程
描述弹性体内应力与应变关系的方程，我们称之为物理方程，也叫材料的本构方程。
弹性力学通常研究的是各向同性材料，在三维应力状态下的应力应变关系。
当弹性体处于小变形条件下，正应力只会引起微元体各棱边的伸长或缩短，而不会影响棱边之间角度的变化，剪应力只会引起角度的变化而不会引起各棱边的伸长或缩短。
因此运用力的叠加原理、单向虎克定律和材料的横向效应（泊松效应），我们就可以很容易的推导出材料在三向应力状态下的虎克定律，也就是通常所说的广义虎克定律。
式中，E——材料线弹性模量
G——材料剪切弹性模量
μ——材料横向收缩系数，即泊松系数。
三者不是独立的。
具有以下关系：
这些参数都是材料的固有属性系数，可以通过查材料手册获得。
如：
钢材的弹性模量E＝196～206GPa之间，通常取2.1×105MPa，μ＝0.24～0.28之间，也取为0.3进行计算，G＝79GPa。
将以上空间问题的物理方程运用到平面问题，其形式如下：
1）平面应力问题的物理方程
前面分析已知，平面应力问题有σz=0、τzx=0、τzy=0。
所以：
从以上物理方程，也论证了我们前面说的εz≠0的结论，但由于它是由x和y方向应力产生的附加无约束变形，所以通常不予以考虑。
在有限元分析中更多的是运用应变表示的应力关系，所以我们将上式变形一下：
以上方程的矩阵表达形式为：
简记为：
式中：
{σ}、{ε}为该问题的应力、应变向量。
[D]为弹性矩阵。
它是一个对称矩阵，且只与材料的弹性常数有关。
2）平面应变问题的物理方程
因为εz=0所以由空间物理方程的第三式得：
，代入
（1）、
（2）式得
同理变型为应变表示应力的形式：
矩阵形式：
也可简记为
平面应变问题的弹性矩阵不同于平面应力问题的弹性矩阵，比较可以发现只需将平面应力问题弹性矩阵[D]中的材料常数E换为E/（1-μ2），μ换为μ/（1-μ）就得到了平面应变问题的弹性矩阵。
其实弹性矩阵的这种转换方法，是弹性力学中将平面应力结果，转换到平面应变问题结论的一般方法。
因为在两种平面问题的描述方程中（平衡微分方程、几何方程和物理方程），只有物理方程是不同的。
§2.6边界条件
求解弹性力学问题实际就是在确定边界条件下，求解8个基本方程（平面问题而言），以确定8个未知变量。
所以从数学的角度看，就是求解偏微分方程的边值问题。
边界条件的给出通常是各式各样的。
大体可以分为三类：
1，第一类边值问题
给定物体的体力和面力条件，确定弹性体的应力场和位移场。
此类问题边界以力的形式给出，所以也称为应力边界条件。
我们可以来考察一下应力边界的一般形式：
是在Sσ面上给出的力的分量。
平面问题如左图所示，设阴影部分的微元体弧长为ds，厚度为单元厚度“1”，其法线与X轴的夹角为θ，由阴影部分微元体的平衡条件可以推出：
化简后得：
此即为平面问题应力边界方程。
2，第二类边值问题
给出弹性体的体力和物体表面各点得位移条件，确定弹性体得应力场和位移场。
由于以位移给出已知得边界条件，所以也称为位移边界问题。
一般得位移边界条件为：
在Su面上。
3，第三类边值问题
给定弹性体得体力和一定边界上得面力，其余边界上的位移，确定其应力场和位移场。
由于边界以力和位移两种形式给出，所以也称为混合边界问题。
针对不同的边界条件，弹性力学求解的方法也有所不同。
§2.7弹性力学的解题方法（解析法）
1，应力法
由于第一类问题的边界条件以应力形式给出，所以以应力作为基本的未知量的求解过程，就是人们通常所说的应力法。
由于平衡方程中有三个未知量，而只有两个平衡微分方程，必须找出另外一个包含应力分量的方程，才能求得方程解。
考虑到弹性体变形前是一个连续体，变形后也应是连续体的基本假设，所以要求微元体的变形一定要协调，才能使变形前、后，不会发生裂缝、重叠等现象。
要使变形协调，就要研究几何方程。
前面介绍的平面问题几何方程如下：
分别对εx、εy、求y、x的二阶偏导，然后相加：
上式表明三个应变分量之间应满足的连续性条件，我们称之为形变协调方程（相容方程）。
通过物理方程，使上述的形变协调方程换成应力表示的形式，使之与平衡微分方程就构成了应力法中需要求解的方程组。
具体我们来看看：
1）利用物理方程消去相容方程中的形变分量（以平面应力为例）
2）利用平衡微分方程，消去上述公式中的剪应力
式对X求偏导，式对y求偏导，然后两者相加
代入相容方程，化简
对于平面应变而言，运用前面讲过的物理方程的转换方法，只需将上式中的μ代以μ/（1-μ）就可以了。
3）最终求解的方程组
平面应力问题
平面应变问题
三个微分方程，三个未知变量，再考虑边界条件，即可求得。
如果是单连体（只具有唯一的封闭边界）的对象，满足了以上方程组后就是实际的解。
但对于多连体（具有多个封闭边界）的对象，中还包含有待定系数，这些待定系数会导致位移的解出现多值性。
所以对于多连体的问题，还应考虑位移的单值条件，才能最终确定。
该部分的内容可以参见徐芝纶编的《简明弹性力学教程》中圆环受均布压应力的情况（P87）
2，位移法
位移法主要针对第二类边界条件问题求解。
解题步骤：
1）改写物理方程使之成为应变表示应力的形式
2）应用几何方程表示以上得到公式中的应变
3）将它们代入平衡微分方程
经整理最后得到的位移法求解平面应力问题方程为：
两个未知量，两个方程，再加以边界条件即可求得问题的解。
以上介绍的解析法中，应力法和位移法是求解弹性力学问题的基本方法。
但都需要解联立的偏微分方程组。
求解过程中的数学难度，常常导致这种求解是无法进行的。
由于应力法在体力为常量的情况下可以进一步简化为求解一个单独的微分方程的问题，所以应力法在解析法中相对应用较多。
但即使这样，在应力法中，也常常采用逆解法或半逆解法。
§2.8常体力情况下应力法的简化、应力函数及实例分析
我们前面讲述了弹性力学的三大方程，及应用这三大方程的应力法和位移法解题步骤。
但是也说了要求解这些联立的偏微分方程在数学上是存在很大难度的。
很多情况下，根本无法进行。
那么弹性力学如何在实际中进行应用，它们和我们前面学过的材料力学区别究竟在哪里？
我们将通过这一节的学习，一方面了解如何应用这些弹性力学的方程求解问题，另一方面加深对力学概念的理解，建立力学分析问题的直观感觉，为建立有限元模型打好基础。
我们知道在大多数情况下，我们分析的对象，体力是常数，它不随x、y坐标变化。
如此以来，前面讲解的第三个方程（应力表示的相容方程），就可以简化为了：
简记为：
以上方程称为拉普拉斯微分方程，数学上也称之为调和方程，满足调和方程的函数称之为调和函数，及这里的。
是拉普拉斯算子。
这样以来常体力情况下的应力法方程就是：
以上方程都不含有材料常数E、μ，所以平面应力和平面应变两类问题具有相同的方程，这表明：
在单连体问题中，只要边界相同、受同样的分布外力，应力分布与材料无关；也与是平面应力还是平面应变的状态无关。
以上结论的意义：
弹性力学平面解答的应用范围加宽。
为实验应力分析提供了理论依据（光弹实验）
下面我们考察平衡方程：
其解由其次
升级会员 