CN101334803B - 基于变形修正的动力减缩算法及超单元构造方法 - Google Patents

Abstract

本发明公开了一种基于变形修正的动力减缩算法和一种超单元构造方法，在运动同步性假设的前提下，采用适当的位移模式去逼近各同步性区域的位移，得到整体位移与局部位移模式间的转换关系，进而用总体应变能最小的算法去优化整体位移与局部位移模式间的转换关系，可得到精确度极高的简化模型。基于这种简化模型，可构造出计算量大大降低的、各种各样的高效动态超单元。

Description

基于变形修正的动力减缩算法及超单元构造方法 

技术领域

本发明涉及结构动力学领域，尤其是基于变形修正的动力减缩算法及超单元构造方法。 

背景技术

动态系统模型简化方法的研究和应用，一直是大型动态系统理论研究和结构设计中的重要课题。同样，在结构动力学领域，模型简化也是结构动力分析中非常关键的技术。结构动力学模型简化的根本目的是为了获得一个满足工程精度要求的低阶、有效的计算模型，从而可以用简化模型对原始复杂模型进行性能分析以及模拟仿真。现有的结构动力模型简化方法主要可分成以下三类。 

(1)自由度减缩方法，基本思路是从结构一般运动方程或特征方程出发，用保留自由度表示缩聚掉的自由度，从而实现对模型的简化。典型的此类方法有Guyan-Irons法、Kuhar法、IRS法、和模态缩聚法等。 

(2)动态子结构方法，这类方法是直接得到低阶模型的方法。首先得到各个子结构的低阶动力特性，然后通过子结构间位移和力的双协调条件得到整体结构以低阶模态坐标表示的综合振动方程。 

(3)结构等效方法，此类方法从结构力学分析出发，针对某种特定的结构，用简单的结构等效复杂结构的主要特征，得到简化的力学模型。 

随着结构动力学研究领域的进一步深入，近些年来，又出现了参数识别、载荷反演、振动控制等很多新问题。这些新问题的求解对动力学系统的模型简化提出了更高的要求。首先，虽然现代计算机技术得到了很大的发展，但是计算量在某些问题中依然是一个瓶颈。其次，随着结构自由度的剧增，计算中引入的误差项可能导致解的不稳定，而无法得到解的最优极值。再次，大型复杂结构的控制、工况检测及故障检测等问题中，如果结构系统的数学模型的自由度太大，就无法实现“实时”。最后，如果问题求解中涉及非线性，同样会有计算量过大的困难。对于 以上诸多方面对模型简化提出的要求，现有的模型简化方法有时远远不能解决。 

发明内容

针对现有技术存在的问题，本发明的首要目的是提供一种基于变形修正的动力减缩算法，该算法以复杂工程结构的有限元模型为基础，根据运动同步性假设将结构分成若干同步性区域，采用适当的位移模式去逼近各同步性区域的位移，得到整体位移与局部位移模式间的转换关系，进而用总体应变能最小的算法去优化整体位移与局部位移模式间的转换关系，使得整体位移与位移真值最逼近。本发明进一步的目的是提供一种在上述模型简化基础上的超单元构造方法。 

为实现上述目的，本发明基于变形修正的动力减缩算法，具体为： 

1)根据运动同步性的假设将结构的有限元模型划分为若干个同步性区域，采用适当的位移模式去逼近各个同步性区域的位移，得到整体位移u与局部位移模式q间的转换关系： 

u＝Rq    (3) 

其中一种比较简单而有效的位移模式就是刚体位移模式。 

2)用总体应变能最小的算法来优化整体位移u与局部位移模式q间的转换关系，使得整体位移u与位移真值 

最逼近，由关系式： 

$U\left(\mathrm{\δ\ϵ}\right)=\frac{1}{2}{\left(\mathrm{\δu}\right)}^{T}K\left(\mathrm{\δu}\right)---\left(8\right)$

$=\frac{1}{2}{(\tilde{u}-\mathrm{Rq})}^{T}K(\tilde{u}-\mathrm{Rq})$

$\frac{\∂U}{\∂q}=0---\left(9\right)$

