有限元6-动力有限元

有限单元法第六章结构动力分析有限元法1此前述及的问题属于静力分析问题，即作用在结构上的荷载是与时间无关的静力。由此求得的位移、应力等均与时间无关。实际工程中的大部分都可简化成静力问题。当动载与静载相比不容忽略时，一般应进行动力分析。如地震作用下的房屋建筑，风荷载作用下的高层建筑等，爆炸作用下的建筑结构，船舶受海浪冲击，都应计算动荷载作用下的动力反应。研究课题中以动力问题为主。解决动力问题有两大工作要做：一是动荷载的模拟和计算，二是结构反应分析。本章将讨论如何用有限元来解决动力计算问题。26.1结构动力方程一、单元的位移、速度和加速度函数设单元的位移函数为式中：单元位移函数列阵{f}，结点位移函数列阵{d}e均是时间t的函数。由(6-1-1)可求得单元的速度、加速度函数:(6-1-1)(6-1-2)(6-1-3)3二．单元的受力分析设图示三角形单元，当它处于运动状态时，其上的荷载一般应包括：1）其他单元对其结点的作用力2）单元内部的单位体积的惯性力和阻尼力（设正比于运动速度）：惯性力：

阻尼力：3）干扰力（已知的荷载）：(6-1-4)(6-1-5)根据达朗贝尔原理，这些力构成一瞬时平衡力系，使单元处于动平衡状态。4三．结点力与结点位移、速度和加速度之间的关系用虚功原理推导，令单元结点发生任意可能的虚位移{d*}，它满足单元所定义的位移场，即虚位移场成立。作用在单元上的外力所作的外力虚功：单元内部应力在由于虚位移所引起的虚应变上所做的内力虚功：虚实5根据虚功原理（T=W），若将惯性力{Fm}，阻尼力{Fc}用上面的(6-1-4)，(6-1-5)代替，得：由于虚位移的任意性，可从等式两边各项中消去{d*}T，得：简写为：6式中：单元阻尼矩阵（第二项为阻尼力)包括由作用在单元上的干扰力转化成的等效结点荷载(6-1-6)即为单元结点力之间的关系式。(6-1-6)单刚（第一项为弹性恢复力）质量矩阵（第三项为惯性力）7四、结构的动力方程有了上述单元力关系式，象在静力问题中对每个结点建立平衡方程一样，根据达朗贝尔原理，对每个结点建立动平衡方程后，即可得到结构的动力方程组：(6-1-7)式中：[M]、[C]、[K]分别为总质量、总阻尼、总刚度矩阵。[R]为外力（结点力），实为干扰力（当不考虑静载时）当不计阻尼影响时，上式成为：(6-1-8)8若干扰力为零，得：(6-1-9)即结构的无阻尼自由振动微分方程组。由此可求得结构的自振特性（频率，振型）。由上可见，动力问题首先要解决如下问题1)[M]、[C]、[K]及[R]形成2)方程组的求解时域解法、频域解法、振型分解法96.2单元质量矩阵和单元阻尼矩阵单元质量矩阵可分为：1）协调质量矩阵[m]C（或一致质量矩阵）。单元质量矩阵表达式为：2）集中质量矩阵[m]L认为质量集中在节点，因此单元质量矩阵为对角矩阵。也有的采用上述两种质量矩阵的线性组合（h一般取0.5）。10一、协调质量矩阵（Consistent

mass）由于单元质量矩阵表达式为

式中包含了推导单元刚度矩阵时相同的形函数[N]，因此常将按此式形成的[m]称为协调质量矩阵（或一致质量矩阵）。几种常见的协调质量矩阵1.梁单元设梁单元位移函数：11式中形函数设单元的质量沿梁的长度方向均匀分布，则有：(6-2-1)W=mg单元重量，g为重力加速度12考虑轴力的梁单元协调质量矩阵式中G为梁单元的重量132.平面三结点三角形单元形函数设单元的质量均匀分布，则有单元位移函数式中：为三角形单元的重量。14二、集中质量矩阵（Lumped

Mass

