第3有限元分析的数学求解原理

1、第3有限元分析的数学求解原理第1页，共62页，2022年，5月20日，7点49分，星期六前一章针对任意形状变形体，基于物体内的微小体元dxdydz定义了描述弹性变形体的所有基本力学信息（ui，ij，ij）、基本方程（平衡、几何、物理）及边界条件。接下来的任务就是对这些方程在具体的条件下进行求解，也就是说在已知边界条件下，由基本方程求出相应的位移场、应力场和应变场。一般来说，求解方程的途径有两大类：（1）直接针对原始方程进行求解，方法有：解析法(analytical method)、半逆解法(semi-inverse method)、有限差分法(finite difference method)

2、等；（2）间接针对原始方程进行求解，方法有：加权残值法、虚功原理、最小势能原理、变分方法等第2页，共62页，2022年，5月20日，7点49分，星期六主要内容3.1 简单问题的解析求解3.2 虚功原理3.3 应用举例3.4 基本步骤第3页，共62页，2022年，5月20日，7点49分，星期六3.1 简单问题的解析求解1D拉杆问题有一个左端固定的拉杆，其右端承受一外力P。该拉杆的长度为l，横街面积为A，弹性模量为E，如图所示：第4页，共62页，2022年，5月20日，7点49分，星期六(1)基本变量由于该问题视为沿x方向的一维问题，因此只有沿x方向的变量，而其它变量为零。即： 位移：u(x) 应

3、变：x(x) 应力：x(x)3.1 简单问题的解析求解第5页，共62页，2022年，5月20日，7点49分，星期六(2)基本方程对原三维问题的所有基本方程进行简化，只保留沿x方向的方程，得到该问题的三大类基本方程和边界条件平衡方程（无体力）几何方程物理方程3.1 简单问题的解析求解第6页，共62页，2022年，5月20日，7点49分，星期六3.1 简单问题的解析求解边界条件（BC）上述方程中，力的边界条件为一种近似，因为在x=l的端面，x(x)不应是均匀分布的。由圣维南原理（Saint-Venant principle），在远离x=l的截面，力的边界条件才较好的满足。第7页，共62页，2022

4、年，5月20日，7点49分，星期六(3)求解对上述的方程直接求解，可以得到以下的结果：3.1 简单问题的解析求解其中c及c1为待定系数，由边界条件可以求出上式中的常数为c1=0，c=P/A，因此有最后的结果：第8页，共62页，2022年，5月20日，7点49分，星期六(4) 讨论1上述问题若用经验方法求解（如材料力学的方法），则需要先作平面假设，即假设x为均匀分布，这样可以得到再由Hooke定律算出：再计算右端的伸长量为：3.1 简单问题的解析求解第9页，共62页，2022年，5月20日，7点49分，星期六通过比较可以看出，经验方法求解的结果和弹性力学的解析结果完全一致。比较以上解析方法和经验

5、可以看出：解析方法的求解过程严谨，可以得到物体内各点力学变量的表达，是场变量。经验方法的求解过程比较简单，但需要事先进行假定，往往只能得到一些特定位置的力学变量表达，而且只能应用于一些简单情形。3.1 简单问题的解析求解第10页，共62页，2022年，5月20日，7点49分，星期六(5) 讨论2根据计算能量的方法，得到：应变能外力功势能3.1 简单问题的解析求解第11页，共62页，2022年，5月20日，7点49分，星期六平面梁的弯曲问题 假设有一个受分布载荷作用的简支梁如图所示，由于简支梁的厚度较小，外载沿厚度方向无变化，那么该问题可以认为是一个oxy平面内的问题3.1 简单问题的解析求解第

6、12页，共62页，2022年，5月20日，7点49分，星期六1.基本方程 有两种方法来建立基本方程。 方法一：采用一般建模及分析方法，即从对象取出dxdy微元体进行分析，建立最一般的基于（ui，ij，ij）描述的方程，类似于2D问题的基本变量及方程，这样，所用的变量较多，方程复杂。 方法二：采用特征建模（characterized modeling）的简化方法来推导的三大方程，其基本思想是采用工程宏观特征量来进行问题的描述。3.1 简单问题的解析求解第13页，共62页，2022年，5月20日，7点49分，星期六应此简支梁问题的特征为：梁为细长梁(long beam)，因此可只用x坐标来刻画；主