可得： $\stackrel{\‾}{u}=\mathrm{Tq}---\left(13\right)$

其中， 为位移真值的近似值，T＝K^-1((R^TR)^-1R^T)^TR^TKR。 

一种在上述模型简化基础上的超单元构造方法，具体为： 

1)根据原模型的几何形状，把结构划分为若干部件，每一部件作为一个超单元。再根据各部件的运动特征划分为若干运动同步性区域。保留部件边界节点并且在部件各同步性区域特征点处增设超节点； 

2)利用基于变形修正的动力减缩算法，得到超单元区域简化模型位移向量Q与相应有限元结构节点位移向量U的转换矩阵 

； 

3)把区域边界节点位移关系回代到简化模型中保证了超单元边界位移的协调性，得到超单元位移向量 

与相应有限元结构节点位移向量U的转换矩阵Δ；进而得到：超单元的质量矩阵M_s、超单元的阻尼矩阵C_s、超单元的刚度矩阵K_s、超单元的节点力向量F_s； 

4)得到超单元的单元矩阵后，根据总体结构超单元节点的总体编号，把各个超单元的单元矩阵集成为总的矩阵； 

5)在求解精度要求较低的情况下，可以对超单元的质量进行凝聚，只保留超单元的内部超节点，从而使得计算量大大降低。 

附图说明

图1为基于变形修正的动力减缩算法实施例的有限元网格图； 

图2为图1中保留边界节点简化模型示意图； 

图3为图1中不保留边界节点简化模型示意图； 

图4为图1实施例的两种简化模型的频率误差曲线； 

图5为超单元示意图； 

图6为图1中超单元模型示意图； 

图7为图6实施例的频率误差曲线； 

图8为图6中超单元质量凝聚模型； 

图9为图8实施例的频率误差曲线； 

图10为图1实施例的简化系统的频率结果与ANSYS求解的结果对比； 

图11为图6实施例的超单元模型频率结果与ANSYS求解的结果对比； 

图12为图8实施例的频率结果与ANSYS求解的结果对比； 

具体实施方式

基于变形修正的动力减缩算法： 

在一定的位移模式下，结构某一区域运动可以用某一模式进行描述，即所谓运动同步性。刚体运动自然是运动同步性的一种特殊情况，我们在每个运动同步性区域都用准刚体模态去简化，由于采用了刚体模态，故称这样的区域为刚体化区域。 

针对相应结构的有限元模型，先根据运动同步性的假设将结构划分为若干个同步性区域，然后用准刚体模态位移去等效每个区域的位移，得到最初的局部刚体化简化模型。于是，每个刚体化区域内任一点j的运动可以用该区域特征点的6个基本的准刚体模态进行叠加： 

u_j＝R_jq_j    (1) 

其展开形式为： 

$\left\{\begin{array}{c}{u}_{\mathrm{jx}}\\ u_{\mathrm{jy}}\\ u_{\mathrm{jz}}\end{array}\right\}=\left[\begin{array}{cccccc}1& 0& 0& 0& z_{\mathrm{cj}}& -y_{\mathrm{cj}}\\ 0& 1& 0& {-z}_{\mathrm{cj}}& 0& x_{\mathrm{cj}}\\ 0& 0& 1& y_{\mathrm{cj}}& {-x}_{\mathrm{cj}}& 0\end{array}\right]\left\{\begin{array}{c}{q}_{\mathrm{jx}}\\ q_{\mathrm{jy}}\\ q_{\mathrm{jz}}\\ q_{\mathrm{j\θx}}\\ q_{\mathrm{j\θy}}\\ q_{\mathrm{j\θz}}\end{array}\right\}---\left(2\right)$

式中，u_j表示该区域中第j点的位移向量，q_j表示该区域的刚体模态，R_j表示6个刚体模态的叠加矩阵或称位移变换矩阵，x_cj，y_cj，z_cj表示该区域中任一点j的相应坐标与该区域特征点相应坐标的差值。 

把整个模型的所有节点分成边界节点和内部节点两部分，内部节点的自由度全部用刚体化区域特征点的刚体模态凝聚，所以在刚体化区域特征点处增加内部超节点，即： 