）单元的协调质量矩阵和单刚具有相同的阶数，因此，总质量矩阵[M]的阶数也与总刚[K]相同。或者说采用协调质量矩阵后，结构的振动自由度和结构的静力自由度是相同的。动力问题的这种做法，其求解是很费时的：(1)形成质量矩阵的工作量等同于总刚(2)特征值的求解工程实际和试验证明，在某种干扰力作用下，结构的动力反应是有明显主次之分的。工程上通常把单元的分布质量集中到各结点而成为集中堆聚质量，这样可使问题得到很大简化，且计算经验表明，二者给出的计算精度相差无几。15如果非结构构件的质量占较大的比例，则两种质量矩阵的差别将更小质量集中按静力等效原则，且常忽略转动惯量的影响，上述各单元的质量矩阵简化为：梁单元：常应变三角单元16矩形薄板单元：由此可见，采用集中质量时，集中质量即采用“就近堆积”的原则。17对于具有中间节点的高次单元，质量集中的方式常有1）按静力等效原则根据等效节点力公式得集中质量矩阵各对角元素为与协调质量矩阵[m]C的关系为缺点：角节点可能出现负值2）按协调质量矩阵的对角线元素的比例18三、阻尼矩阵在土木结构的动力问题中，结构振动耗能是客观存在着的，根据结构振动耗能的原因，阻尼可被分为外部阻尼和内部阻尼。

外部阻尼是指结构周围介质对结构振动能量的消耗。它可分为介质阻尼、外摩擦阻尼和辐射阻尼。

内部阻尼是指结构自身对振动能量的耗散性能，可分为材料阻尼和结构阻尼。19土木结构动力分析中，需要将结构的阻尼作用抽象为既能够合理反映结构主要阻尼机制，又便于计算的简单数学模型。由于实际结构的振动阻尼是多方面因素共同形成的，其综合作用机理相当复杂，迄今尚未有既具方便数学处理，又能很好体现其物理意义的计算方法。一百多年来，人们已提出了多种阻尼模型，目前在工程结构动力计算中常用的阻尼计算模型有：粘滞阻尼模型、滞变阻尼模型、库伦阻尼模型、比例阻尼模型和非比例阻尼模型。其中用得最多的粘滞阻尼模型和比例阻尼模型。20比例阻尼在结构动力分析中，较多采用的是粘滞阻尼理论，即假定阻尼力与速度成正比，由此得到单元的阻尼矩阵：这样做，可给方程组的求解带来方便。但这个假定并不能很好的符合结构的实际情况。因此在实际应用中也常采用[M]和[K]的线性组合。或(1)(2)(3)即瑞利(Rayleigh)阻尼模型。21也可写成更一般的形式:(4)即柯西(Caughey)阻尼模型。这样处理的阻尼矩阵可以通过模态向量正交化完全解耦为对角阵。由于阻尼矩阵依赖于质量矩阵和刚度矩阵，故可通过[M]、[K]而获得[C]第i阶振型的阻尼比

xi与（4）式中各系数的关系为22当取S=0、1两项时，即为式（2）如果已知某2阶振型的频率与阻尼比，则式中系数由下式确定：ωi─第i个固有频率,ξi─阻尼比,ξ是实际阻尼与临界阻尼的比值,当ξi=ξj=ξ时,上式简化为236.3结构的自由振动和特征值问题

一、特征值方程这是一个大家都很熟悉的问题，故着重讨论程序设计上的一些处理方法。结构作无阻尼自由振动的微分方程(6-1-9)：设结构作简谐振动：将其代回(6-1-9)得：求解特征方程即可获得几个运动自由度所对应的频率和振型。(6-3-1)24另一种求特征值的频率方程：由结构力学可知，如果从质点的位移方程出发，建立运动方程（柔度法），则n阶自振方程：相应频率方程：(6-3-2)式中：[F]—柔度矩阵；[E]—n阶单位矩阵。解此方程亦得到上述相同结果。

25结构固有频率与振型的特点1）结构的固有频率都是正实数如果存在刚体位移，则固有频率为“0”2）结构的固有频率的分离性，即