7、要形变为垂直于x的挠度，可只用挠度（deflection）来描述位移场。3.1 简单问题的解析求解补充概念：挠度（ deflection ），弯曲变形时横截面形心沿与轴线垂直方向的线位移称之为挠度。简言之，就是指梁、桁架等受弯构件在载荷作用下的最大变形，通常指竖向，就是构件的竖向变形第14页，共62页，2022年，5月20日，7点49分，星期六针对这两个特征，可以做出以下假定：直法线假定小变形与平面假定该问题的三类基本变量：位移： （中层性挠度）应力：（采用 x，其他应力分量很小，不考虑），该变量对应于梁截面上的弯矩M应变：（采用x，满足直线假设）3.1 简单问题的解析求解第15页，共62页，

8、2022年，5月20日，7点49分，星期六下面取具有全高度梁的dx“微段”来推导三大类方程3.1 简单问题的解析求解平衡方程首先是x方向的合力平衡然后是y方向的合力平衡y为距梁中性层的坐标最后是弯矩平衡第16页，共62页，2022年，5月20日，7点49分，星期六几何方程 由变形后的几何关系，可得到其中，y为距中性层的坐标，k为梁挠度的曲率，即：3.1 简单问题的解析求解物理方程 由Hooke定律有：第17页，共62页，2022年，5月20日，7点49分，星期六对上述方程整理，就得到了平面简支梁弯曲问题的基本方程：3.1 简单问题的解析求解式中， 为梁截面的惯性矩（moment of iner

9、tia）。可以看出：将原始基本变量定位中性层的挠度v(x),则可以求出其他参数第18页，共62页，2022年，5月20日，7点49分，星期六边界条件 该简支梁的边界为梁的两端，作用在梁上的q(x)已在平衡方程中考虑，因此不作为力的边界条件。两端位移：两端力（弯矩）：将弯矩以挠度的二阶导数来表示，即：3.1 简单问题的解析求解第19页，共62页，2022年，5月20日，7点49分，星期六2.求解 若用基于dxdy微体所建立的原始方程（即原平面应力问题中的三大类方程）进行直接求解，比较麻烦，并且很困难，若用基于简化的“特征建模”方法所得到的基本方程进行直接求解则比较简单，对简支梁问题求解，其方程为

10、：3.1 简单问题的解析求解第20页，共62页，2022年，5月20日，7点49分，星期六对上述的常微分方程，其解的形式有：其中c0c3为待定系数，可由四个边界条件BC求出，有结果：3.1 简单问题的解析求解第21页，共62页，2022年，5月20日，7点49分，星期六3.2 虚功原理虚位移与虚功原理如图所示的平衡力系，由于该系统处于平衡状态，则有：假想在该平衡力系上作用有微小的扰动（不影响原平衡条件），且外力所作用的位置产生了微小的位移变化，即A ，B。该假想的位移如果不影响原平衡条件，应满足以下几何关系：第22页，共62页，2022年，5月20日，7点49分，星期六 进一步分析。当杠杆处于

11、平衡状态时，A和B这两个位移是不存在的，但是如果某种原因，例如人为地振一下让它倾斜，一定满足这种关系。 将这个客观存在的关系抽象成一个普遍的原理，去指导分析和计算结构。 对于在力的作用下处于平衡状态的任何物体，不用考虑它是否真正发生了位移，而假想它发生了位移，(由于是假想，故称为虚位移)，那么，物体上所有的力在这个虚位移上的总功必定等于零。这就叫做虚位移原理，也称虚功原理。 在图a中的PA和PB所作的功就不是发生在它本身(状态a)的位移上，(因为它本身是平衡的，不存在位移)，而是在状态(b)的位移上作的功。可见，这个位移对于状态(a)来说就是虚位移，亦即是状态(a)假象的位移。3.2 虚功原理

