一 简介

飞机地面振动试验(GVT)之一为纯模态试验 ^[1] 。在试验过程中，若干个激振器固定在某阶模态的频率同时激励。 通过调节每个激振器的大小和方向，可激励出纯模态，直接得到模态振型。这套机理可借鉴到试验模态分析(EMA) 和运行状态模态分析(OMA)的参数识别，对具有共同分母的频响函数FRF或半谱HSP进行提纯^[2~4]，得到参数识别的提纯算法。

对单输入多输出SIMO或多输入单输出MISO模态试验，FRF函数矩阵的一行或一列已知，这些FRF具有相同的分母，决定了模态的频率和阻尼。FRF由多阶模态叠加而成。只要为每个FRF找到合适的权系数，计权相加后可得到只含一阶模态的新FRF，此即所谓的提纯过程。由此再识别得到模态的频率和阻尼。改变权系数，可得到其它阶的模态频率和阻尼。

对于多输入多输出MIMO模态试验，多行或多列FRF已知，把一行或一列FRF看成一组，为每组FRF找到一个合适的权系数，同测点不同组频响函数计权相加后可得到一组只含一阶模态的FRF，此即所谓的MIMO提纯过程。由这组提纯后的FRF，再识别得到模态的频率，阻尼和振型。改变权系数，可得到其它阶的提纯FRF，识别模态频率，阻尼和振型。

对于EMA的MIMO试验，通过多变量指示函数(MvMIF) ^[5,6]，可计算得到权系数，进行提纯。当MvMIF的数值小于0.1，可通过FRF虚部直接得到模态振型。

对于EMA和OMA的MIMO试验，通过复模态指示函数(CMIF) ^[2,3]，可计算得到权系数，进行提纯。

当模态密集时，需要适当增加MIMO试验的参考点数目。

二 理论推导

2.1 SIMO或MISO EMA模态试验

对于SIMO或MISO的EMA试验，FRF的一行或一列已知。

激励点和响应点之间的FRF可写成：

其中，i为激励点，j为响应点，r为模态阶数，n为频率分析区间内总的模态阶数，m_r为广义模态质量。

(1)式常写成

此处●*为复共扼， Ar为复留数，且

当Ar为纯虚数，Ar的实部为0，Ar*= -Ar 。由式(1)

Ψ_irΨ_jr为实数。此处开始，本小节模态振型假定为实数。

假定采用质量归一提取振型，式(1)可重写为：

所有FRF有共同的分母。在对某阶模态进行参数识别时，需要考虑其它阶模态的影响。

假定选中m个FRF通过权系数进行叠加，满足m≥n，合成得到新的FRF：

即：

假定采用质量归一提取振型，式(1)可重写为：

此处r是模态阶数，r ≠ k。只在一阶模态r=k，

为非零。归一化处理，令。

这样，由式 (7)，可得到只含一阶模态的频响函数。不用考虑其它阶模态的影响，通过新的FRF可识别得到精确的模态频率和阻尼。

当所有模态的频率和阻尼已知后，由原始的FRF可逐点提取模态振型。

为了得到只含单阶模态的FRF，需要先计算得到实系数α_j（j=1，2，…m）。

单阶模态加速度或位移的频响函数如图1所示：

本小节此后的频响函数类型都看作加速度或位移。对于速度类型的频响函数，实部和虚部需要交换。对于OMA的半谱，可看成速度类型。

在图1中，以模态频率为中心，选择一个频率区间。在此频率区间外，虚部和实部之比的绝对值为：

确保μ<0.01,频率范围是

当ξ_r<0.005，外部频率区间为ω<0.62ω_r，ω>1.62ω_r

当ξ_r<0.001，外部频率区间为ω<0.90ω_r，ω>1.11ω_r

在所有FRF中，选m个虚部绝对值较大的FRF或全体FRF。

假设在频率区间为FRF谱线条数为N，因此求实系数α_k的方程为