其中，u_o代表边界节点的位移向量，u_in代表内部第n个同步性区域节点的位移向量，R_n代表内部第n个同步性区域节点的位移变换矩阵，q_n代表内部第n个同步性区域超节点的刚体模态，q代表内部自由度缩聚后的节点位移向量，R表示总体位移变换矩阵。 

在各个运动同步性区域用准刚体模态来假定位移模式，而实际情况可能不是，这就带来了误差。基于能量法和变分法，用总体应变能最小的算法来对简化模型进行变形修正，以达到总体优化的目的。 

假设位移真值为 

，依据误差位移引起的系统总误差应变能需达到最 小，推导最优的q，使得简化模型的位移u与位移真值 

最接近。 

位移误差向量为： 

$\mathrm{\δu}=\tilde{u}-u=\tilde{u}-\mathrm{Rq}---\left(4\right)$

其中，u代表(3)式所示的局部刚体化简化模型的位移向量。 

δu对应的应变误差为： 

δε＝Bδu    (5) 

B是一般有限元理论中的应变矩阵。应变误差δε对应的总体误差应变能为： 

$U\left(\mathrm{\δ\ϵ}\right)={\∫}_{\mathrm{\Ω}}\frac{1}{2}{\left(\mathrm{\δ\ϵ}\right)}^{T}D\left(\mathrm{\δ\ϵ}\right)\mathrm{dV}$

$={\∫}_{\mathrm{\Ω}}\frac{1}{2}{\left(\mathrm{B\δu}\right)}^{T}D\left(\mathrm{B\δu}\right)\mathrm{dV}---\left(6\right)$

$=\frac{1}{2}{\left(\mathrm{\δu}\right)}^T{\∫}_{\mathrm{\Ω}}{B}^T\mathrm{DBdV}\left(\mathrm{\δu}\right)$

而： 

K＝∫_ΩB^TDBdV (7) 

K表示系统总刚度矩阵，D为一般有限元理论中的弹性矩阵。 

可得： 

$U\left(\mathrm{\δ\ϵ}\right)=\frac{1}{2}{\left(\mathrm{\δu}\right)}^{T}K\left(\mathrm{\δu}\right)---\left(8\right)$

$=\frac{1}{2}{(\tilde{u}-\mathrm{Rq})}^{T}K(\tilde{u}-\mathrm{Rq})$

为使系统总误差应变能最小，需使： 

$\frac{\∂U}{\∂q}=-R^{T}K(\tilde{u}-\mathrm{Rq})=R^T\mathrm{KRq}-R^{T}K\tilde{u}=0---\left(9\right)$

即： $R^T\mathrm{KRq}=R^{T}K\tilde{u}---\left(10\right)$

用 

近似表示位移真值 

由(10)式可得： 

$K\stackrel{\‾}{u}={\left(R^T\right)}^{+}{R}^T\mathrm{KRq}$

$={\left({\left(R^{T}R\right)}^{-1}{R}^T\right)}^T{R}^T\mathrm{KRq}---\left(11\right)$

其中，(R^T)⁺代表R^T矩阵的广义逆矩阵。因为广义逆是一种基于最小二乘的近似，所以 

与位移真值 

还是不完全一致的。 

由(11)式可得： 

$\stackrel{\‾}{u}=K^{-1}{\left({\left(R^{T}R\right)}^{-1}{R}^T\right)}^T{R}^T\mathrm{KRq}---\left(12\right)$

令：T＝K^-1((R^TR)^-1R^T)^TR^TKR，则得简化模型优化后的位移变换表达式为： 

$\stackrel{\‾}{u}=\mathrm{Tq}---\left(13\right)$

利用(13)式的变换关系代入到系统的动力学方程： 

$M\stackrel{\·\·}{U}+C\stackrel{\·}{U}+\mathrm{KU}=F---\left(14\right)$

再在方程的两边同时左乘T^T，得到简化后的动力学方程： 

$M_R\stackrel{\·\·}{q}+C_R\stackrel{\·}{q}+K_{R}q=F_R---\left(15\right)$

