大规模动力系统改进的快速精细积分方法.docVIP

大规模动力系统改进的快速精细积分方法高强，张洪武*，钟万勰 (大连理工大学，工业装备结构分析国家重点实验室，工程力学系，辽宁大连 116024） 1 引言 对于线弹性结构的动力学问题，比较成熟和常用的时域积分方法是Newmark方法[1]和Runge-Kutta方法[2-4]，由于计算稳定性和精度的要求，这两种方法的积分步长一般都取得比较小，计算量比较大。钟万勰提出了精细积分方法[5-7]，这种方法计算精度非常高，稳定性好，允许采用很大的积分步长，特别是在处理刚性具有明显优势。精细积分方法提出后，得到了很多发展[8-10]，但是这种方法的一个缺点是在处理规模很大的系统时，由于计算矩阵指数的计算量比较大，效率是一个主要问题。 本文针对大规模动力系统提出一种计算其动力响应的高效率算法。本文以精细积分方法为基础[5-7]，利用动力问题的物理特性，一个关键思想，也就是结构动力问题能量传递的有限性，来降低矩阵指数的计算量。利用上述思想，本文分析了大规模动力结构对应矩阵指数的稀疏性质，并基于此给出一种计算矩阵指数的高效率方法。在高效和精确计算矩阵指数的基础上，给出了大规模动力系统响应的高效率和高精度算法。 2动力系统的精细积分方法 假设系统的刚度矩阵、阻尼矩阵和质量矩阵分别为，和，那么结构动力学方程为 其中为外力。方程可以写为状态空间形式，即 其中 其中为动量。 数值积分时，将积分区间等分成步长为的积分区间，即 若记 则方程的解可以表示为 其中 通过方程计算结构的动力响应，需要解决两个主要问题，一是要准确高效的计算方程定义的矩阵指数，二是要计算方程中的积分。对于常见的外载荷，方程中积分大多可解析积分，例如 如果外力在一个积分时间步长内多项式，即，则 其中 上式中和分别为矩阵指数对应的分块矩阵，即 如果外力在一个积分时间步长内为简谐变化，即，则 因此，上述计算过程的关键是计算矩阵指数。目前，计算矩阵指数最好的方法是精细积分方法[5-7]。但是，精细积分法计算矩阵指数的计算量比较大，即使对于为稀疏矩阵的情况，也很难利用其稀疏性，得到的矩阵指数可能是满阵。 3 改进的快速精细积分方法 精细积分方法基于两个要点，一个是加法定理，另一个是增量计算和存储。对于给定矩阵，它对应的矩阵指数有如下性质，即 如果足够大，则矩阵比较小，那么矩阵的矩阵指数可用如下的阶Taylor级数近似，即 精细积分方法[5-7]将的矩阵指数分为两部分，即 然后对增量部分应用加法定理，即 通过执行次方程，然后将加上单位矩阵，则得到对应的矩阵指数。 上面简要地介绍了精细积分方法，此方法编写程序，计算精度非常高但是对于，由于系统的自由度数目很大，计算矩阵指数将非常耗和内存。也是稀疏矩阵，但是直接通过以上的精细积分方法计算其矩阵指数，在计算过程，特别是方程的加法定理的合并过程中，矩阵将逐渐变为满矩阵或非常稠密的矩阵，很难利用矩阵的稀疏性质，从而计算量很大。 本文利用结构动力响应的特点，给出一种计算矩阵指数的算法，且可节省存储要求。 众所周知的事实是能量在一维杆中的传播速度是有限值，同理，在离散的结构中，虽然其能量的传播速度很难确切得到，但其动力学响应的能量传递速度肯定是有限的，这将对矩阵指数的结构产生很大影响。根据方程，如果外力为零，则方程可写为如下分块形式，即 由上述方程可以得到矩阵指数元素的物理意义，即：的第行第列元素表示第个自由度上给定单位位移，而其他自由度位移为零且所有自由度为零时，经过一个积分步长后，第个自由度的位移响应；的第行第列元素表示第个自由度上给定单位，而其他自由度为零且所有自由度位移为零时，经过一个积分步长后，第个自由度的位移响应；的第行第列元素表示第个自由度上给定单位位移，而其他自由度位移为零且所有自由度为零时，经过一个积分步长后，第个自由度的响应；的第行第列元素表示第个自由度上给定单位，而其他自由度为零且所有自由度位移为零时，经过一个积分步长后，第个自由度的响应。 由于能量在结构中的传递速度是有限的，假设某个自由度上有扰动，在，必然只能传播到有限的自由度，而不可能传播到所有自由度。根据上面给出的的物理含义，则它们一定是稀疏矩阵。这样，既可以将矩阵指数作为稀疏矩阵计算和存储，从而节约计算量和存储空间。对于给定矩阵和积分步长，按如下过程计算矩阵指数：，然后令，其次按照方程计算，然后执行N次方程，最后将与单位矩阵相加得到矩阵指数。 若采用上述的计算过程，虽然矩阵为稀疏矩阵，但是如果观察根据方程计算得到的矩阵，可以发现它可能是一个满矩阵或很不稀疏的矩阵，但仔细观察会发现，此矩阵有很多很小的元素。另一方面，此时的矩阵相当于很小的时间步长