12、第23页，共62页，2022年，5月20日，7点49分，星期六 必须指出，虚功原理的应用范围是有条件的，它所涉及到的两个方面，力和位移并不是随意的。对于力来讲，它必须是在位移过程中处于平衡的力系；对于位移来讲，虽然是虚位移，但并不是可以任意发生的。它必须是和约束条件相符合的微小的刚体位移。 还要注意，当位移是在某个约束条件下发生时，则在该约束力方向的位移应为零，因而该约束力所作的虚功也应为零。这时该约束力叫做被动力。(如图中的反力Rc由于支点C没有位移，故Rc所作的虚功对于零)。反之，如图中的PA和PB是在位移过程中作功的力，称为主动力。因此，在平衡力系中应当分清楚哪些是主动力，哪些是被动力，

13、而在写虚功方程时，只有主动力作虚功，而被动力是不作虚功的。3.2 虚功原理第24页，共62页，2022年，5月20日，7点49分，星期六虚功原理表述如下： 在力的作用下处于平衡状态的体系，当发生与约束条件相符合的任意微小的刚体位移时，体系上所有的主动力在位移上所作的总功(各力所作的功的代数和)恒对于零。虚功原理用公式表示为：这就是虚功方程，其中P和相应的代表力和虚位移。3.2 虚功原理第25页，共62页，2022年，5月20日，7点49分，星期六虚功方程是按刚体的情况得出的，即假设图示杠杆是绝对刚性，没有任何的变形，因而在方程中没有内功项出现，而只有外功项。将虚功原理用于弹性变形时，总功W要包

14、括外力功(T)和内力功(U)两部分，即： W = T - U ；内力功(-U)前面有一负号，是由于弹性体在变形过程中，内力是克服变形而产生的，所有内力的方向总是与变形的方向相反，所以内力功取负值。根据虚功原理，总功等于零得： T - U = 0 外力虚功 T = 内力虚功 U 弹性力学中的虚功原理可表达为：在外力作用下处于平衡状态的弹性体，如果发生了虚位移，那么所有的外力在虚位移上的虚功(外力功)等于整个弹性体内应力在虚应变上的虚功(内力功)。3.2 虚功原理第26页，共62页，2022年，5月20日，7点49分，星期六虚应变能虚应变分量外力虚功内力虚功即应力在虚应变上做的的虚功，也称虚应变能

15、外力虚功即作用于弹性体上的外力在虚位移上做的功由于虚位移是微小的，可认为在虚位移发生过程中外力保持为常量，则上式的变分符号可提到积分号外。3.2 虚功原理第27页，共62页，2022年，5月20日，7点49分，星期六最小势能原理： 表明在满足位移边界条件的所有可能位移中，实际发生的位移使弹性体的势能最小。即对于稳定平衡状态，实际发生的位移使弹性体总势能取极小值。显然，最小势能原理与虚功原理完全等价。最小势能原理 在有限元的理论中，最小势能原理是在所有满足给定边界条件的位移时，满足平衡微分方程的位移使得势能取得最小值。 最小势能原理就是说当一个体系的势能最小时，系统会处于稳定平衡状态。或者说在所

16、有几何可能位移中，真实位移使得总势能取最小值 3.2 虚功原理第28页，共62页，2022年，5月20日，7点49分，星期六虚功原理的矩阵表示i点外力分量j点外力分量外力分量用 表示；引起的应力分量用 表示3.2 虚功原理第29页，共62页，2022年，5月20日，7点49分，星期六假设发生了虚位移虚位移分量为用 表示；引起的虚应变分量 用 表示虚功原理的矩阵表示3.2 虚功原理第30页，共62页，2022年，5月20日，7点49分，星期六虚功原理的矩阵表示在虚位移发生时，外力在虚位移上的虚功是：式中 是 的转置矩阵。 同样，在虚位移发生时，在弹性体单位体积内，应力在虚应变上的虚功是：因此，在

