随机子空间模态参数识别精度的参数影响研究★

范倪圣，陈稳杰，周振纲

(嘉兴学院建筑工程学院，浙江 嘉兴 314001)

随机子空间方法[1-2]只需要根据测点响应数据即可识别结构工作模态参数，但该方法进行模态识别时可能存在模态遗漏和出现虚假模态问题[3-5]，对于此，国内外学者通过采用稳定图方法或模糊聚类方法等对其进行了改进，在一定程度上避免了出现虚假模态问题，较好地解决了随机子空间方法系统阶数n取值问题。然而，该方法的识别精度实际还会受到方法中相关参数(如Hankel矩阵行数i、测点数L，识别用到的测点响应数据量s)取值的影响。国内个别学者就Hankel矩阵行数i对随机子空间方法模态参数识别精度问题做了初步研究[6-7]。但同时考虑参数i和s不同取值组合对随机子空间方法模态参数识别精度的影响规律研究还尚未见相关报道。因此，为了将随机子空间方法更好地应用于环境激励下实际工程结构的模态参数识别中，有必要对参数s和i不同取值组合对随机子空间方法模态参数识别精度的影响规律问题进行研究。

1 随机子空间方法模态参数识别方法

采用基于协方差的随机子空间方法进行结构模态参数识别时，首先根据各测点的加速度响应测试数据构造如式(1)所示的Hankel矩阵。

(1)

其中,H为Hankel矩阵；Hp为Hankel矩阵过去分块矩阵；Hf为Hankel矩阵未来分块矩阵；yk(k=0,…,2i+j-2)为根据测点响应数据构造的列向量，上述构造的列向量的行数为测点数L；i和j分别为Hankel矩阵分块矩阵的行数和列数，Hankel矩阵的分块矩阵的行数和列数应满足关系式s=2i+j-1，s为计算时选取的测点响应数据数量。

根据生成的Hankel矩阵过去分块矩阵和未来分块矩阵按式(2)计算得到Toeplitz矩阵T：

(2)

对T矩阵按式(3)进行奇异值分解：

(3)

其中,U,S,V均为T矩阵经奇异值分解得到的矩阵，U和V为正交矩阵，S为对角矩阵；S1为S矩阵的分块矩阵，其为对角矩阵，对角元素T为矩阵异值分解得到的正奇异值，S1∈Rn×n；U1和U2为U矩阵的分块矩阵，U1∈RLi×n，U2∈RLi×(Li-n)；V1和V2为V矩阵的分块矩阵，V1∈RLi×n，V2∈RLi×(Li-n)。

系统的观测矩阵Γ及其逆矩阵Γ-1可表示为：

(4)

(5)

其中，Γ∈RLi×n；Γ-1∈Rn×Li。

根据式(6)，式(7)计算得到矩阵P1和P2。

P1=Γ-1(1∶n,1∶L(i-1))

(6)

P2=Γ(L+1∶Li,1∶n)

(7)

其中，P1∈Rn×(Li-L+1)；P2∈R(Li-L)×n；Γ-1(1∶n,1∶L(i-1))为取Γ-1矩阵第1行至第n行和第1列至第L(i-1)列范围内元素的构成的矩阵；Γ(L+1∶Li,1∶n)为取Γ矩阵第L+1行至第Li行和第1列至第n列范围内元素的构成的矩阵。

系统的状态矩阵A和输出矩阵Q可表示为：

A=P1P2

(8)

Q=Γ(1∶L,1∶n)

(9)

其中，A∈Rn×n;Q∈RL×n;Γ(1∶L,1∶n)为取Γ矩阵第1行至第L行和第1列至第n列范围内元素的构成的矩阵。

对计算得到的状态矩阵A按式(10)进行特征值分解：

A=EDE-1

(10)

其中，E为特征向量构成的矩阵，E∈Rn×n；D为特征值构成的对角矩阵，D∈Rn×n。矩阵E和矩阵D的元素都为复数。

提取矩阵D的对角线元素构成列向量Z，即：

Z={z1z2…zn}T

(11)

其中，zk为Z向量的第k行元素值，其大小等于D矩阵第k行的对角元素值；Z∈Rn×1。

根据式(12)计算得到识别模态频率列向量Ft。

Ft=abs(ln(Z)·fs)/2π

(12)

其中，fs为测点响应测量的采样频率，Hz；π为圆周率；ln(·)函数为对复向量元素进行自然对数运算；abs(·)函数为对复向量元素进行取模运算；Ft∈Rn×1。

根据输出矩阵Q和状态矩阵A特征值分解得到的特征向量构成的矩阵E按式(13)计算复矩阵G。

