《第十章-动力问题的有限元法(共34页).doc》由会员分享，可在线阅读，更多相关《第十章-动力问题的有限元法(共34页).doc（34页珍藏版）》请在taowenge.com淘文阁网|工程机械CAD图纸|机械工程制图|CAD装配图下载|SolidWorks_CaTia_CAD_UG_PROE_设计图分享下载上搜索。
1、精选优质文档-倾情为你奉上第十章 动力问题的有限元法动力问题：结构受载荷作用没有达到平衡状态，或由于结构的弹性和惯性而围绕平衡位置振动时，其位移、应力、应变都是时间的函数，除各点有位移外，还有速度和加速度。10-1 结构的动力学方程在动态情况下，结构承受的载荷可以随时间而变化，是时间的函数。按有限元方法，将此载荷以作功等效的原则，移置到结点上，得到的结点载荷列阵也应该是时间的函数，用表示。此外，结构在运动中，各点除位移外，尚有速度和加速度，分别用、和表示。显然，、和不仅是坐标的函数，还是时间的函数，即 采用有限元法对结构进行离散，对任一单元，若单元的结点位移为、单元的结点速度为、单元的结点加速
2、度向量为，则单元内的位移可由结点位移插值得到，于是 （10-1） 按照达朗贝尔原理：一个运动中的物体，若加上惯性力，则可看为静力平衡状态。设材料密度为，则结构内单位体积的惯性力为，这对结构来说，相当于受到另外一种体力分布载荷，大小与质点的加速度成正比，而方向与加速度方向相反。另外，在介质中运动的质点总是会受到阻力作用，如阻力与速度成正比，则单位体积的阻力为，表达式中为阻力系数，负号表明阻力方向与速度相反。其次，设结构还受到体力和面积力作用，则结构所受总的体力为： 根据达朗贝尔原理，结构受到所有力（包括惯性力）作用处于静力平衡，于是 （10-2）上式等式右边的第一积分相当于体力的等效结点力。将惯
3、性力和阻力项代入上式，得 （10-3）上式右边第一、第二项即为一般静力分析时的等效结点力，令 （10-4） （10-5） （10-6）、分别为惯性力和阻力对应的等效结点力，考虑到积分与结点位移无关，于是 （10-7） （10-8）式中 （10-9）称为单元质量矩阵，而 （10-10）称为单位阻尼矩阵。于是，单元动力方程化为 （10-11）若结构系统有个单元，则对单元动力方程扩大后相加，得 （10-12）式中称为整个结构的质量矩阵，而称为整个结构的阻尼矩阵。2质量矩阵式（10-9）所表示的单元质量矩阵称为一致质量矩阵，因为形成质量矩阵所用的形函数与位移插值函数所用的形函数相同，单元的动能和位能是
4、相互协调的。此外，还可采用一种较为简单的质量矩阵，它假定单元的质量团聚在单元的结点处，因此称为团聚质量矩阵。单元的一致质量矩阵是一满阵，而团聚质量矩阵为一对角阵。以常应变三角形单元为例，一致质量矩阵团聚质量矩阵为 3阻尼矩阵式（10-10）所示的单元的阻尼矩阵还可进一步表示为由此可见，当阻力正比于速度时，阻尼矩阵正比于单元的质量矩阵，这种阻尼称为粘滞阻尼。此外，还可假定阻尼力正比于应变速度，此时可见此时阻尼矩阵正比于单元的刚度矩阵。由材料内摩檫引起的阻尼通常可假设为这种情况。但结构的阻尼往往是多重因素造成的，为了能更一般地描述结构的阻尼，有关文献建议将单元阻尼表示为式中为常系数，实际计算中，通
5、常取前两项，于是得 （10-13）即认为阻尼矩阵是质量矩阵和刚度矩阵的线性组合，这种阻尼通常称为比例阻尼。10-2 无阻尼自由振动为了计算动力方程，首先讨论无阻尼自由振动，即考虑 的情况于是，动力学方程（10-12）变为 （10-14）这便是无阻尼自由振动方程，该方程组具有下列形式的解： （10-15）将上式代入（10-14）式，得 （10-16）记 可将式（10-16）改写为 （10-17）求解方程组（10-17）的问题称为广义特征值问题，满足方程（10-17）的解及其相应的向量分别称为特征值和特征向量。由此求得的就是结构的固有频率，而为相应的固有振型。我们已经知道，当结构上的约束足以消除刚
6、体位移时，其总体刚度矩阵是对称正定的。由于材料的质量密度，不难证明，一致质量矩阵是正定的，而团聚质量矩阵是正定或半正定的对角阵，例如，对梁、板单元，当考虑转动自由度而忽略转动惯量时，团聚质量矩阵因与转动自由度相对应的对角元为零而成为半正定（不为零的元素全为正值）。自由振动时，结构中各点的振幅不可能全为零，因而齐次方程组必有其系数行列式等于零，即它称为广义特征方程，如果矩阵的阶数为，则广义特征方程是的次代数方程，由此可解出个广义特征值，从而可得到个固有频率。对每个固有频率，由式（10-17）可确定维列向量，它代表个质点的振幅构成的振型。个特征值和相应的个特征响量常称为个特征对。求标量和非零向量，
7、使满足方程组的问题，称为标准特征值问题，方程（10-19）的解和称为矩阵的特征值和特征向量。广义特征值问题可以化为标准特征值问题，如果质量矩阵是正定的，可以将进行三角分解其中为下三角阵，方程（10-17）可写为 用前乘上式，并注意到，则上式变为 若记 ， （10-20）就得到了式（9-17）所对应的标准特征值问题。如果是正定的，则也是对称正定的，且当为对角阵时，是一个与相同的带状矩阵，在求得标准特征值问题的特征值和特征向量后，即可求得固有频率，并利用（10-20）求出相应的振型。 （10-21） （10-22）在总体刚度矩阵是正定的情况下，可以将进行三角分解采用同样的方法，也可将（10-15）
8、化为如下的标准特征值问题 （10-23）其中 若是对称正定的，则也是对称正定的，但它一般不是带状的。在求得特征值和特征向量后，由式（10-24）可求得固有频率和和相应的振型。 （10-25） （10-26）目前标准特征值问题已有充分的研究，有许多现成的算法和程序可供使用。10-3 特征值和特征向量的性质在研究特征值问题求解之前，首先介绍特征值和特征向量的性质。1 实对称阵的特征值问题 可以证明：对于和均为实对称阵的广义特征值问题，其所有的特征值都为实数，所有的特征向量也为实向量。 事实上，阶对称阵和的元素都是实数，且有 ，（）对广义特征值问题，应有即 用的共轭数乘上式的两端，并对求和，得 （1
9、0-28）由于为对称阵，因此即式（10-28）的左端是实数，同理，式（10-28）的右端也是实数。因此，由式（10-28）可知必为实数，从而由（10-27）解得的特征向量也比为实向量。显然，若为特征向量，则也必为特征向量。2 特征向量的正则化由上述可知，若为特征向量，则也必为特征向量，即相差一常数因子的一系列平行线均为特征向量，为确定起见，我们规定： （10-29）由此确定的特征向量称为正则化特征向量。式（10-29）就为正则化条件。3 特征向量的正交性设和为广义特征值问题的两个特征对，则有分别用和前乘上两式，得 （） （）（）式转置减（），考虑到和均为实对称阵，得若，且，则 （10-30）这
10、就是特征向量关于的正交性，将式（10-29）和（10-30）合并，可写成 （10-31）其中为kronecker符号。将式（10-31）代入（）得 这说明特征向量对矩阵也具有正交性。假设是由个特征向量组成的特征向量矩阵，该矩阵称为模态矩阵，即 （10-33）设是由特征值组成的对角阵，即 （10-34）则特征向量的正交性可表示为4 雷利（Rayleigh）商关于特征值的性质还可以从雷利商出发去研究它，对于广义特征值问题，雷利商定义为 （10-36）式中为任一向量，设它可以由特征向量的线性组合表示，即 （10-37）式中是由个系数组成的向量。将（10-37）代入（10-36），注意到特征向量的正交
11、性，得展开后可以写成 （10-38）将特征值由小到大排列：可见，雷利商在最小特征值和最大特征值之间，即并且，当时这说明，当取第阶特征向量时，雷利商达到它的一个极值，该极值就是对应的特征值，特别是这就是雷利商的极小值原理，次原理常被工程界作为预估振型的依据，用以计算特征值，从而求得基频的近似值。如果，即向量与前阶特征向量正交，则由式（10-38）看出：且当时这就是说，是雷利商在向量与前阶特征向量正交条件的极小值。5 移位 对于广义特征值问题，可以使作一移位10-4 逆迭代法求解代数特征值问题的方法很多，这里介绍迭代法。对于阶广义特征值问题当和均为正定阵时，我们已经证明了它有个特征对。由个特征向量
12、可以构成模态矩阵同时，为了说明迭代过程，令， （10-39）注意式（10-39）与式（10-34）中的区别，此时广义特征方程变为用模态矩阵表示为 （10-40）式（10-40）两边右乘任意向量得 （10-41）若已知，将看为某一未知量，则由代入试式（10-41），求得 （10-42）即，经过一次迭代后，对应的特征向量的分量扩大倍。如此经过次迭代后得 （10-43）由于，当迭代次数较大时，有 （10-44）而前后两次向量的模之比的平方为 （10-45）由此可见，迭代格式使向量逐渐与平行，前后两次向量的模之比趋于，向量的收敛速率为，特征值的收敛速率为。实际上，这样的迭代是无法进行的，因为和正是需要
13、求解的。为了克服这一矛盾，并考虑到初向量的任意性，只要把作为初向量，把作为迭代后的向量，上述迭代过程就可进行了。此时 （10-46）实际进行迭代时还需要解决一个问题，从（10-43）可以看到，若，则的模趋于无穷大；若，则的模趋于零。不加处理在计算上是无法实施上述迭代的。为解决这一问题，只要在迭代过程中不断使初向量的模取为1即可。逆迭代法的实施步骤概述如下：（1） 选取初向量，计算；（2） 由方程 解出；（3） 计算 （4） 计算近似特征值注意：此处，若满足精度要求则转向第五步，计算特征向量，否则计算下一次迭代的初向量：上式表明，不断使迭代向量取正则化，然后转向第（2）步进行下一次迭代。（5）
14、计算特征向量取特征值为结束迭代。这是基本迭代格式，当时，收敛非常快。这种迭代过程称为逆迭代法。如果迭代方向反过来，即 也可求得最大特征值，这种迭代称为幂法。逆迭代法除求最小特征值和特征向量外，还可以扩大运用范围。和克莱姆-史密特（cram-schmidt）正交化过程相结合，可以用来求取最低的几阶特征对；也可求得重特征值对应的特征向量；和移频法相结合可以对奇异的刚度矩阵进行迭代，可以改善某一近似特征值的精度和求取对应的特征向量。现逐一叙述如下：由式（10-43）可见，假如，而，那么迭代结果应该得到第阶特征对。当求得前阶特征向量以后，要做到这一点很容易，只要选择初向量与前阶特征向量正交。虽然任意初
15、向量不一定与正交，但可以构造一个新的初向量： （10-47）与前阶特征向量正交，只要取 （10-48）这就是克莱姆-史密特（cram-schmidt）正交化过程。从迭代公式看，似乎这种正交化过程可以一劳永逸，一旦初向量与前阶特征向量正交，则所有的迭代向量都与前阶特征向量正交。但实际上并非如此，由于计算机误差的原因，不可避免产生低阶特征向量，迭代到最后还是得到最低阶的特征向量。为了克服这一点，必须每迭代一次都进行正交化处理，不断把前阶特征向量从迭代向量中清除掉。对于重特征值，用带正交化过程的迭代，可得出重特征值对应的一组正交特征向量。用逆迭代法求得第一特征对后，利用正交化过程的迭代法，可以依次求
16、取各阶特征对。为保证计算精度，前面的特征向量必须计算精确些。逆迭代一般用来计算少数前几阶特征对。移频逆迭代可以用来改善某几阶近似特征值的精度和求取相应的特征向量。设已经知道第阶近似特征值，则对广义特征值问题作以下移频 （10-49）得到新的特征值问题 （10-50）式中新特征值问题与原特征值问题的特征向量相同，各阶特征值移动一个距离。新特征值问题绝对值最小的特征值是，与它邻近的两个特征值是和。由此可知，对（10-50）进行逆迭代，将得到和，从而得到。其收敛率为和中较大的一个，即和中较大的一个，所以只要比较接近收敛便很快。最后还要提一下初向量怎么选取的问题。初向量选取得好坏，对迭代效率影响较大，
17、显然。要求初向量选得应尽量接近所求的振型。一个简单的选取办法是，把结构简单地考虑为各个自由度是互不耦合的，认为特征方程便为非耦合形式所以，将从小到大排序，得各阶近似特征向量。设排在第位，记，。一般取第一阶初向量为 取第阶初向量为10-5 振型叠加法在10-5中我们已经导出了系统的动力学方程（10-12）在一定的初始条件下求解上述方程，得出任一时刻的位移、速度和加速度响量，这就是结构动力响应的计算。数学上方程（10-51）是一个二阶线性微分方程，原则上可用解常系数微分方程组的标准方法求界。但是如果矩阵的阶数很高，求解将要付出高的代价，而有限元分析中，有的方法是行之有效的。这里我们将介绍两种动力响
18、应的计算方法：振型叠加法和逐步积分法。本节先介绍振型叠加法。振型叠加法的基本思想是，利用系统无阻尼自由振动的振型将动力学方程转换为一组互不耦合的微分方程，然后分别求解这些方程，并将结果叠加而得到方程（10-51）的解。首先解出方程（10-51）的所对应的无阻尼自由振动的固有频率和振型。设结构被激发起个振型，其相应的特征对为。其次，将位移矢量看成是各振型的线性组合，即引入变换 （10-52）式中为振型矩阵，将此变换方程代入（10-51），两端前乘，由的正交性，得 （10-53）式中 当阻尼取比例阻尼时，即则为方便起见，设 （10-54）此时，方程（10-53）就变为互不耦合的二阶常微分方程，即
19、（10-56）此时，初始条件化为 其分量形式为由式（10-56）及其初始条件，利用杜哈美尔(Duhamel)积分得式中，而可由初始条件确定。于是方程（10-51）的解为 （10-59）关于比例阻尼系数的确定：待定系数的一般由实验确定，设和分别为第和第个固有频率，和分别为第和第个振型对应的阻尼比，即实际阻尼与临界阻尼的比值，可由实验测定。根据式（10-54）得 由此解得 （10-60a） （10-60b）上式中一般取第一和第二阶特征值和相应的阻尼比进行计算。在求得后，其它阻尼比可由下式求得 由上式可以看出，比例阻尼中高阶阻尼很高，这实际上就等于在计算中将高阶振型的响应消除了。事实上，高阶振型的响
20、应也是很小的，忽略它是合理的。所以振型叠加法不必考虑全部振型，只需取少数最低几阶振型即可。 振型叠加法适用于象地震等只激发起较少振型，所需计算时间较长这类问题。例：考虑一个两个自由度系统，其运动方程为 （a）初始条件是 当，试用振型叠加法求系统的位移响应。解：该系统对应的广义特征值问题是 （b）方程存在非零解的条件为系数行列式等于零，即 （c）上式是关于特征值的特征方程，令，则特征方程化为 (d)由上式求得 (e) (f)将(e)代入（b）得显然方程的系数行列式为零，设代入，则得，于是为系数因子，无论为何值，方程均满足，为了唯一确定特征响量，根据正则化条件：于是，得 (g)将(f)代入（b），
21、并由特征向量的正则化条件，同理可得 (h)利用(g)和(h)将原问题转换为互不耦合的运动方程 (i) (j)初始条件化为 应用杜哈美尔积分式，可求得方程(i)和(j)的解为最后将两个振型叠加，10-6 逐步积分法逐步积分法根据动力学方程（10-51），引进某些假设：任一时刻动力学方程都成立；于是建立由时刻结构状态响量到时刻状态响量的递推关系，从而从时刻的初始状态向量出发，一步一步第求出各时刻的状态向量。在时刻之间的位移、速度和加速度采用某种假设，根据假设的不同，可以有各种个样的方法。我们这里介绍纽马克(Newmark)法和Wilson-法。1 Newmark法 建立逐步积分法的关键是建立由时刻
22、的状态向量的递推关系。时刻有三组未知量，它们满足动力学方程（10-51），即 （10-62）方程有三个未知量，显然还应补充二组方程才能求解。这两个补充方程由速度和位移的泰勒级数展开采用某种近似得到。Newmark法取速度的一次展开式式中是在区间中某点的值。Newmark法取近似假设 于是有 （10-63）位移的二次展开式并对取类似的假设 于是又有 （10-64）式（10-62）、（10-63）、（10-64）便是Newmark法的基本公式，由此可以建立从时刻到时刻状态向量的递推关系。由（10-64）可得 （10-65）将（10-63）和（10-65）代入（10-62）得到关于的线性方程组，从而
23、可以解出，然后将其代入（10-64）求出加速度向量；将代入（10-63）求出。参数的选择对算法的影响很大，算法稳定性分析指出，当时，Newmark法是无条件稳定的，这时可只根据精度的要求选择步长。取时，Newmark法即为平均加速度法。根据以上分析，综合计算步骤如下：1 初步计算（1） 形成刚度矩阵，质量矩阵和阻尼矩阵；（2） 获得初始状态向量；（3） 选择步长，以及参数，（），并计算下列常数， ，（4） 形成“刚度”矩阵 （5） 分解矩阵 2 对每个时间步计算（1） 计算时刻的“载荷”向量 （2）求时刻的位移（2） 计算时刻的加速度、速度 二、Wilson-法Wilson-法是线性加速度法的
24、推广，它把线性加速度假定的范围从扩充到。Wilson-法作如下的假设 （10-66）在此假设下，积分便有速度和位移的表达式 （10-67） （10-68）在式（10-67）和（10-68）中，令时刻的速度和位移 （10-69） （10-70）而时刻的动力学方程 （10-71）就可以建立从到的递推关系。但这不是我们的目的，因为这就是普通的线性加速度法，无非把步长扩大为，我们的目的是建立从到的递推关系。 从（10-69）和（10-70）写出加速度、速度用位移的表达式 （10-72） （10-73）将（10-72）和（10-73）代入（10-71），得的方程：解之得。然后求时刻的状态向量。令式（10
25、-66）中，得将（10-72）代入上式，得 （10-74）由式（10-69）及（10-70），令，得 （10-75） （10-76）式（10-74）、（10-75）和（10-76）就是Wilson-法的递推关系式。Wilson-法的特点：）首先由线性加速度法时刻求时刻的；）由求得；）由时刻加速度及时刻的加速度插值得到；）再由时刻到时刻间的线性加速度关系式的积分式求和，从而可以进行下一步递推。由算法的稳定分析指出，Wilson-法当时是无条件稳定的，计算表明一般取为宜。计算还表明Wilson-法有较高的“算法阻尼”，即由于算法引进的所谓阻尼，对于结构的高阶振型算法阻尼尤高，这就取消了高阶振型对响应的贡献。由于有限元离散必然引起较大的高阶振型的误差，考虑高阶振型没有实际意义，通过积分过程的算法阻尼，把高阶振型的响应消除掉是合理的。不难看出，当，Wilson-法就是普通的线性加速度法。现将Wilson-法的计算步骤综合如下：1 初始计算（1）形成刚度矩阵，质量矩阵和阻尼矩阵；（2）获得初始状态向量；（3）选择计算步长，并计算积分常数：， ，（4）形成“刚度”矩阵 （5）分解矩阵 2对每个时间步计算（1）计算时刻的“载荷”向量 （2）求时刻的位移（3） 计算时刻的加速度、速度和位移 专心-专注-专业