w1<w2<...<wn3）振型的正交性26二、特征值方程求解方法

特征方程的解法很多，可分为1）迭代法2）雅可比法3）能量法瑞利法、里兹法3）子空间迭代法无论哪种方法，也不管是6-3-1或6-3-2，对一个大型结构系统来说，如果它的运动自由度与静力自由度一致的话，其求解是很费时的。工程实际应用中，所要求的运动自由度远小于静力自由度。27自由度的缩减如图示刚架，在水平地震力作用下，我们通常只取各楼层的水平位移为运动自由度，即对该结构来说静力自由度为72，而运动自由度仅6(不考虑轴向变形)。在运动方程6-1-7里，如果只取与运动自由度有关的位移分量，则可使求解大为简化。此时,可将n阶自由度体系的运动方程改为(6-3-3)式中下标E代表只与振动自由度有关的位移分量。28从包含运动全部自由度的运动方程6-1-7中分离出运动自由度的运动方程6-3-3，常称为凝聚或缩减。2930从包含运动全部自由度的运动方程6-1-7中分离出运动自由度的运动方程6-3-3，常称为凝聚或缩减。凝聚可采用多种办法：1.自由度编号分块法（王志楷《高等结构力学》)在形成总刚时，人为地将与运动自由度相应的位移分量排在一起，而将其余自由度另放一起，设与质点运动自由度有关的位移分量为{D2}，则总刚可写成如下分块形式：31(1)(2)因P1＝0，由（1）得(3)将（3）代入（2）得简写为式中在房屋建筑抗震问题中，[K]E又称为抗侧刚度矩阵（只含楼层的水平位移），此法缺点：带宽大，存多矩阵、多个矩阵乘、[K11]-1甚大，只限于研究性、小问题。322.根据质量矩阵[M]中主元为零的信息，在解方程中直接加以处理（特殊的解频率方程的方法，通常只用于自由振动）3．模态综合法（子结构法用于动力计算）4．直接形成与[K]E相应的柔度矩阵[F]这是较好的一种方法，在结构抗震的弹性分析中较多采用。其思路是，根据柔度中柔度系数的物理意义（单位力作用下的位移）可分别在各质点加水平单位力，逐一求出其与运动自由度对应的位移，从而获得[F]。33式中fij表示当质点j作用一单位水平力时，质点i所在结点的水平位移。由于求[F]时，总刚已经分解，将每个单位力Pi=1作为一组独立工况，可很快解出其位移列阵{D}。从中挑出与质点水平位移有关的分量即可组成[F]。计算工作量小，不需附加数组，程序简单。可由[F]直接代入频率方程求解，也可对[F]求逆得出抗侧总刚[K]E成柔度矩阵的程序段：DO30I=1,MA动力自由度数DO10J=1,NN 未知量总数10P（J）＝0 右端项P（JW（I））＝1只分解加相应单位力CALLJFC（NN，LD，2）解方程DO20J＝1，MA20A（J，I）＝P（JW（I））柔度矩阵30CONTINUE杆系层模型35层剪切模型抗侧刚度矩阵[K]E的形成方法层间刚度由各个竖向构件的抗侧刚度合成366.4地震反应分析的计算机方法-振型分解反应谱法

求解结构自振特性是结构动力分析中很重要的一步。但在工程应用中常需动力反应分析，如地震反应分析。地震反应分析通常采用两种方法：反应谱法和时程分析法（直接积分法）。我国抗震规范提供了三种计算动力反应的方法：底部剪力法、振型分解反应谱法和时程分析法。底部剪力法是一种简化的手算方法。37－、振型分解法概念振型分解法是利用振型的正交性，将原来的n阶联立方程分解成相互独立的振动方程，把结构的复杂振动，分解成按各个振型的独立振动的叠加，也可理解成把结构结点位移列矢量{D(t)}表达成：式中为上节求得的振型矩阵，称为振型坐标，或广义坐标。(1)(6-4-1)(2)(6-4-2)383940将(2)代入(1)并利用振型的正交条件：（设：阻尼亦满足正交条件）便可将式(1)化为n个独立非耦联的微分方程：(3)(6-4-3)两个不同的主振型的对应位置上的位移相乘，再乘以该质点的质量，然后将各质点所求出的上述乘积作代数和，其值等于零。41第j振型的广义刚度系数第j振型的广义阻尼系数则可化成一般形式的单质点振动力方程：(4)(6-4-4)若引入记号：第j振型的广义质量42省去上式中的下标j，并将右端项换成P(t)，便得常见单质点振动微分方程形式。它的解可用杜哈美积分表示：(6-4-5)式中为第j振型的振型参与系数