式中，M，C，K为原始有限元结构的质量、阻尼和刚度矩阵；M_R＝T^TMT，C_R＝T^TCT，K_R＝T^TKT为优化后简化结构的质量、阻尼和刚度矩阵，F_R＝T^TF为优化后简化模型的载荷向量。 

基于变形修正的动力减缩算法的实施例： 

本实施例为一平面问题，一种材料。ANSYS单元类型：Plane42，单元尺寸为0.1m，几何尺寸0.5×5m，材料参数：E＝1.0×10¹⁰Pa，μ＝0.25，ρ＝4500kg/m³，左端固定。把模型均匀划分为十个同步性区域，在每个同步性区域质心处各放置一个超节点。有限元网格图如图1所示，保留边界节点简化模型示意图如图2所示，不保留边界节点简化模型示意图如图3所示。 

得到简化结构的质量矩阵和刚度矩阵之后便可以求得系统的频率，简化系统的频率结果与ANSYS求解的结果对比见图10所示。把ANSYS求解的结果作为频率真值，简化系统的频率误差见图4所示。由此可见简化模型的结果是相当精确的，已经完全可以满足工程计算需要。 

本发明基于变形修正的动力减缩算法，并且通过具体算例的应用得到以下主要结论： 

(1)本发明基于变形修正的动力减缩算法，对于一些复杂结构动力模型进行简化时，计算过程简单，而且简化模型和原模型之间始终保持着直接简单的联系。 

(2)本发明提出的方法对单元类型和形函数没有限制，只需要知道系统的一些基本信息，如：节点坐标、柔度矩阵等，便可对模型进行简化，操作非常简单、步骤清晰。 

(3)本发明提出的模型简化方法在保留模型边界节点的基础上把内部节点的自由度全部用各个同步性区域特征点的自由度来替代，当模型很庞大的时候可以对原模型进行超大规模降阶，计算效率大大提高。 

(4)本发明方法不同于其它简化方法，根据运动同步性的假设将结构的有限元模型划分为若干个同步性区域，采用适当的位移模式去逼近各个同步性区域的位移，得到整体位移与局部位移模式间的转换关系。其中一种比较简单而有效的位移模式就是刚体位移模九用总体应变能最小的算法来优化整体位移与局部位移模式间的转换关系，使得整体位移与位移真值最逼近。其实质是一种模型修正的概念。通过具体算例的分析，证明本文提出的方法可以很好的应用到工程计算中。 

超单元构造方法： 

首先得到超单元区域简化模型位移向量Q与相应的有限元结构节点位移向量U的转换关系，即：用基于总体应变能最小的方法得到转换矩阵 

。它们之间的关系如下： 

$U=\left\{\begin{array}{c}{u}_I\\ u_B\end{array}\right\}=\widetilde{\mathrm{\Λ}}Q=\left[\begin{array}{cc}{\widetilde{\mathrm{\Λ}}}_{\mathrm{II}}& {\widetilde{\mathrm{\Λ}}}_{\mathrm{IB}}\\ {\widetilde{\mathrm{\Λ}}}_{\mathrm{BI}}& {\widetilde{\mathrm{\Λ}}}_{\mathrm{BB}}\end{array}\right]\left\{\begin{array}{c}{q}_I\\ q_B\end{array}\right\}---\left(16\right)$

其中，u_I为有限元结构内部节点位移向量，u_B为有限元结构边界节点位移向量，q_I为超单元区域简化模型内部节点位移向量，q_B为超单元区域简化模型边界节点位移向量。 

式(16)可以写成如下两个式子： 

$u_B={\widetilde{\mathrm{\Λ}}}_{\mathrm{BB}}{q}_B+{\widetilde{\mathrm{\Λ}}}_{\mathrm{BI}}{q}_I---\left(17\right)$

$u_I={\widetilde{\mathrm{\Λ}}}_{\mathrm{IB}}{q}_B+{\widetilde{\mathrm{\Λ}}}_{\mathrm{II}}{q}_I---\left(18\right)$

由式(17)可求得： 

$q_B={{\widetilde{\mathrm{\Λ}}}_{\mathrm{BB}}}^{-1}(u_B-{\widetilde{\mathrm{\Λ}}}_{\mathrm{BI}}{q}_I)---\left(19\right)$