G=QE

(13)

按式(14)计算得到识别模态振型向量构成的振型矩阵Mt。

Mt=abs(G)*sign(180·(|angle(G)|-π/2)/π)

(14)

其中，angle(·)函数为对复矩阵元素进行取相位运算；abs(·)函数为对复矩阵元素进行取模运算；sign(·)函数为符号函数，当因变量分别为正值、零和负值时，函数值分别取1，0和-1。

2 随机子空间方法模态参数识别精度评价方法

本文主要就参数i和s取值对基于协方差的随机子空间方法的模态频率和模态振型的识别精度影响规律进行研究。以下对本文模态频率和模态振型识别精度计算方法进行简介。

2.1 模态频率识别精度评价方法

前述随机子空间方法识别得到的结构的各阶模态频率构成的列向量Ft可表示为：

(15)

其中，ftk为根据试验测量数据由随机子空间方法识别得到的第k阶模态的模态频率，k=1,2,…,m，m为识别模态数。

与随机子空间方法识别得到的各阶模态对应的结构真实的模态频率构成的列向量Fa可表示为：

(16)

其中，Fa∈Rn×1；fak为与随机子空间方法识别得到的第k阶模态对应的结构真实的模态频率，k=1,2,…,m，m为识别模态数。

根据式(17)可计算得到随机子空间方法识别得到的各阶模态的模态频率识别精度。

(17)

其中，eFk为与随机子空间方法识别得到的第k阶模态的模态频率识别精度。

参数i和s取值对随机子空间方法不同阶模态的模态频率识别精度的影响有一定差异，故本文采用按式(18)计算得到的模态频率平均识别精度指标EM来对随机子空间方法模态频率识别精度进行评价。EM值越大则随机子空间法的模态频率整体识别精度越高。

(18)

其中，m为识别模态数。

2.2 模态振型识别精度评价方法

前述随机子空间方法识别得到的结构的各阶模态振型构成的振型矩阵Mt可表示为：

(19)

其中，φtk为根据试验测量数据由随机子空间方法识别得到的第k阶模态的模态振型列向量，k=1,2,…,m。

与随机子空间方法识别得到的各阶模态对应的结构真实的模态振型向量构成的振型矩阵Ma可表示为：

(20)

其中，φak为与随机子空间方法识别得到的第k阶模态对应的结构真实的模态振型列向量，k=1,2,…,m。

模态振型识别的精度可采用识别模态振型与真实模态振型的模态保证准则MAC值来评价。反映随机子空间方法各阶模态的模态振型识别精度的MAC值可按式(21)计算得到。

(21)

其中，MACk为随机子空间方法识别得到的第k阶模态模态振型的MAC值。

参数i和s取值对随机子空间方法不同阶模态的模态振型识别精度的影响有一定差异，故本文采用按式(22)计算得到的模态振型平均模态模态保证准则MACM值对随机子空间方法模态振型的识别精度进行评价。MACM值越大则随机子空间法的模态振型整体识别精度越高。

(22)

3 五质点的弹簧-质量-阻尼系统模态识别算例

五质点的弹簧-质量-阻尼系统算例模型如图1所示。算例模型的质量、刚度和阻尼参数如表1所示。

表1 算例模型参数

在该算例模型的节点1和节点2处施加竖向白噪声激振力f1(t)和f2(t)来模拟算例模型系统受到的环境激励，通过计算获得各节点竖向加速度响应数据作为本文随机子空间方法模态识别的测点响应数据。

根据各节点的动力平衡关系可建立算例模型受上述白噪声激励时的动力方程[8]，如式(23)所示。

(23)

(24)

(25)

(26)

(27)

(28)

(29)

(30)

在节点1和节点2处均输入相同的高斯白噪声激励力，高斯白噪声激振力时程及功率谱曲线如图2所示，激振力时长为10 s，步长为0.001 s。采用Newmark-β法用Matlab软件编程求解算例模型动力响应时程。限于篇幅，此处仅给出节点1加速度响应时程及功率谱曲线，如图3所示。

算例模型的模态频率和归一化模态振型的真实值(理论值)可通过对式(23)忽略阻尼项和无激振力时得到的无阻尼自由振动方程进行特征值和特征向量分析计算得到。算例模型的模态频率和模态振型真实值分别如表2，图4～图8所示。

表2 i和s不同取值时随机子空间法模态频率和模态振型识别平均精度初步比较

在测点数L=5和系统阶次n=12取值一定时，对i=38，s=6 542和i=50，s=2 838两种参数取值情况根据前述基于协方差的随机子空间模态参数识别方法对算例模型进行了模态参数识别计算，识别得到的各阶模态频率和模态振型如表2，图4～图8所示。同时，给出了两种参数取值情况时模态频率平均识别精度指标EM值和模态振型平均识别精度指标MACM值，如表2所示。