为第j振型阻尼比43二、地震反应谱概念（思路）反应谱理论创立于1940年代，它考虑了结构动力特性与地震动特性之间的动力关系。通过反应谱来计算由结构动力特性（周期、振型和阻尼）所产生的动力效应。反应谱方法基于以下假定(水平振动)（1）假定地震时建筑物的地基只作平行于地面的刚体运动；（2）假定地面运动过程可以用强震仪记录来表示；（3）假定建筑物是弹性体系。44利用仪器所记录的地面加速度时程曲线，按振型分解的方法将建筑物在地震作用下的强迫振动化为一系列单自由度体系问题。由杜哈美积分(6-4-5式)，通过数值积分方法便可计算出任何一个频率为ω，阻尼比为ξ的单质点体系在给定加速度影响下的相对位移反应曲线。45在相同地震加速度记录作用下，不同ω的单质点反应各不相同，即各自反应曲线的峰值和频率特性也各不相同。设计中最关心的往往是最大反应。阻尼比ξ一定时，对于不同的自振周期T（频率ω）都可以找出相应最大位移反应ymax。于是，对于每一个地震加速度记录，都可以算出一组以ξ为参数的ymax与周期T之间的关系曲线，这就是我们常说的位移反应谱。此外，还有加速度反应谱等。4647不同地震波输入计算出来的加速度反应谱48工程设计中，一般都采用地震作用（荷载）的概念也就是通过荷载来计算结构的内力。为了计算地震荷载，还要涉及到速度和加速度反应谱，以及这些反应谱之间的关系。最后可以计算出地震作用时在单质点体系上的最大惯性力，即地震作用力。式中：地震系数，即以重力加速度g为单位的地面运动最大加速度。动力系数，即以地面最大加速度为单位的加速度反应谱。地震影响系数49我国抗震规范提供了地震影响系数α曲线。规范规定，α与自振周期、场地类别等有关。对于多自由度体系，即可利用反应谱，查得对应于各振型的绝对最大值的地震反应，得到简化的地震作用计算公式。50水平地震影响系数αmax地震影响6度7度8度9度多遇地震0.040.08（0.12）0.16（0.24）0.32罕遇地震—0.50（0.72）0.90（1.20）1.40设计地震分组场地类别ⅠⅡⅢⅣ第一组0.250.350.450.65第二组0.300.400.550.75第三组0.350.450.650.90特征周期Tg应根据场地类别和设计地震分组按下表采用：51关于阻尼调整系数和形状参数应符合下列规定：(1)曲线下降段的衰减指数应按下式确定：

式中：γ──曲线下降段的衰减指数;ζ──阻尼比。(2)直线下降段的下降斜率调整系数应按下式确定：

式中：η1──直线下降段的下降斜率调整系数，小于0时取0。52(3)阻尼调整系数应按下式确定：式中：η2──阻尼调整系数，当小于0.55时，应取0.55。除有专门规定外，建筑结构的阻尼比应取0.05。这时。对于周期大于6.0秒的结构，地震影响系数应专门研究。设计分组参见附录《我国主要城镇抗震设防烈度、设计基本地震加速度和设计地震分组》。53三、地震作用计算对于多自由度体系，可利用反应谱，查得对应于各振型的绝对最大值的地震反应，从而得到简化的地震作用计算公式（规范公式）。式中：αj---地震影响系数；

γj---振型参与系数；1.不进行扭转耦连计算的结构规范给出第j振型i质点上的水平地震作用标准值54xji－振型向量