将式(19)代入式(18)中得： 

$u_I={\widetilde{\mathrm{\Λ}}}_{\mathrm{IB}}{{\widetilde{\mathrm{\Λ}}}_{\mathrm{BB}}}^{-1}(u_B-{\widetilde{\mathrm{\Λ}}}_{\mathrm{BI}}{q}_I)+{\widetilde{\mathrm{\Λ}}}_{\mathrm{II}}{q}_I$

$={\widetilde{\mathrm{\Λ}}}_{\mathrm{IB}}{{\widetilde{\mathrm{\Λ}}}_{\mathrm{BB}}}^{-1}{u}_B+({\widetilde{\mathrm{\Λ}}}_{\mathrm{II}}-{\widetilde{\mathrm{\Λ}}}_{\mathrm{IB}}{{\widetilde{\mathrm{\Λ}}}_{\mathrm{BB}}}^{-1}{\widetilde{\mathrm{\Λ}}}_{\mathrm{BI}})q_I---\left(20\right)$

为保证相邻超单元位移的协调性，引入： 

u_B＝u_B               (21) 

式(20)、(21)可合并为： 

$U=\left\{\begin{array}{c}{u}_I\\ u_B\end{array}\right\}=\left[\begin{array}{cc}{\widetilde{\mathrm{\Λ}}}_{\mathrm{II}}-{\widetilde{\mathrm{\Λ}}}_{\mathrm{IB}}{{\widetilde{\mathrm{\Λ}}}_{\mathrm{BB}}}^{-1}{\widetilde{\mathrm{\Λ}}}_{\mathrm{BI}}& {\widetilde{\mathrm{\Λ}}}_{\mathrm{IB}}{{\widetilde{\mathrm{\Λ}}}_{\mathrm{BB}}}^{-1}\\ 0& I\end{array}\right]\left\{\begin{array}{c}{q}_I\\ u_B\end{array}\right\}---\left(22\right)$

令： $\mathrm{\Δ}=\left[\begin{array}{cc}{\widetilde{\mathrm{\Λ}}}_{\mathrm{II}}-{\widetilde{\mathrm{\Λ}}}_{\mathrm{IB}}{{\widetilde{\mathrm{\Λ}}}_{\mathrm{BB}}}^{-1}{\widetilde{\mathrm{\Λ}}}_{\mathrm{BI}}& {\widetilde{\mathrm{\Λ}}}_{\mathrm{IB}}{{\widetilde{\mathrm{\Λ}}}_{\mathrm{BB}}}^{-1}\\ 0& I\end{array}\right],$ $\tilde{Q}=\left\{\begin{array}{c}{q}_I\\ u_B\end{array}\right\},$ 式(22)便可表示为： 

$U=\mathrm{\Δ}\tilde{Q}---\left(23\right)$

式(23)即为构造的超单元位移向量 

与相应有限元节点位移向量U之间的转换关系。 

一个超单元区域的一般动力学方程： 

$M_e\stackrel{\·\·}{u}+C_e\stackrel{\·}{u}+K_{e}u=F_e---\left(24\right)$

把式(23)代入式(24)中，并在等式两边同时左乘Δ^T，得： 

${\mathrm{\Δ}}^T{M}_e\mathrm{\Δ}\stackrel{\stackrel{\·\·}{~}}{Q}+{\mathrm{\Δ}}^T{C}_e\mathrm{\Δ}\stackrel{\·}{\tilde{Q}}+{\mathrm{\Δ}}^T{K}_e\mathrm{\Δ}\tilde{Q}={\mathrm{\Δ}}^T{F}_e---\left(25\right)$

令：M_s＝Δ^TM_eΔ，C_s＝Δ^TC_eΔ，K_s＝Δ^TK_eΔ，F_s＝Δ^TF_e，式(25)可以表示为： 

$M_s\stackrel{\·\·}{\tilde{Q}}+C_s\stackrel{\·}{\tilde{Q}}+K_s\tilde{Q}=F_s---\left(26\right)$

