基于节点里兹势能主自由度的结构动态缩减方法

毛虎平，高鹏飞

（中北大学 能源动力工程学院,太原 030051）

复杂结构动力学分析需要对巨大自由度的模型进行计算，当高频激励力作用时，要求计算步长非常小，这将造成计算耗时指数级增加。为了提高计算效率，可在保证一定精度的情况下，用少量自由度模型代替大量自由度模型，即模型缩减。所谓模型缩减是通过一定的变换，将对总体结构动力学分析影响较小的次自由度用对总体结构动力学分析影响较大的少量自由度表示，达到减少自由度的目的，其中少量自由度就是主自由度。然而如何从庞大的自由度中选择主自由度，目前在结构动力学领域仍属极具挑战性的问题。不过学术界目前提出一些选择主自由度的原则，最具代表性的有：

(1)将结构振动方向定为主自由度；

(2)在质量或转动惯量相对较大而刚度又相对较小的位置选择主自由度；

(3)在施加力或非零位移的位置选择主自由度。

这些原则在具体选择主自由度时，仅仅是指导思想，随意性较大。主自由度的位置和数目直接影响模态分析缩减质量矩阵的精度。针对无阻尼系统的各种有效合理PDOFS选择方法不能直接推广到阻尼系统中的问题，提出基于自由度能量分布比值的阻尼系统PDOFs选择方法，并通过算例说明了该方法的有效性和可靠性[1]。由于结构的动力行为受低阶模态控制，不需要计算系统的所有特征值，提出了一种单元级能量估计方法，建立了简化的有限元模型，该方法有效地节省了计算时间，并能准确预测全局系统的特征值[2]。罗虹等[3]提出2种主自由度选择方法，以单层悬臂梁为对象分析了这2种方法的特点和适用范围。刘孝保等[4]对多种主自由度选择方案进行了分析研究，表明结构静态缩减模态分析中误差的最大影响因素是主自由度数量和分布，特别是主自由度分布。包学海等[5]以转向架为分析对象，提出了选取主自由度的部分准则。

在选择主自由度后，面临的问题就是模型缩减。最早的模型缩减法是Guyan[6]提出的，称为Guyan缩减法，也称为静态缩减法。该方法忽略了自由度相关的惯性项和阻尼项，质量矩阵、刚度矩阵、状态向量和载荷向量均被分为主自由度和次自由度两个部分，然后经过矩阵变换，用包含主自由度的部分表示包含次自由度的部分。考虑惯性项，就产生了改进的缩减方法，通过Guyan缩减法获得惯性项，其结果与整体结构模态更加接近[7-8]。通过能量估计选择候选元素来构建简化系统，在所选元素关联的候选自由度中通过依次删除来选择主自由度，该方法有效地节省了计算成本，并能从最低频率到截尾频率范围内恢复整个系统的高精度特征值[9]。在模态叠加法中，所求的需要叠加的模态与载荷完全没有关系，实质上有些模态可能贡献很小，因此可以考虑用里兹向量叠加，因为里兹向量更能体现结构动态特性[10]。将简化系统与子域方案相结合，对每个子域进行尺寸与形状优化，由于在每个子域中都采用了约简方案，因此该方案对大规模问题的设计优化问题非常有效[11]。针对Guyan缩减法的不足，提出了改进的思路并推导了相关公式[12]。将里兹向量叠加与静态子结构法结合起来，形成了动力分析的子结构方法[13]。在工程结构损伤识别中应用逐级近似模型缩聚法，采用改进Guyan递推缩聚法的一级缩聚模型获得识别精度最高的结果[14]。提出了一种求解特征问题迭代凝聚的加速方法，采用序列消除法或能量法选择主自由度，并寻求了不同缩聚过程中方程的系统推导和比较，矩阵更新不仅包含了逆迭代，还隐含了子空间变换[15]。用结构各阶模态的DC增益作为其价值判断的准则，进行模态截断并实现模型缩减[16]。