中的元素i（j振型i质点的水平相对位移）Gi

－质点i的重量（重力荷载代表值）

对于第j振型对结构的地震作用效应，可将当作一组静载施加于结构，计算其内力、位移（地震作用效应）。对于各振型的总作用效应：上式称为平方和开方法（SRSS法），适用于平面问题的计算，当用振型分解反应谱法对结构物作空间分析，也即考虑结构物的空间扭转效应时，采用完整二次项组合法（CQC法）552.按扭转耦连振型分解法计算1）单向水平地震的扭转效应：式中：SEK：地震作用标准值的扭转效应；

Sj、Sk：分别为j、k振型地震作用标准值的效应，可取前9～15个振型；ζj、ζk分别为j、k振型的阻尼比；

ρjk第j振型与k振型的耦联系数；λT：k振型与j振型自振周期比。56按扭转耦联振型分解法计算，可取两个正交水平位移和一个转角，共三个自由度计算。j振型i层水平地震作用标准值式中：Fxji、Fyji、Ftji—分别为j振型i层的x方向、y方向和转角方向的地震作用标准值；Xji、Yji—分别为j振型i层质心在x、y方向的水平相对位移；2）双向水平地震的扭转效应：57Φji—j振型i层的相对扭转角；ri—i层转动半径，可取i层绕质心的转动惯量除以该层质量的商的正二次方根；γtj—计入扭转的j振型的参与系数。586.5时程分析法

（直接积分法、步步积分、逐步积分）反应谱理论尽管考虑了结构的动力特性，然而它仍然是把地震惯性力当作静力来对待，所求的是一个综合的最大值，所以还只是一种准动力方法。时程分析法是采用地震加速度时程曲线作为输入参数，对结构进行地震反应的时间历程分析，是一种完整的动力分析方法。可以用于:（1）地震反应分析（2）拟动力试验59优点：（1）全面考虑了强震中震动的振幅、频谱和持时三要素。（2）全面考虑了长周期分量对高层建筑的不利影响。（3）可了解到地震过程中结构物的屈服机制和较准确地找出结构的薄弱部位。6061塑性铰发展过程62(子结构)拟动力试验计算子结构试验子结构63基本思想：将本来应在任何时刻都应满足的运动方程的位移矢量{D}简化为只要在时间离散点上满足方程。而在一个时间间隔内，对位移、速度、加速度的关系则采用某种假定。依所取假定不同而有各种不同的积分方法。如假定在一个时间间隔内加速度按线性变化的线加速度法，以线加速度法为基础的wilson-θ法，以及Newmark-β法等。64一、线加速度法的基本原理

因为wilson-θ法是以线加速度法为基础的，故首先介绍其基本原理。在时刻ti到ti+1（ti+1=ti+Δt）之间，假定ti时刻的质点的地震响应()是已知的，ti+1时刻的质点的地震响应（）是待求量。线加速度法假定在每个时间间隔Δt内的加速度呈线性变化。因此，位移对时间的三阶导数应为常数，即：(a)在时段Δt内，刚度、阻尼都不变化65对τ求导得：当τ=Δt时，即时程由ti变化到ti+1时，由式（c）得：(d)或写成：(e)将τ时刻的位移反应yτ（即前面的D）在时刻ti开始时展开成泰勒级数：(b)将式（a）代入式（b）得：(c)66由式（e）得：(g)将式（g）代入式（f）可得：(h)由式（d）得：或写成：(f)加速度增量和位移增量的关系速度增量和位移增量的关系67在t+Δt时，式（i）同样成立：(j)将速度增量的表达式（h）和加速度增量表达式（g）代入上式并整理得：(k)(l)从式（g）和式（h）可看出，速度增量、加速度增量都可以通过位移增量Δy算出，为了求Δy，我们先讨论单自由度的情况，然后再扩展到多自由度。单自由度体系的运动方程为：(i)68式（l）叫增量拟静力平衡方程，据此可算出Δy（因为步长Δt及、等均为已知），再由式（