式(26)即为超单元的一般动力学方程，M_s为超单元的质量矩阵，C_s为超单元的阻尼矩阵，K_s为超单元的刚度矩阵，F_s为超单元的节点力向量。构造的超单元示意图如图5所示。当然超单元比较大的时候内部超节点也可取多个，超单元的几何形状也不仅限于矩形，与具体结构的几何形状有很大关系。 

得到超单元的单元矩阵后就可以根据总体结构超单元节点的总体编号，把各个超单元的单元矩阵集成为总的矩阵，即： 

$\stackrel{\‾}{K}=\underset{i}{\mathrm{\Σ}}{K_s}^i$ $\stackrel{\‾}{M}=\underset{i}{\mathrm{\Σ}}{M_s}^i$ $\stackrel{\‾}{C}=\underset{i}{\mathrm{\Σ}}{C_s}^i$ $\stackrel{\‾}{F}=\underset{i}{\mathrm{\Σ}}{F_s}^i---\left(27\right)$

超单元的实施例： 

问题还是前面的平面悬臂梁。该悬臂梁平均分成10个超单元，如图6所示。频率求解结果如图11所示，频率误差如图7所示。 

超单元简化方法： 

在超单元的基础上做进一步的简化，使求解方程的规模进一步降低，简化方法的出发点就是假定超单元的质量集中于超单元的内部超节点。 

把模型划分为几个超单元后，系统的无阻尼动力学方程为： 

$\stackrel{\‾}{M}\stackrel{\·\·}{\stackrel{\‾}{Q}}+\stackrel{\‾}{K}\stackrel{\‾}{Q}=\stackrel{\‾}{F}---\left(28\right)$

其中，Q为各个超单元位移向量集成后的超节点位移向量，M为各个超单元质量阵集成后的质量阵，K为各个超单元刚度阵集成后的刚度阵，F为各个超单元载荷向量集成后的载荷向量。 

把Q重新调整次序，写成(29)式的形式。质量矩阵和刚度矩阵也做相应的调整，以使原来的动力学方程保持不变。 

$\stackrel{\‾}{Q}=\left\{\begin{array}{c}{\stackrel{\‾}{Q}}_B\\ {\stackrel{\‾}{Q}}_I\end{array}\right\}---\left(29\right)$

其中，Q_B为边界超节点位移向量，Q_I为内部超节点位移向量。 

现在把各个超单元的质量都集中在超单元内部超节点处，刚度矩阵不发生变化。此时的系统动力学方程将变为(30)式所示。 

$\left[\begin{array}{cc}0& 0\\ 0& {\stackrel{\‾}{M}}_{\mathrm{II}}\end{array}\right]\left\{\begin{array}{c}{\stackrel{\·\·}{\stackrel{\‾}{Q}}}_B\\ {\stackrel{\·\·}{\stackrel{\‾}{Q}}}_I\end{array}\right\}+\left[\begin{array}{cc}{\stackrel{\‾}{K}}_{\mathrm{BB}}& {\stackrel{\‾}{K}}_{\mathrm{BI}}\\ {\stackrel{\‾}{K}}_{\mathrm{IB}}& {\stackrel{\‾}{K}}_{\mathrm{II}}\end{array}\right]\left\{\begin{array}{c}{\stackrel{\‾}{Q}}_B\\ {\stackrel{\‾}{Q}}_I\end{array}\right\}=\left\{\begin{array}{c}{\stackrel{\‾}{F}}_B\\ {\stackrel{\‾}{F}}_I\end{array}\right\}---\left(30\right)$

(30)式可拆分为两个等式，即： 

K_BBQ_B+K_BIQ_I＝F_B                (31) 

${\stackrel{\‾}{M}}_{\mathrm{II}}{\stackrel{\·\·}{\stackrel{\‾}{Q}}}_I+{\stackrel{\‾}{K}}_{\mathrm{IB}}{\stackrel{\‾}{Q}}_B+{\stackrel{\‾}{K}}_{\mathrm{II}}{\stackrel{\‾}{Q}}_I={\stackrel{\‾}{F}}_I---\left(32\right)$

由(31)式可得： 

Q_B＝K_BB ^-1(F_B-K_BIQ_I)            (33) 

将(33)式代入(32)式中，可得： 