通过文献分析可以看出，人们对结构动态缩减法以及对其影响最大的主自由度的选择方法均进行了大量研究，然而主自由度的选择一直以来没有一种确切的精度高的方法，鉴于此，本文提出一种基于节点里兹势能主自由度的结构动态缩减方法，利用里兹向量与结构自身动态特性和结构所承受的载荷分布形态相关联的特点，定义了节点里兹势能的概念，并在此基础上给出加权系数的公式，将两者点乘获得了节点里兹势能向量，以其为依据选择主自由度，最后用改进的动态缩减方法获得小规模的结构动力方程并求解。

1 节点里兹势能计算及主自由度选择

模态叠加法是计算结构动态响应的一种有效方法，其利用振型的正交性将动力学方程解耦，分别求解每一个方程再叠加起来。然而在参与计算的模态中，有些对响应影响较小，这是因为振型与结构所受载荷分布方式无任何关系。里兹向量是一组正交的并与载荷空间分布有关的向量，其是通过将质量矩阵归一化并正交化获得的，能够反映惯性力的影响。里兹向量构造过程如下：

（1）求初始向量{x1}。其可通过对[M] 归一化处理获得，即

其中：

（2）构造迭代式。

里兹向量与模态向量是相对应的，在前k阶里兹向量里分别取p1,p2,…,pk个最大里兹向量分量对应的自由度，将其组合并删除重合项，获得最终的主自由度，这称为基于里兹向量的主自由度选择方法。

节点里兹势能是指将模态空间转换到里兹向量空间，用里兹向量与节点自由度质量向量点乘而得到的势能向量。节点质量是将质量矩阵每一行元素求和作为结构有限元节点自由度的质量。选取节点里兹势能向量较大的分量对应的自由度作为主自由度。然而这样做会过分强调低阶频率，并忽略高阶频率，因此可以通过定义加权系数来减低计算结果过分集中于低阶频率的问题。节点里兹势能计算公式如下

2 构造缩减系统

其中：[K] 为n阶刚度矩阵，[M] 为n阶质量矩阵，ω

其中：φp为主自由度对应的特征向量；φs为次自由度对应的特征向量；下标p表示主自由度对应的量，下标s表示次自由度对应的量。

结构稳态简谐响应可表示为

其中：Ω为简谐激励力频率，f为激励力。

如果fs=0，根据式（12）可以获得精确的缩减关系

应用二项式定理将式(13)展开并省略2阶Ω以上的项可得

当Ω=0时

那么，静态缩减转换矩阵为

通过以下近似来消除Ω

3 基于节点里兹能量主自由度的结构动态缩减实施

图1说明了节点里兹能量主自由度的结构动态缩减流程，其实施步骤如下：

(1)构造里兹向量。通过式（1）至式（5）构造里兹向量，第1个里兹向量可以通过质量矩阵对角元素与刚度矩阵逆矩阵相乘，并与质量矩阵归一化获得。

(2)计算节点里兹势能。节点里兹势能由节点质量与里兹向量对应的分量乘积获得，其表示结构动态特性的贡献率，较大者说明贡献大，较小者贡献小，因此将此作为选择主自由度的依据。然而，这样做过分强调低阶频率，通过定义加权系数可以提高高频的精度。

(3)主自由度选择。通过式（6）计算获得节点里兹势能向量，选择其分量较大者作为主自由度。

(4)构造缩减系统。通过IRS方法，在静态缩减法的基础上，考虑结构惯性力，该惯性项能使模态结果更加逼近完整模型的模态，它采用静态缩减法获得。式（20）是最重要的转换矩阵，式（8）至式（19）为其推导过程。

(5)采用广义Schur分解法[17]求解该矩阵。对于任意n阶矩阵A，存在一个酉矩阵U，使得U′AU成为上三角矩阵，且该上三角矩阵的对角线元素为A矩阵的特征值，利用该性质进行缩减系统求解。

4 实例分析

4.1 圆柱曲板

圆柱形曲板半径为100 mm，高为100 mm，两侧端固定，弹性模量为0.3 MPa，泊松比为0.3，密度为0.01 kg/m3，采用SHELL63单元，共有216单元，247个节点，1326个自由度，如图2所示。

根据里兹向量和节点里兹势能选择主自由度，并将与主自由度相连的单元显示出来（其中颜色较深的为主自由度对应的单元），图3为选择90个主自由度时2种方法的不同效果。

图1 基于节点里兹能量主自由度的结构动态缩减流程

