结构动力方程的精细时程积分法

第34眷第2期 1目91年4月 大 连 理 工 大 学 学 报 Journal 0f~ lian University of Technology Vo1．34，Not 2 Apr - 1 9 9 4 G I 卜 结构动力方程的精细时程积分法。 。3 。 钟 万 勰” ， — — — — — — — — 一 (工程力学研究所) 摘 要 对线性定常结构动力系统提供了精细时程积分法．它的解在积分点处， 在数值上逼近于精确解的强值结果．数值例题验证了方法高度精确的特点． 关量词：结构 力学{积分法／时程积分法 分类号： 鬲 —— r 1 问题的提 出 工程结构在倒如突加荷载或冲击荷载的作用下，往往要求作瞬态历程的分析．由于结构几 何形状的复杂性，在对空间坐标采用有限元法离散后，得到的往往是数十、百、千 自由度的系 统．采用本征向量展开法普通只选用低频的若干个本征解，这是不够的．此时就要采用在时域 中逐步积分的方法． 逐步积分法已经有了多年的研究[1 a]，有许多种方法可供选用．从分类上说有显式积分与 隐式积分两大类．显式积分对于每一时阿步效率较高，但时间步必须取得很小方能保证其稳定 性 隐式积分法则可以通过恰当的参散选择，保证积分的数值稳定性，因此时间步长 delt长一 些也可适用．熟知的 Newmark法及 W'~,on-9法都是隐式积分格式．然而这些积分法也有其弱 点}由于delt选得大一些，一些高频振动的分量不能正确她反映出来．又如果将其用于一个保 守体系，则经过若干步的计算之后，系统的能量等不能保持守恒．这种情况被称为“人工阻尼 或。算法阻尼 I这是积分法选用所带来的，是无可奈何的． 对于最常见的动力体系方程 + G + 五 ；r(1)． (0)． (O)给定 (1) 其中： 、G、x为 × 矩阵， 待求，r(f)为给定外力．应当尽力设法寻求一类积分方法，对于 定常体系，可使系统积分的结果高度准确．当村 为对称正定，G为反对称， 为对称且 r=0 时，即系统为保守时，积分结果保持系统的守恒量不变．这类积分方法虽不是提供解析解公式， 但其数值结果却是高度准确的，因此可称之为精细分析}这也就是本文的目标． ．：蛊累皇羹莩差善鑫 国皇家学会联台资助项 收蔫日期 -1993-07—2．5 维普资讯 http://www.cqvip.com 大 连 理 工 大 学 学 报 3‘眷 2 方程的变换与时程积分 在研究偏微分方程与计算结构力学对最优控制理论的模拟关系时[ ，对于哈密顿体系 的积分推出了一套有效的算法．虽然方程(1)已经不是哈密顿体系，但其基本方法与思路仍可 用于当前课题． 仿效哈密顿体系对偶变量的引入，令口 p— M + Gx／2 或 = M～P—M Gx／2 (2) 此时，式(1)成为 一 一 ( 一 GM—G／4)x— GM p／z+ r (3) 以上方程可以写成为线性体系的一般形式 口= + Dp+ rq ‘ (4a) 一 脚 +cp+ rp (4b) 其中各个量对比于方程(2)及(3)为 口：x，P—P， A一一M G／2， B=GM G／4一K 1 ⋯ C= 一 GM一／2， D — M一， ，p一 ，， = 0 J 当M 阵对称正定，G反对称， 为对称，且 =， =O时，体系成为保守，即成为自由陀螺系统， 此时B、D为对称矩阵，^一一c ．对于陀螺系统的分析推导及本征问题的计算可见文献[7]． 保守体系可运用变分原理来推演，并充分利用哈密顿体系的已有成果．但当有阻尼时，方程(4) 并不是保守体系I本征函数的频域法用起来就感到不便 因此，采用时程积分是很自然的． 方程组(4)是非齐次的．从线性常教分方程的理论知，应当先拽出齐次方程组的通鼻，然后 再用叠加原理或变量代换法找出其非齐次方程的特解．因此，应首先着重研讨齐次线性方程 口一 由 + D， (6a) j=脚 +cp (6b) 常微分方程组常用的效值 计算方法 煤矿单位产值综合能耗的计算方法营养成分理论值计算方法电缆末端电压降计算方法初中24点计算方法与技巧答案24点计算方法与技巧下载 往往是将微分方程化成差分方程的积分法．然而如文献[8] 所指出的，这些差分方程往往会破坏保守体系的守恒性质{这是难于接受的· 将式(6)合并写成为 ， 一 ， Ⅳ = A ， 一 }：) ㈣ 对于定常系统，H是常矩阵，其通解可形式地写为 U— exp(H ·t) ． (8) 令时间步长 dt／t=fv，则 u(r)： exp(H ·r)·uo T · (8 ) 其中： 一 exp(H ·r) (9) 问题就归结到 阵的计算I只要很精细地算出了 阵，则时程积分就成为 一 Tuo； uz一 ，⋯ ； 一 丁 一l'．-· (10) 维普资讯 http://www.cqvip.com 2期 钟万勰 ：结构动力方程 的精细时程分法 l33 一 系列的矩阵乘法了 3 指数矩阵的精细计算 文献Es]已经给出了指数矩阵的精细计算方法．其要点是利用加法定理，有 exp(H ·f)一 Eexp(H ·v／m)]- (11) 其中可选用 一 2 }例如 Ⅳ一20，则m=1 048 576．由问题看出r=de~本来是不大的时间区 段，从而At=r／m将是非常小的一个时间区段了．对于 的时间区段，有 exp(H ·At)≈ ，+ H · + (H ·At)0／2= ，+ (12) 其 中： 一 (H ·At)·(J+ (H ·At)／2) (13) 在数值计算时至关重要的一点是，矩阵的存储只能是式(13)的 ．而不可以是(12)的』+ ． 因为 是非常小的值；当它与 1相加时．就成为尾数f在计算机的舍入操作时其精度将丧失殆 尽．注意倍精度数的表示范围是 lO ，所以 的值不会消失．另外应注意倍精度数的有效位 数为 l6位十进数，式(12)、(13)的展开式中的项取得太多已无太大作用，然有益无害．究竟用 几项还应当依课题而定．本文的例题计算都是用式(13)的． 为了计算 丁阵先应当将式(11)作 T=[ + z =[ +』 X[ +J] (14) 的分解．其次应当注意 (J+ )X (J+ )一 J+ + + × (15) 因此式(14)相当于执行语句 for(iUr= 0；衙 < N~／t,T++) 一 2T,+ L X f (16) 当循环结柬时有 T — I+ (17) 式(13)、(16)、(17)便是指数矩阵的精细 计算公式 六西格玛计算公式下载结构力学静力计算公式下载重复性计算公式下载六西格玛计算公式下载年假计算公式 ． 还应当对于此算法的精度作一分析．设 阵有 HY=yrvj或 y 肿 = J (18) 其中ty为本征向量所组成的矩阵， 为相应的本征值．于是就可以导出 exp(H ·At )一 Yexp(rpjAt )Y一 一 Yrexp(p ·At )J Y一 (19) 于是式(12)、(13)的近似相当于 exp(p·At)≈ 1+ · + ( ·At) ／2 因此计算精度的需要应当使 · 很小．当然应当选用最大的abs( )· 来验证．由于 Ⅳ一 20的选择，即使有 ddt×Bbs( )一i00X 2x，即de／t有 100个周期，仍有 At~"-lO 的周期．划分 仍是很细的． 式(16)、(17)对于迭代后有 一一J的情况就丧失了其相对精度}出现这种情况意味着该 体系的本征解是高度阻尼的，在 ddt的时间区段内已将近衰减掉了．这类解本来也是不太关 心的，因此不必在意．如果关心很短时间内的系统反应，则ddt本来就应选得更小的． 维普资讯 http://www.cqvip.com 大 连 理 工 大 学 学 报 34卷 4 非齐次方程 可以认为非齐次项是在时间步 (tj，“+ )内线性的，即方程为 = Hu+ r0+^( 一 tDI 当 t= 时 ，u一 (20) 式中： 、r 是 2 维给定向量．该方程可以用叠加原理求解{令 0一“)是齐次方程的解，即 一 H ， (。)= J (21) 于是可以写出 u一 (t一 )[ + H一 (to+ H一 r1)]一 H一 ( + H一 r1)一 (爿叫，、)(t= tD (22) 代入式(20)，利用(21)，可知(20)的微分方程与初始条件都已满足- 数值计算中虽没有 0一t。)的解析式，然而逐步积分要求提供的是 tt+l=t·+r时刻的向 量 +1，此时 ( Hl— tD 一 (r)一 T (Z3) 而 是已算得的．因此 +1一 丁[ + H一 (to+ H一 ，、)]一 H一 [ro+ H叫r1+ ，l‘f] (24) 这就是所耍推导的时程积分公式． 5 数 例 第 1个数倒选用文献[2]中的倒 8．2．方程为 ㈥㈡+[ _42]{ ㈦ 。0) 汹， delt 0．28．经过精细时程积分的计算，见表 1．在有效位数内与精确解完全一致．这里应当指 出，文献[2]的精确解的数值列表是不精确的．对比是根据解析式重算的- 表1 精细时程积分计算结果 003 0 038 0 1 0 486 0 996 1 657 338 2 1 052 805 ， m ． ． 76 ． - · 2· ·85 3· Z- z。l l l。l ， 。 o．382 1．412 z_181 4"．094 ‘996 11 塑 !：!!! ． ： !! !：!竺 !： 坠 !：! ! 第2个数例见图 1，质量阵M：diag(8)，而 G = 0 0 O K = 4 — 4 — 4 8 — 4 - 一 4 8 — 4 — 4 8 — 4 — 4 8 — 4 O 一4 8 — 4 o_ 维普资讯 http://www.cqvip.com 2期 钟万勰 t结构动力方程 的精细 时程分法 韧始条件为aT 一 一 z 一0，两一10} = _．．·= =0，当t一0时．这个课题的特点是质 心守恒．选择 r=delt=1以及 r一4分别作计算，得到的相同时刻的反应完全相同．并且还可以 看出，其质心完全不动，即满足了守恒律．课题中有一阻尼器，所以能量并不守恒；其振动是逐 步衰减的．表 2列举了某些时刻的位移}它表明在合理的范围内选择时间步长r：delt，对于结 果没有影响．这表明了精细积分法的可靠性，其数值结果当然应当与精确的数值结果相符合． 表 2每行的(为+⋯l r 都等于 10．n．表明了质心守恒律． 磊LJ_；； ．r；； J L j ； J ～ t 了 k 鬲 k 瓣 }k k )； l嬲 }l k } l — I_ 靠 捌 L ， b kg． J_、 n1，g— j N nlJ 表 2 时程积分的输出结果(部分) 6 结束语 结构动力学时程积分是应用中的重要环节．本文对定常体系给出了精细时程积分法，在数 值结果上可与精确解的数值结果相比拟．其计算高教的原因是充分利用了定常幂统对时间的 均匀性，即动力方程对于时间坐标的平移群是不变的性质．因此适用 2 的算法}因而可将时间 区段划分得足够小，从而使数值结果足够精确． 这种方法对于分段为定常的体系也可运用，类同的方法在LQ控制中已有论述嘲． 作者对孙焕纯、林家浩教授有益的讨论深表感谢 维普资讯 http://www.cqvip.com 大 连 理 工 大 学 学 报 3‘卷 参 考 文 献 1 孙焕纯．曲乃滔．林束浩．高等计算结构黄力学．大连 大连理工大学出版社，1991· 2 Blithe K J，Wilson E L．Nasmsrkd Methods如 Fibre捌 A~ ysis．New York=Premice-Hall·1976· 3 subba"j K．Dokalnish M A．A survey of db~ct time integ~atlon methods in computetiou~t structural dynamics：1 Explicit methods，1 Implicit methods．Computers＆ sIruHt·w ，1989，32(6)：1371～1386f1387 ～ 1401 ’ 4 钟万艇．钟期期．计算结构力学，曩优控制及壤截舟方程牛解析法．计事结构力学及其应用，1990·7(1)：1～ l5 5 钟万勰，钟朔期．LQ拉削区段混合能矩阵的微分方程及其应用．自神化学报，1992，18(3)；325~331 6 钟万蕊，橱再石．连壤时间LQ控制主要丰征对的算法．应用数学和力学·1991，12(1)：45~50 7 钟万艇．林索浩．疙螺系统与反对称矩阵辛丰征鬻的计算．计算结构力学及其应用·1993·10(3)：237~253 8 冯 康，橐盂兆．Hamilton神力体幕的HamilUm算法．茸热科学进展，1991，1(2)=110~120 On precise time-integration method for structural dynamics Zhong W amde" (Research Institute Engir~ecdng Mechanics，DUT) Abstract The precise time-integration method is proposed for linear time—invariant structura1．dynamic systems．Its numerical results at the integration polat~are a／m ost equal to the results of the exact solution．Numer~ examples show the highly precise~hlu'acteri~ics 0f the present method· Key Wordsl structural dynamics；integratlons／t／me-lntegration method ．：HJoi唧nt pcoryjectPld tehse∞Nra ， t u ura n l iv ~ e 哼 F0fouwnd a a I t 瞄 imx c d o1] C h ge ina 0f an ( d 二a t rd he m Roy s0dety o‘u K 维普资讯 http://www.cqvip.com