图1 单自由度FRF实部和虚部

简写为

在式 (13)中，●^T表示转置矩阵。HI为FRF虚部。元素HI_ij，当i不为0，为频率区间外的谱线序号；i=0表示为频率区间内谱峰的位置，允许有小的误差，不会对识别参数的精度造成影响。j为FRF的点号。式(12)中，右边矩阵的第一个元素为1，也可假设成其它的常数。

从式 (13)，得到最小二乘解为

即

对于SIMO或MISO试验，以上算法的明显缺点是不适合密集模态或大阻尼模态。

如两个靠近的模态频率位置大于2条谱线，即频率差大于2Δf，此处 Δf=SF/N，SF为采样频率，N为进行FRF分析的FFT点数。这两阶的模态参数可识别出来。

注意到图1，模态频率位置处频响函数的实部为0，因此在式(12)中增加一行，新的方程为

两阶位置靠近的模态可以分开。新增的一行是HR，FRF实部，位于FRF主峰值的位置，和第一行元素所在位置相同。

对于包含只含一阶模态的FRF，在由式(10)和(11)决定的频率区间内，需要识别三个参数ω_r，ξ_r，Φ_r。方程为

此处HR_k为H的实部，HI_k为H的虚部，k为谱线序号，ω_k为第k条谱线的ω，N为频率区间内谱线总数。

简化成

方程的最小二程解为

此处X=[X₁  X₂  X3]^T，因此

得到了模态频率，阻尼和参考点的振型。根据这些识别出的参数，可得到合成的理论FRF。合成的理论FRF和合成的纯FRF的吻合程度，反映了识别出的模态参数的精度。当吻合程度不好时，其可能的原因有：在选定的频带区间，模态不够纯，即相邻模态的影响没有完全消除；或模态阻尼很大而选定的频率区间宽度不够。

2.2 SIMO或 MISO OMA模态试验

在式 (16)中, 矢量α_m×l为复数，式(16) 改写为

2.3 MIMO试验通过MvMIF计算

对MIMO试验，多变量指示函数MvMIF可直接用来计算系数，从而构造出一组纯FRF。

已知的FRF矩阵为H_q×p，q为总的响应点，p为激励点个数。

实系数可看成激励力矢量H_q×p，‖F‖=1，一组新的FRF为（HF)_q×l。因此，问题和飞机GVT纯模态试验相同，可定义为寻找激励力矢量，使得 

在每条谱线上，定义矩阵

从式(22), 得到

因此

在每条谱线，对下面的实对称矩阵进行奇异值分解

此处S为从大到小排列的实对角阵，U为归一化的实矩阵，有UU^T=I。I为单位阵。

多变量指示函数MvMIF曲线可以画出，如图2所示。

有p个激励点，就有p条MvMIF曲线，每条曲线由对角阵S相同序号的元素构成。

由MvMIF曲线的最小值点，可找到模态的位置。实系数为对应此最小奇异值的矩阵U的相应的一列。

当系数知道后，可得到一组提纯后的新频响函数FRF (HF)_q×l 。当MvMIF的值小于0.1，模态振型可通过相应谱线位置FRF的虚部直接得到。

图2  MvMIF曲线

MvMIF最小值依赖于激励的位置和数量，当有4阶模态密集相邻时，至少需要3个激励点。只要激励点数量足够，而且位置分布适当，即使重根也能轻易解耦。

当一组提纯后的FRF得到后，式(17) 扩展为

此处在元素ⁱHR和ⁱHI中, i表示一组FRF中的第i个FRF，一组共有q个。

直接求解方程(28)很费时，从等式 (28)，可得到等式 (29) 

可先求出模态频率和阻尼。

然后，再由式 (27)，求出模态振型

对于质量归一，存在等式

此处Ψ_v可看成单输入多输出SIMO系统的虚拟激励点，满足

Ψ₁，Ψ₂，…，Ψ_p 为激励点对应的振型。

由式(31)，推导出 