17、整个弹性体内，应力在虚应变上的虚功是：根据虚功原理得到： 这就是弹性变形体的虚功方程，它通过虚位移和虚应变表明外力与应力之间的关系。这是以后推导有限元方程的基础。第31页，共62页，2022年，5月20日，7点49分，星期六3.3 应用实例A1, l1A2, l2R31. 离散化1232. 位移函数A, lijF i uiF j uj第32页，共62页，2022年，5月20日，7点49分，星期六xx2. 位移函数A, lijF i uiF j uj（单元内位移线性分布） 形函数矩阵 3.3 应用实例第33页，共62页，2022年，5月20日，7点49分，星期六343. 单元刚度矩阵方程A 虚功

18、原理 外力虚功虚应变能应变应力应变矩阵弹性矩阵3.3 应用实例第34页，共62页，2022年，5月20日，7点49分，星期六单元刚度矩阵单元的刚度方程单元刚度矩阵ElementElement3.3 应用实例第35页，共62页，2022年，5月20日，7点49分，星期六36B 最小势能定理 外力虚功虚应变能由势能变分原理（势能最小原理）得势能变分，整理得平衡方程3.3 应用实例第36页，共62页，2022年，5月20日，7点49分，星期六4 整体分析 整体分析就是建立整个离散结构所有节点位移与外力之间的关系，实现未知节点位移的求解整体平衡方程整体刚度方程可基于势能变分原理建立，也可根据节点的静力

19、平衡来实现（即每个节点静力平衡）。节点i的平衡为三个节点三个自由度，即3.3 应用实例第37页，共62页，2022年，5月20日，7点49分，星期六（A）扩充单元刚度方程法位移协调性载荷的叠加性A1, l1A2, l2R31233.3 应用实例第38页，共62页，2022年，5月20日，7点49分，星期六整体刚度方程3.3 应用实例第39页，共62页，2022年，5月20日，7点49分，星期六A1, l1A2, l2P123（B）“对号入座”法 （方便编程） i j 1 2 2 3Total1 2 31 2 33.3 应用实例第40页，共62页，2022年，5月20日，7点49分，星期六5.

20、引入边界条件求解A1, l1A2, l2R3123边界条件支反力3.3 应用实例第41页，共62页，2022年，5月20日，7点49分，星期六结构离散单元分析整体分析3.4 基本步骤为三大步骤第42页，共62页，2022年，5月20日，7点49分，星期六1、结构离散:就是用假想的线或面将连续物体分割成有限个单元组成的集合体且单元之间仅在节点处连接，单元之间的作用仅由节点传递。（基本要求）注意的问题单元：满足一定几何特性和物理特性的最小结构域节点：单元与单元间的连接点节点力：单元与单元间通过节点的相互作用力节点载荷：作用于节点上的外载3.4 基本步骤为三大步骤第43页，共62页，2022年，5月

21、20日，7点49分，星期六2单元分析:1），选择插值（位移）函数；2），构造位移函数。插值函数：用以表示单元内物理量变化（如位移或位移场）的近似函数。由于该近似函数常由单元节点物理量值插值构成，故称为插值函数，如单元内物理量为位移，则该函数称为位移函数。选择位移函数的一般原则位移函数在单元节点的值应等于节点位移（即单元内部是连续的）；所选位移函数必须保证有限元的解收敛于真实解。位移函数一般采用多项式形式，在单元内选适当阶次的多项式可得到与真实解接近的近似解3.4 基本步骤为三大步骤第44页，共62页，2022年，5月20日，7点49分，星期六构造位移函数：如平面问题位移函数的一般形式为1.多项

22、式项数越多，则逼近真实位移的精度越高，项数的多少由单元的自由度数决定。2多项式选取应由低阶到高阶，尽量选择完全多项式以提高单元精度。3.选取多项式时，还应使所选取的多项式具有坐标的对称性，即按Pascal（帕斯卡）三角形来选择3.4 基本步骤为三大步骤第45页，共62页，2022年，5月20日，7点49分，星期六位移函数构造方法：1.广义坐标法：2.插值函数法：即将位移函数表示为各个节点位移与已知插值基函数积的和3.4 基本步骤为三大步骤第46页，共62页，2022年，5月20日，7点49分，星期六3.单元特性分析：单元特性分析的基本任务就是建立单元的平衡方程，也称为刚度方程。在选择了单元类型