${\stackrel{\‾}{M}}_{\mathrm{II}}{\stackrel{\·\·}{\stackrel{\‾}{Q}}}_I+{\stackrel{\‾}{K}}_{\mathrm{IB}}{{\stackrel{\‾}{K}}_{\mathrm{BB}}}^{-1}({\stackrel{\‾}{F}}_B-{\stackrel{\‾}{K}}_{\mathrm{BI}}{\stackrel{\‾}{Q}}_I)+{\stackrel{\‾}{K}}_{\mathrm{II}}{\stackrel{\‾}{Q}}_I={\stackrel{\‾}{F}}_I---\left(34\right)$

整理后可得： 

${\stackrel{\‾}{M}}_{\mathrm{II}}{\stackrel{\·\·}{\stackrel{\‾}{Q}}}_I+({\stackrel{\‾}{K}}_{\mathrm{II}}-{\stackrel{\‾}{K}}_{\mathrm{IB}}{{\stackrel{\‾}{K}}_{\mathrm{BB}}}^{-1}{\stackrel{\‾}{K}}_{\mathrm{BI}}){\stackrel{\‾}{Q}}_I={\stackrel{\‾}{F}}_I-{\stackrel{\‾}{K}}_{\mathrm{IB}}{{\stackrel{\‾}{K}}_{\mathrm{BB}}}^{-1}{\stackrel{\‾}{F}}_B---\left(35\right)$

令：M_l＝M_II，K_l＝K_II-K_IBK_BB ^-1K_BI，F_l＝F_I-K_IBK_BB ^-1F_B，则(35)式可写成： 

$M_l{\stackrel{\·\·}{\stackrel{\‾}{Q}}}_I+K_l{\stackrel{\‾}{Q}}_I=F_l---\left(36\right)$

(36)式即为基于超单元质量集中的动力学运动方程。 

超单元简化的实施例： 

问题还是前面的平面悬臂梁。各个超单元的质量都集中在内部超节点处，模型示意图如图8所示。频率求解结果如图12所示，频率误差如图9所示。 

Claims (1)

1.一种在模型简化基础上的超单元构造方法，具体为：

1)根据原模型的几何形状，把结构划分为若干部件，每一部件作为一个超单元，再根据各部件的运动特征划分为若干运动同步性区域，以得到简化模型，保留部件边界节点并且在部件各同步性区域特征点处增设超节点；

2)利用下述的基于变形修正的动力减缩算法进行简化模型优化：采用位移模式去逼近各个同步性区域的位移，得到整体位移u与局部位移模式q间的转换关系：

u＝Rq    (3)

用总体应变能最小的算法来优化整体位移u与局部位移模式q间的转换关系，使得整体位移u与位移真值

最逼近，

由关系式：

$\begin{array}{c}U\left(\mathrm{\δ\ϵ}\right)=\frac{1}{2}{\left(\mathrm{\δu}\right)}^{T}K\left(\mathrm{\δu}\right)\\ =\frac{1}{2}{(\tilde{u}-\mathrm{Rq})}^{T}K(\tilde{u}-\mathrm{Rq})\end{array}---\left(8\right)$

$\frac{\∂U}{\∂q}=0---\left(9\right)$

可得： $\stackrel{\‾}{u}=\mathrm{Tq}---\left(13\right)$

其中，

为位移真值的近似值，R表示总体位移变换矩阵，K表示系统总刚度矩阵，U(δε)为总体误差应变能，δε为应变误差，T＝K^-1((R^TR)^-1R^T)^TR^TKR；利用上述基于变形修正的动力减缩算法，得到超单元区域简化模型的位移向量Q与相应有限元结构节点位移向量U的转换矩阵

3)把区域边界节点位移关系回代到简化模型中保证了超单元边界位移的协调性，得到超单元位移向量

与相应有限元结构节点位移向量U的转换矩阵Δ；进而得到：超单元的质量矩阵M_s、超单元的阻尼矩阵C_s、超单元的刚度矩阵K_s、超单元的节点力向量F_s；

4)根据总体结构超单元节点的总体编号，把各个超单元的单元矩阵集成为总的矩阵。

Images