由式(33)，可计算得到Ψ_v。第r阶模态的实际振型为

2.4 MIMO试验通过CMIF计算

对于OMA试验以及FRF相干较差的EMA试验，MvMIF曲线无法清楚指定模态所在的位置，而复模态指示函数CMIF却可以。对某一频率处的FRF矩阵进行H_q×p进行复奇异值分解

此处Σ 为正的实对角阵，元素按从大到小的顺序排列。U为复矩阵， UU^H=I_q×q，I为单位阵。V为复矩阵，VV^H=I_p×p， I为单位阵。

各组的加权系数由下面复数方程得到

写成实数的方式为

将不同频率处Σ相同秩的元素连成曲线，得到了复模态指示函数，如图3所示。在曲线的最大值点，对应一阶模态。绝大部分模态在第一条曲线上，U矩阵和V矩阵的第一列和第一条曲线对应。U矩阵的一列和模态振型成比例，V矩阵的一列和模态参与因子成比例。当模态比较密集或有重根时，模态也有可能在CMIF的第二条或第三条曲线上。

图3   CMIF曲线

         首先从CMIF曲线确定模态频率所在的位置和在第几条CMIF曲线上，并得到U矩阵中和谱线的秩相应的一列矢量u_q×l 。

β为待定的复系数。此处，振型Ψ_q×l为复数  

         假定从式(35) CMIF算法得到的模态振型是可靠的。一组复系数ν_i(i=1,2,…p）为矩阵V的一列ν_p×l 。有

此处σ为对角阵Σ的一个和CMIF曲线对应的元素，即 CMIF曲线的值。

         对于提纯后的FRF，有

此处ω₀为CMIF最大值点对应的频率。

         对于提纯后的FRF，只有一阶模态，由公式 (6)，可得

此处最前的p点振型和参考点位置相同。

         由式(40)，可得到下面等式

或改写为

此处

式(44) 和(45)可写成矩阵的形式

此处元素ⁱHR和ⁱHI, i 表示一组提纯后FRF的第i个FRF。

对于q个提纯后的FRF， N为半功率带宽内谱线的条数，用于构建等式(48)。通过最小二乘法求解式(48)。 再由(46)和(47)，得到模态频率和阻尼。

当由式 (48)，求出γ。再由式 (42)，求出β。由式(38)，得到模态振型。

三 工程实例

实例1：这是在北航实验室完成的一个带方孔板的模态试验，只测和板垂直的方向，边界条件为自由。将4个传感器布在四个角，板在面内方向用钢丝绳悬挂，对垂直方向不造成约束，如图4所示。用锤子敲击所有的网格点，网格划分如图5所示。得到四组FRF，多变量指示函数MvMIF曲线如图4所示。

图4  带方孔板的模态试验

图5  模态试验测点的网格划分

将提纯算法得到的模态分析结果和特征系统实现算法ERA得到的分析结果进行比较，如表1所示：得到的模态频率和阻尼的结果非常一致。图6为对两种算法振型进行比较的Cross MAC结果，所得振型也非常一致。

表1 两种算法的分析结果比较范例

图6  提纯算法对ERA的Cross MAC

在以上的分析结果中，第6阶和第7阶模态非常接近。

以上算例表明，提纯算法所得分析结果可靠，能识别密集模态。

实例2：这是一个超大阻尼隔振台的模态试验，台子由工信部十院进行设计。台子长40米，宽7米，高0.7米，质量为1000公斤。前3阶刚体模态阻尼很大，阻尼比在20%到40%。这三阶大阻尼模态耦合在一起，造成模态参数用普通的参数识别算法很难识别，尤其是很难得到协调的模态振型。

考虑到结构的对称性，MIMO试验时选择4个角的垂直方向进行激励，得到了4组FRF。振动台面均匀分成了10*4的网格，共有55个测点。台面下和台面上的振动认为是相同的，只测量垂直方向。传感器采用国家地震局力学所设计生产的941B，测量速度。频率分析区间为0到8Hz。

图7 模态试验测点的网格划分

表2为分析结果。

表2  分析结果

图8为4阶模态的振型。

图8  4阶模态振型

前3阶模态严重耦合在一起，集总显示的FRF如图9所示。

图9  集总FRF显示

图10是ERA参数识别算法得到的稳定图，图中曲线为CMIF复模态指示函数曲线。

图10  ERA算法得到的稳定图

图11为PolyLSCF[3]算法得到的稳定图。

图11  PolyLSCF算法得到的稳定图

图12为PolyIIR[4]算法得到的稳定图。

图12  PolyIIR算法得到的稳定图

由于耦合严重，三种算法的稳定图都不够清晰。ERA算法稍好一点，但得到的刚体模态仍有个别点振型不够协调。

为得到提纯后的FRF，各个激励点的系数如表3：

表3  构造提纯FRF的4个激励点的系数

四组FRF乘以表3中的系数，每次得到一组新的提纯后的FRF。图13为第1阶和第4阶提纯FRF集总显示，参数识别时可选不同的频率区间分别识别，得到的模态振型非常协调，如图8所示。图14为第2阶提纯FRF集总显示，图15为第3阶提纯FRF集总显示。通过式(46)求得的模态振型非常协调。

图13  第1阶和第4阶提纯FRF集总显示

图14  第2阶提纯FRF集总显示

图15  第3阶提纯FRF集总显示

四   结论

1. 在SIMO和MISO试验，提纯算法能识别小阻尼以及频率相差大于2Δf的模态参数。

2. 对MIMO试验, MvMIF理论可用来得到一组提纯的FRF。识别参数的精度依赖于MvMIF的最小值，由激励点个数和位置决定。当激励点数量足够且位置合适，即使重根的模态也容易识别。复模态指示函数CMIF的提纯算法，既可用于EMA试验也可用于OMA试验。

3. 提纯FRF得到后，参数识别过程因为只需考虑一阶模态，非常简单，可看作实时完成。MIMO试验的FRF提纯因可以利用MvMIF或CMIF的数据，也可看作实时完成。提纯算法是目前所有模态识别算法中，唯一对密集模态参数识别能达到实时级的算法。

4. 从工程实例可以看出，对于耦合严重或者超大阻尼的模态试验，提纯后的FRF已完成了模态解耦，能直观显示模态的频率和阻尼。

5. 新算法进一步发展了多变量指示MvMIF和复变量指示函数CMIF的理论。以前这两种指示函数只能用来指示模态的位置，提纯算法利用它们可进一步直接识别出所有的模态参数，即模态频率、模态阻尼和模态振型。MvMIF指示函数值的大小反映了激振器的数量是否足够以及位置是否适当。

参考文献：

[1]       J. M. Liu,Q.H.Lu,H.Q. Ying, The Best Force Design of Pure Modal Test Based Upon a Singular Value Decomposition Approach, Advanced Aerospace Applications, Volume 1,Conference Proceedings of the Society for Experimental Mechanics Series, 2011, Volume 4, 119-129, DOI: 10.1007/978-1-4419-9302-1_11

[2]       W.Heylen, S.Lammeans, P.Sas, Modal Analysis Theory and Testing, K. U. Leuven, 1998

[3]       Robert E. Coleman, Ranall J.Allemang, Experimental Structural Dynamics, 2004, Author House

[4]       J. M. Liu, S. W. Dong, M. Ying, S. Shen.Dynamic Parameters Identification Technology in Bridge Health Monitoring, Proceeding of the 14th ASIA PACIFIC VIBRATION CONFERENCE, 2011, Dec. HongKong

[5]       William R. the Multivariate Mode Indicator Function in Modal Analysis Proc, 3th IMAC 1985

[6]       Nash M. Use of the Multivariate Mode Indicator Function for Normal Mode Determination. Proc, 6th IMAC 1988