从图3可以看出，基于里兹向量法选择的主自由度过于向结构中心集中，这是过分强调低频的一个重要表现；而基于节点里兹势能法选择的主自由度大部分也集中在结构中心，但有一部分主自由度向两侧扩散，这可能是提高高阶固有频率的一个重要现象。如图4所示。

图2 圆柱曲板几何参数

图3 圆柱曲板上选择主自由度对应的单元

图4 根据不同主自由度选择方法所得圆柱曲板计算结果相对误差比较

从图4可以看出，在选择相同主自由度数的情况下，任意选择主自由度误差太大，而基于里兹向量法选择主自由度误差非常小，在前30阶模态中，最大相对误差不超过10%，然而基于节点里兹势能法选择主自由度误差最小，在前30阶模态中，最大误差3%。

表1是对应的固有频率数值，从中可以看出本文方法的优势。如图5所示。

从图5可以看出，本文方法选择不同主自由度数时，结构模态的相对误差变化随着自由度数量的增加越来越小，当主自由度数为450时，前30阶模态相对误差最大不超过2%；另一方面，自由度数量增加相同，而结构模态相对误差减小得越来越慢，450个自由度大约是全部自由度数的1/3，如果继续增加主自由度数，虽然可以继续减小误差，然而对于工程应用来说意义不大，因此在使用本方法时只要误差满足精度要求即可。

表1 根据不同算法所得圆柱曲板模态比较（400个主自由度）

4.2 曲轴

曲轴轴径为6 mm，曲柄直径为4 mm，曲柄长4 mm，曲轴总长34 mm，两端固定，弹性模量为210 GPa，泊松比为0.3，密度为7850 kg/m3，采用十节点四面体结构实体SOLID92单元，共有1022个单元，2127个节点，6231个自由度，如图6所示。

图5 根据不同主自由度数所得圆柱曲板计算结果相对误差比较

图6 曲轴几何参数

图7是根据2种方法选择的500个主自由度对应单元（颜色较深的单元），从中可以看出，基于里兹向量法选择的主自由度集中在结构中部，而基于节点里兹势能法选择的主自由度从中部向两端扩散，而且非常明显。

图7 所选择曲轴主自由度对应的单元

图8是采用随机选择、基于里兹向量法和基于节点里兹势能法分别选择1800个主自由度并用IRS法构造的缩减系统，求解获得的曲轴模态相对误差，可以看出，随机选择不可行，基于节点里兹势能法选择主自由度结果相对误差最小，前30阶模态中，最大相对误差不超过10%，而基于里兹向量法选择主自由度虽然优于随机选择方法，但前30阶的相对误差中，最大误差超过20%。

图9为基于节点里兹势能法选择不同数量的主自由度的计算结果，进一步说明随着主自由度数量的不断增加，结构模态相对误差不断减小，但是减小的幅度越来越小，而主自由度数为2000时，在前30阶中，结构模态相对误差最大不超过5%。

表2提供图8对应的具体数据，更能从数值上看到本文方法的优势。

表2 模态计算比较（2000个主自由度）

图8 根据不同选择主自由度方法所得曲轴计算结果相对误差比较

图9 根据不同主自由度数所得曲轴计算结果相对误差比较

5 结语

本文通过里兹向量法和节点里兹势能法选择主自由度，并采用改进缩减法构造缩减系统，最后应用广义Schur分解法求解，通过圆柱曲板和曲轴分析说明了基于节点里兹势能法选择主自由度方法的有效性，并获得以下结论：

(1)将模态空间转换到里兹向量空间，可以用很少的几个里兹向量捕捉到非常精确的动态特性。由于其过分强调低阶频率，需进一步通过定义加权系数来提高高阶频率的精度。

(2)节点里兹势能法比里兹向量法能捕捉到更加合适的主自由度，从而能更好地反映结构动态特性，在同等条件下，能获得更高的精度，为下一步高效的动态响应优化奠定基础。

(3)在结构缩减中，用大约1/3的主自由度数量比较合适，如圆柱曲面板总自由度数为1326，而取450个主自由度获得了前30阶模态相对误差不超过2%的精度；曲轴总自由度为6231，取2000个主自由度获得前30阶模态相对误差不超过5%的精度。