用于振动分析的有限元方法课件

用于振动分析的有限元方法指导老师：陈益报告人：成志斌韩宗彪何瑜宁鹏用于振动分析的有限元方法指导老师：陈益内容有限元介绍单个元素的运动方程整个系统的运动方程整个系统的边界条件的加载及质量矩阵

MATLAB实例及总结单个元素的质量矩阵、刚度矩阵、力矢量及其转化内容有限元介绍单个元素的运动方程整个系统的运动方程整个系统的有限元法简介有限元法是一种可用于精确地(但近似)解决许多复杂的振动问题的数值方法。对于基本的一维元素进行有限元分析，能得到质量矩阵与刚度矩阵和所需的力矢量，对于二维三维，元素矩阵会转换成相关的更高维的空间。使用一致的和集中质量矩阵的有限元方程并结合边界条件能为复杂系统提供解释。最后,使用MATLAB程序得到在轴向载荷下的指定节点位移,固有振动频率和特征值分析。有限元法简介有限元法是一种可用于精确地(但近似本章目的*认识用于解决不同类型振动问题的刚度和质量矩阵。*将矩阵元素从局部坐标系变换到全球坐标系。*装配单元矩阵和应用边界条件。*对杆、梁元素进行静态分析。*对杆、梁元素进行动态分析来得到固有频率和振型。*在有限元振动分析使用一致的集中质量矩阵。*使用MATLAB解决振动问题。本章目的*认识用于解决不同类型振动问题的刚度和质量矩阵。有限元思想1，实际结构被一些元素所取代,这些元素都是被假定为一个连续的结构部件即有限元，这些元素在特定点即节点上互相关联。2，如果解决方案的各方面都选择得当,那么它可以收敛到精确的解决方案，因为组成总体结构的元素很小，在节点上的力的平衡和元素之间的位移都令人感到满意,这样整个结构(组合的元素)表现为单一实体。3，因为得到准确解很难，所以得到一个方便且逼近的近似解很有价值。有限元思想1，实际结构被一些元素所取代,这些元素都是被假定为元素的运动方程龙门刨铣床有限元模型三角板元素梁元素元素的运动方程龙门刨铣床有限元模型三角板元素梁元素元素的运动方程位移函数形状函数各点对应位移未知节点位移数n动能应变能质量矩阵刚度矩阵位移函数形状函数各点对应位移未知节点位移数n动能应变能元素的运动方程位移函数质量矩阵刚度矩阵位移函数主要内容：一，单元的质量、刚度矩阵、等效节点力二，单元矩阵的坐标变换三，整个系统的运动方程主要内容：一单元的质量、刚度矩阵，等效节点力矢量一杆单元一个杆单元是从杆上划分出的一个小段，如下图所示。由于单元很小，ρ、A均视为常量。现在就以这最简单的杆单元，推导出它的质量、刚度矩阵，等效节点力。图12.1图12.1一单元的质量、刚度矩阵，等效节点力矢量一杆单元一个杆单元(1)求杆单元上任意点的位移u(x,t)本来，杆单元上任意点的位移u(x,t)与节点的位移u1(t)、u2(t)之间的关系是未知的，但是，只要单元划分的足够小，那么其间的关系就无关大局。所以可以假定它们之间有简单的线性关系，即根据节点位移对单元内任意点位移进行插值：（1）式中，Φ1、Φ2称为线性系数，与单元里点的位置有关，是x的函数。此函数与单元的形状有关，又叫形状函数。(1)求杆单元上任意点的位移u(x,t)本来，杆单元上任意形状函数和插值函数一样是任意的，但必须边界条件：只有满足此条件单元才能协调一致运动，而不致破坏系统的完整性，因此这两个条件实际上就是变形协调条件。将式（1）带入（2）中，就可以得到形状函数Φ1(x)、Φ2(x)所满足的边界条件：（2）（3）形状函数和插值函数一样是任意的，但必须边界条件：只有满足此条以上边界条件确定了由于这两个函数的任意性我们可以用简单的线性函数来近似，因此有：（4）代回（1）式中有：（5）为此，我们已经找到了用节点位移表示单元内任意一点位移的表达式。以上边界条件确定了由于这两个函数的任意性我们可以用简单的线性（2）计算此单元的动能和势能杆单元的动能可表示成：（6）（2）计算此单元的动能和势能杆单元的动能可表示成：（6）上式中，ρ是材料的密度，A是杆单元的横截面积。用矩阵形式表示（6）式为：（7）其中，（8）所以，质量矩阵可以认为是：（9）上式中，ρ是材料的密度，A是杆单元的横截面积。用矩阵形式表示杆单元的势能可以写成：（10）式中，E是弹性模量，（10）表示成矩阵形式为：杆单元的势能可以写成：（10）式中，E是弹性模量，（10）表这里，，所以刚度矩阵[k]可以表示成：（11）（12）(3)计算等效节点力设单元上x处作用有分布力f(x,t)，现在要把它等效成节点力遵循等效原则，即原载荷和等效之后的节点载荷在虚位移上所做的虚功相等。这里，，所以刚度矩阵[k]可以表示成：（11）（12）(3)其实，就是对应于广义坐标的广义力，为此，计算所做的虚功：把上式写成矩阵形式：（13）其实，就是对应于广义坐标的广义力，为此，计算所做的虚功：把上所以等效节点力可以写成：（14）所以等效节点力可以写成：（14）二梁单元如下图所示，一个梁单元也是有两个节点，但是有四个自由度，每个节点处，有两种位移形式，一个是线位移，即挠度，一种是角位移。图中，是力，是力矩。是分布载荷是对应的线位移，是对应的转角。是梁单元上任意位移x处的挠度。图12.2二梁单元如下图所示，一个梁单元也是有两个节点，但是有四个自在静载弯曲条件下，梁单元上任意点出的挠度是x的三次方程，可写成：此方程必须满足下面的边界条件：由此可以求解处a(t)、b(t)、c(t)、d(t)，进而挠度方程为：（16）（17）（18）在静载弯曲条件下，梁单元上任意点出的挠度是x的三次方程，可写上式可以写成形状函数的表示：其中，形函数分别为：梁单元的动能、势能、虚功表达式分别为：（19）上式可以写成形状函数的表示：其中，形函数分别为：梁单元的动能式中I是横截面的惯性矩上式中：（20）（21）（22）式中I是横截面的惯性矩上式中：（20）（21）（22）通过上式，可以得到梁单元的质量、刚度矩阵，等效节点力：通过上式，可以得到梁单元的质量、刚度矩阵，等效节点力：二单元矩阵的坐标变换局部坐标系：以各个单元本身的轴线为基准所设立的坐标系。便于计算节点位移。缺点：如果整个系统里各个单元取向各异，各个节点位移方向不一致，如下图。如何使汇交于一个节点的各个杆件的节点位移真正相等?解决方法：进行坐标变换如右图的系统，有四个杆件，u1(t)、u2(t)为局部坐标系的节点位移，Ui

为全局坐标系下的位移图12.3二单元矩阵的坐标变换局部坐标系：以各个单元本身的轴线为基准如下图，节点位移在局部、全局坐标系中的关系：坐标变换矩阵（23）如下图，节点位移在局部、全局坐标系中的关系：坐标变换矩阵（2其中，因为单元的动能、势能与坐标系无关:(24)其中，因为单元的动能、势能与坐标系无关:(24)得到在全局坐标下的单元质量、刚度矩阵为：类似地，根据单元在两个坐标系下的力所做的虚功相等：得到在全局坐标系下的等效节点力：得到在全局坐标下的单元质量、刚度矩阵为：类似地，根据单元在两三全系统运动方程经过坐标变换，各个单元的节点位移方向被统一起来，但是不同的节点有不同的节点位移，为了便于综合出全系统的运动方程，首先要建立全系统的节点位移向量。每个单元的节点位移向量与全系统的节点位移向量之间的关系：长方形矩阵由1和0组成单元节点位移向量单元节点位移向量三全系统运动方程经过坐标变换，各个单元的节

例如，图的12.5中的1单元，方程，变为：把每个单元的动能相加，就得到了整个系统的动能：（把整个系统的动能表示成关于节点速度矢量的形式）例如，图的12.5中的1单元，方程，变为：把每个单元的动能这样就得到了整个系统的质量矩阵:类似地，考虑整个系统的势能，便可以得到整个系统的刚度矩阵：整个系统的广义力向量：最后得到整个系统的运动方程：这样就得到了整个系统的质量矩阵:类似地，考虑整个系统的势能，12.6——添加边界条件前文中，节点没有固定，结构在节点力的作用下会发生刚体位移。也就是说，矩阵[K]是奇异矩阵。通常情况下，我们希望结构的位移为零。

因此，我们需要添加边界条件对矩阵[M]、[K]和向量F进行约束。N：结构中自由节点位移的数目12.6——添加边界条件前文中，节点没有固定，结构例1：杆件分析如图：均质；长0.5m；断面截面积5e-4m^2；杨氏模量200GPa；密度7850Kg/m^3；左端固定。a.节点2处施加1000N静态轴向外力u2，求应力b.求系统固有频率例1：杆件分析如图：均质；长0.5m；断面截解：a：平衡方程：A=5e-4，E=2e11，l=0.5，f2=1000，代入方程得：u1：位移，f1：节点1处应力，添加边界条件：u1=0，解得： u2=5e-10m解：由应力与应变的关系：

表示长度的变化，表示应变；b：

由刚度矩阵和质量矩阵，得特征值方程：

由应力与应变的关

式中w位固有频率，U1、U2分别是节点1、2的振幅，添加边界条件：U1=0；解得：式中w位固有频率，U1、U2分别是节点1

例2：梁的自然频率解：梁被理想化为单一单元，局部和整体的节点位移相同，如图所示：

梁的刚度矩阵：例2：梁的自然频率解：梁被理想化为单一单元，局部

质量矩阵：

节点位移向量：

与端点相关的边界条件：W1=0，W3=0；解得：质量矩阵：节点位移向量：

求解特征值：

乘以l/2EI得：

令系数矩阵的行列式等于0得：

方程的根即梁的自然频率：求解特征值：乘以l/2EI得

结果可以和精确解比较：结果可以和精确解比较：12.7——一致、集中质量矩阵

12.3节中推出的质量矩阵是一致质量矩阵，因为用于推导刚度矩阵的位移模型也用于推导质量矩阵。

一些动态问题可以用形式简单的质量矩阵求解。最简单的质量矩阵——集中质量矩阵，可以通过将质点指定到节点上。

集中质量针对平移和旋转的元素,假设在平均位置两侧的特定位移表现得像个刚体而剩余的元素不参与运动。

因此这种假设不包括元素位移之间存在的动态耦合,因此产生的元素质量矩阵是纯粹的对角矩阵。12.7——一致、集中质量矩阵12.3节中1、杆的集中质量矩阵：2、梁的集中质量矩阵：

旋转自由度的惯性影响被假定为0；若考虑惯性影响，有转动惯量：

集中质量矩阵变为：

1、杆的集中质量矩阵：

对于一般的动态问题，两者谁能得到更精确的解？

两个质量矩阵很相似，他们不考虑各种位移自由度的元素之间的动态耦合。

他们的形状函数也近似，都是用静态的位移模型推导而来。

然而,由于集中质量矩阵对角，在计算时他使用更少的存储空间。

下面的例子说明了在一个简单的振动问题中，集中和一致质量矩阵的应用。集中质量矩阵与一致质量矩阵：对于一般的动态问题，两者谁能得到更精确的解例：杆的一致和集中质量矩阵

用一致和集中质量矩阵求如图所示两端固定杆的固有频率，用两个杆单元建模。解：单元的刚度和质量矩阵分别是：

质量矩阵的下标c和l分别表示一致和

集中质量矩阵。例：杆的一致和集中质量矩阵用一致由于该杆由两个单元建模，组合的刚度和质量矩阵如下：

方框中的部分分别与单元1和2相关。

由于该杆由两个单元建模，组合的刚度和质量矩阵添加边界条件U1=U3=0后，特征值问题为：特征值w由以下方程课解：代入已知条件得：

用一致质量矩阵

用集中质量矩阵

添加边界条件U1=U3=0后，特征值问题为：

事实上，方程的精确解析解是：解得：事实上，方程的精确解析解是：12.8MATLAB应用举例例12.5阶梯轴的有限元分析图12.11中的阶梯轴满足一下条件：A1=16×10-4m2，A2=9×10-4m2，A3=4×10-4m2，Ei=20×1010Pa，i=1,2,3，pi=7.8×103Kg/m3，i=1,2,3，l1=1m，l2=0.5m，l3=0.25m。编写一个MATLAB程序解决以下问题。a，在在载荷p3=1000N下u1，u2，u3的位移b，阶梯轴的固有频率和模态12.8MATLAB应用举例例12.5阶梯轴的有限元分析12.8MATLAB应用举例解决方案：阶梯轴的刚度矩阵和质量矩阵如下所示：12.8MATLAB应用举例解决方案：阶梯轴的刚度矩阵和质在载荷p3的作用下系统的平衡方程如下所示：