23、和相应的位移函数后，即可按弹性力学的几何方程、物理方程导出单元应变与应力的表达式，最后利用虚位移原理或最小势能原理或直接法或加权残值法建立单元的平衡方程，即单元节点力与节点位移间的关系。3.4 基本步骤为三大步骤第47页，共62页，2022年，5月20日，7点49分，星期六3.4 基本步骤为三大步骤第48页，共62页，2022年，5月20日，7点49分，星期六整体分析整体分析的基本任务包括建立整体平衡方程，引入边界条件，完成整体方程求解。整体平衡方程的建立有多种方法，可基于能量原理（势能变分或虚位移原理）推导，也可基于节点力平衡得到。在引入边界条件之前，整体平衡方程是奇异的，这意味着整体方程是

24、不可解的。方程求解包括边界条件引入和数值计算，一旦利用适当的数值方法求出未知的节点位移，则可按前述的应力应变公式计算出各个单元的应变、应力等物理量。3.4 基本步骤为三大步骤第49页，共62页，2022年，5月20日，7点49分，星期六刚度由使其产生单位变形所需的外力值来量度，刚度是指零件在载荷作用下抵抗弹性变形的能力。单元的刚度矩阵：单元刚度矩阵反应的是单元节点力与单元节点位移的关系；总刚度矩阵反应的是整体的节点力与节点位移的关系；刚度矩阵将总体坐标下的节点位移与整个结构的总体力联系在一起。补充实例（刚度矩阵的理解）第50页，共62页，2022年，5月20日，7点49分，星期六单元刚度矩阵行

25、数等于位移向量的分量个数，列数等于为位移的列向量的分量个数，由于两者相等所以单刚是个方阵。结构的总体刚度矩阵即结构的原始刚度矩阵，每1个元素的物理意义就是当其所在列对应的节点位移分量等于单位位移（其余结点位移分量为0）时，其所在行对应的节点力的数值。 表示由于第j个自由度的单位位移dj在第i个自由度需要的力补充实例（刚度矩阵的理解）第51页，共62页，2022年，5月20日，7点49分，星期六单元的刚度矩阵的性质 a. 单元刚度矩阵仅与单元的几何特征和材料性质有关。仅与单元的横截面积A、惯性矩I、单元长度l、单元的弹性模量E有关。 b. 单元刚度矩阵是一个对称阵。在单元刚度矩阵对角线两侧对称位

26、置上的两个元素数值相等，即，根据是反力互等定理。 c. 单元刚度矩阵是一个奇异阵。 d. 单元刚度矩阵可以分块矩阵的形式表示。具有确定的物理意义。第52页，共62页，2022年，5月20日，7点49分，星期六整体刚度矩阵的性质 整体刚度矩阵 中位于主对角线上的子块 ，称为主子块，其余 为副子块。 a. 中主子块 由结点i的各相关单元的主子块扩展之后叠加求得，即 b. 当结点i、 j为单元e的相关结点时， 中副子块 为该单元e相应的副子块，即 。 c. 当结点i、 j为非相关结点时， 中副子块 为零子块，即 。 d. 仅与各单元的几何特性、材料特性，即A、I、l、E等因素有关。 e. 为对称方阵， f. 为奇异矩阵，其逆矩阵不存在，因为建立整体刚度矩阵时没有考虑结构的边界约束条件。第53页，共62页，2022年，5月20日，7点49分，星期六 g. 为稀疏矩阵，整体刚度矩阵中的非零元素分布区域的宽度与结点编号有关，非零元素分布在以对角线为中心的带状区域内，称为带状分布规律，见图a。在包括对角线元素在内的区域中，每行所具有的元素个数叫做把半带宽，以d表示。最大半带宽等于相邻结点号的最大差值加 1 与结点自由度数的乘积，结点号差越大半带宽也就越大。计算机以半带宽方式存储，见图b。半带宽越窄，计算机的存储量就越少，而且可以大幅度减少求解方程所需的运算次数。其效果对大型结构