由表2，图4～图8可知：i=38，s=6 542时随机子空间方法识别得到的算例模型各阶模态频率和模态振型识别值与其真实值非常接近，EM和MACM都接近于1；而i=50，s=2 838时模态频率和模态振型的平均识别精度要低于i=38，s=6 542取值情况。由此可见，参数i和s取值对随机子空间方法模态参数识别精度有一定影响，故有必要对参数i和s取值对随机子空间方法模态参数识别精度的影响规律进行研究。

4 i和s取值对模态参数识别精度影响规律研究

为了进一步研究参数i和s取值对随机子空间方法模态参数识别精度的影响规律，对i和s不同取值组合时采用随机子空间方法对前述算例模型的模态频率和模态振型进行了识别计算，并计算得到了各取值组合时模态频率平均识别精度指标EM值和模态振型平均识别精度指标MACM值。参数i和s不同取值组合构造时参数i取值从30到158，增量为2，参数s取值从2 000到8 998，增量为2，共构造了227 500个i和s不同取值组合。识别计算中测点数L=5和系统阶次n=12保持不变。

图9给出了i和s不同取值组合时采用随机子空间法识别得到的算例模型的模态频率平均识别精度指标EM值变化情况。

图10给出了i和s不同取值组合时采用随机子空间法识别得到的算例模型的模态振型平均识别精度指标MACM值变化情况。

由图9可知：参数i和s不同取值组合时，随机子空间法模态频率的平均识别精度指标EM值均较大，EM值在0.971～0.998范围内变化，表明随机子空间法识别算例模型的模态频率具有很高的精度。参数i和s不同取值对随机子空间法模态频率的平均识别精度指标EM有一定影响，主要表现为：

1)当参数i在42～110取值范围内时，参数s取值对模态频率平均识别精度指标EM的影响相对较明显。

2)参数i取值对模态频率平均识别精度指标值EM的影响相对明显，参数i取值在36～42及110～158取值范围内时，模态频率平均识别精度更高。

由图10可知：参数i和s不同取值组合时，随机子空间法模态振型的平均识别精度指标MACM值变化较大，MACM值在0.283～0.999范围内变化，表明参数i和s取值对随机子空间法识别算例模型的模态振型的精度影响较大。参数i和s不同取值对随机子空间法模态振型的平均识别精度指标MACM的影响主要表现为：

1)当参数i在30～96和130～136取值范围内时，参数s取值对模态振型平均识别精度指标MACM的影响较明显。

2)参数i取值对模态振型平均识别精度指标MACM的影响明显，参数i取值在22～24,32～36,60～62,72～78,82～84,102～120,136～138,142～148及156～158取值范围内时，模态振型平均识别精度较高，MACM值基本在0.7及以上。

综上可知：

1)参数i和s取值对随机子空间法模态频率的平均识别精度指标EM有一定影响，但影响总体相对较小。

2)参数i和s取值对随机子空间法模态振型的平均识别精度指标MACM影响较大。

3)参数i和s在若干取值范围内取值组合时随机子空间法的模型频率及振型均具有较高的识别精度，如i取36～38,s取2 500～7 000时，模态频率的平均识别精度可达0.99以上，模态振型的平均识别精度可达0.8以上。

5 结论及展望

本文以一五质点的弹簧-质量-阻尼系统的模态参数识别算例为例，对参数i和s取值对基于协方差的随机子空间方法的模态频率及模态振型的识别精度影响进行了研究。主要结论如下：

1)参数i和s不同取值组合时随机子空间法模态频率的平均识别精度指标EM值在0.971～0.998范围内变化，变化幅度较小，表明参数i和s取值对随机子空间法模态频率的平均识别精度影响相对较小，但仍有一定影响，同时也表明随机子空间法识别模态频率具有很高的精度。

2)参数i和s不同取值组合时随机子空间法模态振型的平均识别精度指标MACM在0.283～0.999范围内变化，变化幅度较大，表明参数i和s取值对随机子空间法模态振型的平均识别精度影响较大。

3)参数i和s在若干取值范围内取值组合时随机子空间法的模型频率及振型均具有较高的识别精度。

本文对参数i和s不同取值组合时随机子空间法模态频率和模型振型的平均识别精度的影响进行了初步研究。为了更好地将随机子空间方法应用于环境激励下实际工程结构的模态参数识别中，后续将结合相关数学方法对参数i和s合理取值方法进行进一步的研究。