1、第十二章 动力学问题 12.1 动力学问题的有限元方程 12.2 质量矩阵与阻尼矩阵 12.3 直接积分法 12.4 特征值问题及解法 12.5 振型叠加法 12.6 减缩系统自由度的方法 1 12.动力学问题 关键概念 一致质量矩阵 团聚质量矩阵 振型阻尼矩阵 Rayleigh阻尼 显式积分 隐式积分 Guyan 减缩法 动力子结构法 有限元法基础 2 12.动力学问题 12.1 12.1 动力学问题的有限元方程动力学问题的有限元方程 (一)动力学问题的基本方程 平衡方程 几何方程 本构关系 边界条件 初始条件 有限元法基础 3,0V1()V2 V S Sij jiiiiji jj iiji
2、jklkliiuijjifuuuuDuunT?在 中在 中在 中在上在上000(,0)(,)(,0)(,)(,0)(,)iiiiiiu x y zux y zu x y zux y zu x y zux y z?12.动力学问题(二)Galerkin法?平衡方程和力的边界条件的等效积分形式 第一项分部积分 有限元法基础 4?ijijklkliiiiiiVVVSDdVuuu dVu f dVuTdS?,0iij jiiiijjiVSufuu dVunT dS?,iij jiijjijijVVVudVundV?12.动力学问题(三)有限元离散(三)有限元离散?在动力学分析时,物理量是空间(x,y,
3、z)的函数,也是时间(t)的函数,是一个四维问题)的函数,是一个四维问题?有限元离散,或网格剖分是对空间域进行,这一步骤与静力学问题分析时相同步骤与静力学问题分析时相同?时间维的离散使用有限差分法处理时间维的离散使用有限差分法处理 有限元法基础 5 12.动力学问题(四)位移插值函数 只对空间域进行离散,插值函数表示为 写成矩阵形式 有限元法基础 6 111(,)(,)()(,)(,)()(,)(,)()niiiniiiniiiu x y z tN x y z u tv x y z tN x y z v tw x y z tN x y z w t?插值函数 与时间无关 e?uNq1(,)()(
4、,)()(,)()ieiiniu x y z tu tv x y z tv tw x y z tw t?qu=qqq123 3niiN?NNNNNI12.动力学问题(五)有限元方程 将插值函数代入Galerkin 积分表达式,由 的任意性得,系统的求解方程 其中 有限元法基础 7?q()()()()tttt?MqCqKqQeeeeeeee?M=MC=CK=KQ=QeeeeeeTeTVVeTeTTVVSdVdVdVdV+dV?MN NCN NKB DBQN fN T12.动力学问题(六)典型的动力学问题?模态分析(Modal Analysis)确定结构的动力学特征?瞬态分析(Transient
5、Analysis)使用直接积分法或模态叠加法得到结构的瞬态响应?谐分析(Harmonic Analysis)线性结构承受简谐载荷的稳态响应?谱分析(Spectrum Analysis)在响应谱作用下,结构的响应 有限元法基础 8 12.动力学问题 12.2 质量矩阵和阻尼矩阵质量矩阵和阻尼矩阵?动力问题的质量矩阵 它与所使用的有限元列式的原理和位移插值函数保持一致。假定质量集中在节点上,导出的质量 矩阵是对角线矩阵,可提高计算效率。有限元法基础 9 eeTVdV?MN N一致质量矩阵 Consistent Mass 团聚质量矩阵Lumped Mass 12.动力学问题?团聚质量矩阵的计算方法(
6、1)中每一行主元等于 中该行所有元素之和 (2)中每一行主元等于 中该行主元乘以缩放 因子 a 根据平动DOF质量守恒确定,即 有限元法基础 10 elMeM?10eneeikklijijij?MMelMeM?0eeiilijaijij?MM?eeiiViadV?M()i与平动相关的行12.动力学问题?振型阻尼矩阵振型阻尼矩阵 阻尼正比于质点速度 阻尼正比于应变速度 这种阻尼称为比例阻尼或振型阻尼,比例系数与固有频率相关。和 与频率无关,为常数。有限元法基础 11 eeTVdV?CN N阻尼矩阵与质量矩阵或刚度矩阵成比例 eeTVdV?CB DBRayleigh阻尼?CMK?12.动力学问题
7、12.3 12.3 直接积分法直接积分法 半离散的动力学方程的解法分为两类,一是直接进行数值积分,一类是使用固有振型表达动态响应,称为振型叠加法。直接时间积分一般采用差分格式,分为显式时间和隐式时间积分。显式积分式条件稳定的,隐式积分是无条件稳定的,各有优缺点。有限元法基础 12 12.动力学问题 12.3.1 12.3.1 中心差分法 有限差分法的理论依据很简单,以有限增量的比值代替数学上的微分,速度表示为 中心差分格式为 有限元法基础 13 duuudtt?21()21(2)tttttttttuuutuuuut?12.动力学问题?将中心差分格式应用到有限元的半离散方程?整理得递推公式 有限
8、元法基础 14 2221111122ttttttttttt?MC qQKM qM-C q?21122ttttttttttttt?MqqqCqqKqQ12.动力学问题?中心差分法求解运动方程的步骤 1.初始计算 1)形成刚度矩阵K、质量矩阵M和阻尼矩阵C 2)给定 ,和 3)选择时间步长 ,4)计算 5)形成有效质量矩阵 6)三角分解 有限元法基础 15 0q0q0qt?crtt?20002ttt?qqqq211?2tt?MMC?T?MLDL12.动力学问题 2.对每一时间步长 1)计算时间 t 的有效载荷 2)求解时间 的位移 3)如果需要计算时间 t 的加速度和速度 有限元法基础 16(0,
9、2,)ttt?22211?2tttttttt?QQKM qMC qtt?Tttt?LDL qQ21()21(2)ttttttttttttt?qqqqqqq12.动力学问题?特点(1)若已知 和 可直接预测下一步的 ,称为逐步积分法。如果质量矩阵M是对角的,C也是对角或可以忽略,则利用递推公式求解时不需求解方程,直接可得下一时间步的预测值。有限元法基础 17 tt?qtqtt?q显示时间积分(Explicit Time Integral)12.动力学问题(2)当t=0时,需要 和 ,因此必须用专门的起步方法。由速度和加速度的中心差分公式,消去 的量,得 初始加速度可用运动方程求得 有限元法基础
10、18 0qt?q20002ttt?qqqqtt?10000()?qMQ-Cq-Kq12.动力学问题(3)中心差分是条件稳定的,时间步长不能任意取,最大步长与计算的问题相关,以及网格剖分相关。一般步长可取为 为系统的最高阶固有频率,Tn是系统的最小固有振动周期。实际应用中可以用系统中最小尺度单元的最小振动周期代替系统的Tn,因为 。有限元法基础 19 2ncrnTtt?n?minennTT?12.动力学问题(4)时间步长的确定方式 a)网格剖分后,找出尺寸最小的单元,形成单元的特征方程 求出最大特征根 ,得到 。b)网格剖分后,找出尺寸最小的单元的最小边长 L,可以近似地估计 ,由此,得 ,称为
11、Couran,Friedrich 和Lewy条件。有限元法基础 20()2()0ee?KMn?2/nnT?/nTL C?(/)CE?/crtL C?物理解释:时间步长应足够小,以致于在单个时间步内,传播不会超过相邻的两个节点间的距离。12.动力学问题(5)中心差分的显示算法,适合于由冲击、碰撞、爆炸类型的载荷引起的波传播问题的求解。因为这些问题本身就是在初始扰动后,按一定的波速C逐步在介质中传播。对于结构动力学问题,采用显示时间积分不太合适。因为结构的动力响应中低频成分起主要作用,允许大的时间步长。有限元法基础 21 12.动力学问题?例:波的传播例:波的传播 均匀钢杆,无阻尼,开始静止,突然
12、施加轴向端点力。用40个2节点杆单元模拟,材料为线弹性。图中Cn为Courant 数,即实际步长与临界步长的比值。有限元法基础 22 12.动力学问题 有限元法基础 23 12.动力学问题 有限元法基础 24 初始速度为零,开始后在加载。12.动力学问题 12.3.2 12.3.2 Newmark法法 Newmark 积分法假设,在 的时间区域内,有 其中,和 是按积分精度、稳定性和算法阻尼要求决定的参数,取不同的值代表不同的积分方案。有限元法基础 25 ttt?2(1)1()2tttttttttttttttt?qqqqqqqqq?12.动力学问题?几个特例几个特例 1),对应于线性加速度法,
13、即在时间步加速,对应于线性加速度法,即在时间步加速度内线性变化度内线性变化 2),对应于平均加速度法,即在时间步内加,对应于平均加速度法,即在时间步内加速度取平均值速度取平均值 有限元法基础 26 11,62?11,42?()/(0)ttttttt?qqqq1()(0)2ttttt?qqq12.动力学问题?Newmark 法的运动方程?由Newmark 关系式,得?递推公式为 有限元法基础 27 tttttttt?MqCqKqQ2111()12ttttttttt?qqqqq22111112 112ttttttttttttttt?KMqQMqqqCqqtq12.动力学问题?Newmark 法的计
14、算步骤 1.初始计算(1)形成刚度矩阵K,质量矩阵M和阻尼矩阵C(2)给定 ,和(3)选择时间步长 ,以及参数 、和积分常数 (4)形成有效刚度矩阵(5)三角分解 有限元法基础 28 0q0q0qt?012324567111121212ccccttttccctct?01?cK=K+c M+C?TK=LDL12.动力学问题 2.对每一时间步长(1)计算时间 的有效载荷 (2)求解时间 的位移 (3)计算时间 的加速度和速度 有限元法基础 29 tt?Ttttt?LDL qQ(0,2,)ttt?023145?ttttttttttcccccc?QQMqqqCqqqtt?tt?02367()ttttt
15、ttttttttccccc?qqqqqqqqq12.动力学问题?Newmark 法的特点 (1)为隐式积分算法(Implicit Time Integral)每一步都必须求解方程;(2)当 时算法是无条件稳定的,即时间步长得大小不影响解得稳定性;(3)当 时是条件稳定的,;(4)Newmark 法特别适合于时程较长的系统数瞬态 响应分析,而且大时间步长可以滤掉高阶不精确 模态对系统响应的影响。有限元法基础 30 20.5,0.25(0.5)?0.5,0.5?maxcritt?12.动力学问题 有限元法基础 31 12.动力学问题 有限元法基础 32 12.动力学问题 12.4 12.4 特征值
16、问题及其解法 系统的运动方程为 无阻尼自由振动退化为 设方程解的形式为 方程成为 有限元法基础 33()()()()tttt?MqCqKqQ()()0tt?MqKq?()itte?qq2?()0ite?MK q2?()0?KM q广义特征值问题 2?特征根特征向量q 12.动力学问题 四种类型的解法:四种类型的解法:?直接矢量迭代法(幂法)?矩阵变换法?多项式迭代求解法(行列式搜索法)?利用特征多项式的Sturm序列特性求解法 以及 iii?K=M()I()TT?K=M=对角阵单位阵()det()0p?KM()det()p?KM()()()()()()det(),1,2,1rrrrriprn?
17、KM 12.动力学问题 12.4.1逆迭代法(幂法)逆迭代法(幂法)对方程 取近似解 按以下迭代格式求解 则序列 将收敛于相应的特征根 的特征矢量。?K=M1x1(1,2,3)kkk?KxMx12,x x1?12.动力学问题?因为对任一矢量可用特征矢量表示为 代入方程 按迭代方程有 若 ,当 时,110?x Mk?11k?x11 niiix?1111 nniiiiiiix?MMK12111 niiiixx?K M1111111 kkknnkiiiiiiiix?12n?设 12.动力学问题?为了使Xi不受计算的影响,常常需要归一化?正迭代法的计算方案 迭代格式 若 ,当 时,?特征根 的近似解
18、1111/21()kkkTk?xxx M x11111/21()kkkkkTk?MxKxxxxM x10Tn?x Mk?1kn?x11111()kkkTkk?xKxxxM x 12.动力学问题 12.4.2变换法?广义特征值问题化为标准特征值问题 有限元法中的质量矩阵M是对称正定的,则 故有 定义 得到 T?MSST?KMSS?1T?S KST?S?K 12.动力学问题?标准特征值问题标准特征值问题 变换法中有Jacobi 法、法、Givens 法、法、Householder,其实质就,其实质就是通过一系列的变换矩阵,将是通过一系列的变换矩阵,将M变换成单位矩阵,将变换成单位矩阵,将K变换成对
19、变换成对角矩阵。?Jacobi 法法 标准特征值问题的方程 设完成第k步变换成为步变换成为 Pk 是正交矩阵,即,即?K1Tkkkk?KP K PTkk?P PI 12.动力学问题?Pk矩阵的构造矩阵的构造 12.动力学问题?特点?在 时,矩阵K趋于对角阵?由于只能做有限次变换,因此最后的矩阵是对角占优?变换后的矩阵总是对称的,可以减少计算次数?在一次变换使非对角线为零元素,在下次变换中可能成为非零,因此收敛缓慢?需要结合一些其他策略提高计算效率 k?12.动力学问题 12.4.3 子空间迭代法 子空间迭代法是求解大型特征值问题的低阶特征值有效方法,它实际上是Rayleigh-Ritz 法和同
20、时逆迭代法的组合。?子空间迭代法的步骤 1)建立q个初始矢量 (qp,p是要计算的特征根个数,一般q=min(2p,p+8)2)从q个迭代矢量中使用逆迭代法和Ritz分析抽取近似的特征根和特征矢量 3)迭代收敛后,使用Sturm 序列检查验证所得特征根和特征矢量是否符合要求 12.动力学问题?子空间迭代法求解过程 q个初始迭代矢量构成nq阶矩阵 X1 第k步迭代为 形成子空间投影矩阵 求解子空间特征系统 这是Rayleigh Ritz分析,Kk+1 是qq 计算近似特征矢量 Xk+1可作为新的迭代矩阵,当 时,11,2,3,kkk?KX=MX111111,TTkkkkkk?KXKXM=XMX1
21、1111kkkkk?KAMA111Tkkk?X=XAk?11,kk?A 12.动力学问题 12.4.4 Lanczos法 Lanczos方法目前被认为是求解大型矩阵特征值问题的最有效方法,与子空间迭代法相比,其计算量要少得多。?Lanczos变换 选取初始矢量x,并计算 11/2()T?xxx Mx 12.动力学问题?理论上讲,xi(i=1,2,n)是关于M正交的,即?定义矩阵 满足关系 12.动力学问题?经过经过Lanczos 变换后矩阵成为 三对角阵的证明 12.动力学问题?广义特征值方程的变形?使用变换 可得方程?可见Tn特征根是广义特征根问题的倒数 12.动力学问题?由于截断误差Xi并
22、不一定是正交?为了计算效率,而且多数情况下,只需计算一部分低阶特征值,因此变换只需进行q(n)步,这就是截断的Lanczos 变换 这样Tq是原问题的子空间,类似于Rayleigh-Ritz 法、子空间迭代法。12.动力学问题 12.5 振型叠加法(Modal Superposition)(一)固有振型及性质 对于无阻尼的自由振动问题的运动方程为 设 有 求解方程,得n个固有频率和特征向量 其中 有限元法基础 49?0MqKq()i tte?q2()?-0KM q,(1,)iiin?120n?1Tii?M12.动力学问题 有限元法基础 50 根据求特征根的方程,有 两式分别左乘 和 后相减,得
23、 当 不为零时,有 22iiijjj?MKMK22()0Tijji?MTj?Ti?22ij?0Tji?M10Tjiijij?M固有振型关于M正交 12.动力学问题 有限元法基础 51 利用特征向量的正交性,可得 定义 则有 20Tijiijij?K212212200nn?TTM?12.动力学问题 有限元法基础 52(二)系统的动力响应 1.位移基向量的变换 以特征向量表示位移 表达式的意义是将q(t)看成 线性组合,而 看成是广义的位移基向量,xi是广义位移值。代入系统的动力学方程,并利用 的正交性质,得 初始条件为 1()()niiittx?qx?i?i?()()()()()TTttttt?
24、xCxxQR?0000TT?xMqxMq?12.动力学问题 有限元法基础 53 设阻尼为振型阻尼,利用 正交性质 其中 为的i阶振型阻尼比。这样方程解耦,成为?20Tiijiijij?Ci?2()2()()()(1,2,)iiiiiiix tx tx tr tin?每一个方程相当于一个单自由度系统的振动方程每一个方程相当于一个单自由度系统的振动方程 12.动力学问题 有限元法基础 54?特例 1)设Q(t)可分解为空间函数和时间函数表示 如果F(s)与 正交,这表明系统中将不包含 响应成分,也就是说Q(s,t)不能激起与F(s)正交的振型。()()()()()()()Tiiitsttstft?
25、iQFrFi?0,()0iix t?ri?12.动力学问题 有限元法基础 55 2)如果对 作Fourier 分析,可得到所包含的各个频率成分及幅值。根据其中应予考虑的最高阶频率 可以确定进行积分的最高阶 ,例如选择 。综合起来,通常在实际分析时,求解的单自由度方程数远低于系统的自由度数n。()t?p?10p?12.动力学问题 有限元法基础 56 2.求解单自由度系统振动方程 杜哈美积分时将任意激振力分解为为冲量的连续作用,分别求出个系统的响应,然后叠加起来,即 ai和bi由初始条件确定。()01()()sin()(sincos)iiiitttiiiiiiiix tr t etdeatbt?2
26、1iii?一般杜哈美积分需数值积分计算 12.动力学问题 有限元法基础 57 3.振型叠加得到系统响应 获得每个振型的响应后,将它们叠加起来,得到系统的响应,即 a)b)c)振型迭代法不使用于非线性系统。1()()niiiq tx t?12.动力学问题 有限元法基础 58 3.振型叠加得到系统响应 获得每个振型的响应后,将它们叠加起来,得到系统的响应,即 在实际运用中,所取的振型数远小于n,能大大提高计算效率。1()()niiiq tx t?12.动力学问题 有限元法基础 59?特点 a)振型叠加中使用n个单自由度方程求解,应与直接积分的结果一致;b)振型叠加法比直接积分法节省时间,尤其是在选
27、取少量的单自由度方程的情况;c)振型迭代法不使用于非线性系统。12.动力学问题 12.6 减缩系统自由度的方法 12.6.1 Guyan减缩法 Guyan减缩法又称为主从节点法。设节点位移列阵q,分为主自由度qm和从自由度qs两部分,并设 ns和nm分别为qs和qm中的个数,有 有限元法基础 60 sm?qTq*mmms?qIq=qT qqT12.动力学问题 12.6 减缩系统自由度的方法 12.6.1 Guyan减缩法 Guyan减缩法又称为主从节点法。设节点位移列阵q,分为主自由度qm和从自由度qs两部分,并设 ns和nm分别为qs和qm中的个数,有 有限元法基础 61 sm?qTq*mm
28、ms?qIq=qT qqT12.动力学问题 考虑无阻尼自由振动,方程 代入关系式,并右乘 ,得 有限元法基础 62 Mq+Kq=0*TT*()()TTmmTMT q+TKTq=0*()()TT?MTMTKTKT系统方程从n降为nmm 12.动力学问题?Guyan法以静力减缩的方式,导出为主自由度qm和从自由度qs关系式,即 有限元法基础 63 00mmmsmsmsss?KKqKq=KKq1sssm?T=K K10smmsssssssmm?K qK qqKK q*1111*1TmmsssmsmsmsssmsmsssssssmTmmsmsssm?MMK KMMK KK K M K KKKK K K
29、12.动力学问题?特点 1)减缩后的方程,带宽会增加,只有采用较多的从自由度才能给计算带来明显的好处;2)主从自由度关系矩阵使用静力分析中的凝聚原理建立起来;对质量矩阵的减缩实际上是假定将从自由度上的惯性力项按静力等效的原则转移到主自由度上,这只有在从自由度的质量较小、刚度较大,以及频率较低时才合理。有限元法基础 64 12.动力学问题?使用使用GuyanGuyan法法 分析悬臂方板分析悬臂方板 有限元法基础 65 12.动力学问题 12.6.2 12.6.2 动力子结构法动力子结构法 又称为模态综合法。主要分析步骤:1)将总体结构分割为若干子结构 如同静力分析中子结构法,将结构分为若干个子结
30、构,仅在交界面上连接;2)子结构的模态分析 建立子结构的运动方程,建立模态基向量 a)固定界面主模态 在完全固定交界面上的位移条件下,求子结构系 统的固有频率和振型 有限元法基础 66 12.动力学问题 b)约束模态 在界面完全固定条件下,依次释放界面上的每一个 DOF DOF,并令它们取单位值所得的静态位移 c)用模态坐标表示子结构的物理坐标 3 3)集成各个子结构的运动方程 利用界面上的约束模态坐标集成,这样保证界面上位移协调条件 4 4)求解整个结构系统的运动方程 5 5)由模态坐标返回到各个子结构的物理坐标 有限元法基础 67 12.动力学问题 例:空间结构分析 由16根长为30cm,直径8cm 钢杆组成的空间结构,有 限元离散为32各2节点单元,供28各节点,144DOF。采用动力子结构法分析 时分上下两个子结构。方案1:上取8主模态,下16主模态 方案2:上下各取16个主模态 有限元法基础 68 12.动力学问题 空间结构分析的固有频率空间结构分析的固有频率 有限元法基础 69