基于已知激励且同时考虑环境激励影响的贝叶斯模态识别方法与流程

本发明技术方案隶属于结构健康监测的模态识别领域，主要目的是解决已知激励模态识别中如何合理系统地考虑环境激励影响的模态识别问题。

背景技术：

针对于健康监测中的模态识别问题，主要是基于振动测试数据得到结构频率、阻尼比、振型等模态参数进行结构动力性能评估。虽然现有相关技术，但是随着越来越精确的测量技术，以及不断推进的高质量振动数据采集技术，都对模态识别方法提出更高的要求，以求提高结果的准确性，因此需要提出新型的模态识别方法来满足模态识别新形势的需要。另一方面，由于在动力测试整个过程中有一些环境的变化不受人为控制，而且不可能和理论推导所基于的环境状态完全一致，因此不确定性是不可避免的存在，并且不断上升，不确定性的来源很广泛，例如测量噪声，对准误差和建模误差。

因此，新方法的提出是形势的需要，从而更好的解决现存的问题。现有的方法主要有以下两个问题，第一个问题是针对已知外加激励，一般均忽略了环境激励的影响，然而由于结构质量较大，因此即使激励中已知激励为主导，但是从振动幅度的角度，环境激励的响应不容忽视，必将造成较大的精度影响，这个影响对高阶模态尤为显著。第二个问题是对于测试过程中所产生的不确定性没有理论系统的评估，传统技术只能直接利用识别得到的模态参数值本身，而对模态参数存在的误差及不确定性无法获得，从而也就不能利用其对精度的重要衡量价值。即使少数方法对不确定性进行了研究，基本也是通过同类型激励多次重复试验得到，费时费力，而且由于时间及测试环境约束，大多数情况下在实际测试中较难实现。

技术实现要素：

本发明从最开始的基本理论即全面考虑了环境激励的影响，得到了更符合工程实际的理论模型和识别方法，克服了传统方法对环境激励忽视的缺点；并且可以实现对动力参数识别结果精度的不确定性评估，得到的参数区间与传统的单一最优值相比，更充分的利用了所采集数据中包含的结构振动信息，从而更有效的利用大型结构所安装的健康监测系统。

本发明所开发的贝叶斯模态识别方法，最突出特征为将传统方法中为了理论简化而被省略处理的环境振动的影响，在构建理论模型时进行参数化，因为实际工程中这一部分的影响是不容忽视的，因此可以很好的去除传统方法中没有考虑的环境激励的影响；另一方面新方法同时兼顾了识别结果精确度的不确定性评估，从而可以得到更准确的模态识别结果，因此可以得到区间的形式代替传统的单一最优值，描述更符合自然规律，而且所得到的高精度结果以及不确定性参数可以进一步用于结构健康监测、损伤识别等领域。

为此，本发明给出的技术方案：

本发明为基于已知激励且同时考虑环境激励影响的贝叶斯模态识别方法，其特征在于，针对已知激励和环境激励共同激励情况下结构加速度数据进行分析，通过对目标函数的优化，得到结构的固有频率、阻尼比以及振型，并同时计算这些模态参数的不确定性。由于在构建理论模型的时候充分考虑了环境激励，因此从根本上解决了环境激励对结果精度的影响，从而对动力参数进行更准确的评估，是一种更加高效的方法，基于最终优化的程序计算时间只需几秒钟。

基于已知激励且同时考虑环境激励影响的贝叶斯模态识别方法，其特征在于，整个优化过程的流程为：

步骤(1)，对于优化程序所需要输入的结构参数集合θ的初始值，

自振频率f可以通过计算测得加速度的功率谱，并针对频率变量进行绘制，通过拾取曲线的峰值所对应的频率值，作为自振频率的初始值输入优化算法；

阻尼比ζ的初始值，通常会取具有普遍性的1％；

对于质量比r'的初始值，通过式(54)

计算得到；

振型φ的初始值等于式(48)

模态力功率谱s初始值为式(53)

预测误差功率谱se，应用式(51)求取初始值。

步骤(2)，输入初始值之后，接下来需要通过优化算法得到参数{f,ζ,s,se}的最优值，这一步将通过最小化负对数似然函数式(33)

步骤(3)，通过对式(40)进行特征值分解,得到振型的最优值φ，等于最大特征值所对应的特征向量的前一半，之后对其进行归一化。

步骤(4)，对于质量比r，可以先通过式(42)计算得到r'，并基于式(32)r'＝rφ(i)进一步得到r的最优值，从而完成了所有需要参数的最优值求解过程。

本发明是基于结构动力学的基本原理，并结合统计学中的贝叶斯理论，基于实测结构加速度数据和建立的理论模型，针对于结构模态参数(包括固有频率、阻尼比、振型等)，经过推导、构建后验概率密度函数。该函数包含测试中涉及的动力参数以及所测试的数据，从而可以从理论上推导得到目标函数中采集数据与模态参数之间的关系，以及识别模态参数的不确定性等。通过数值验证和实际场地测试验证，发现其可以较准确的识别结构模型参数，并在实际工程应用中具有较好的效果，所得到的结果可以为之后的结构损伤识别和健康监测服务。

本发明技术方案的关键理论包含以下四个部分，第一部分是构建通用的概率密度函数理论框架，同时涵盖已知激励和环境激励两种振动成分；之后为了解决优化的问题，第二部分将问题缩窄至分离模态并进行推导；由于优化过程中需要提供所求参数的初始值，初始值将大大决定之后优化的速度，因此本发明第三部分从理论上得到了各参数初始值的合理推导方式，即低噪音假设下的近似解；第四部分总结了所推导方法的整个流程和步骤，方便对方法的理解和使用。

附图说明

图1为本发明优化过程的流程图。

具体实施方式

以下结合附图对本发明技术方案做进一步介绍。

本发明技术方案的关键理论如下：

1.创建后验概率密度函数理论框架，同时考虑已知激励和环境激励

这部分理论包含两部分，从最开始的广泛的理论框架出发，根据具体解决的问题进行合理细化，从而得到更适应于实际工程应用的理论方法和算法。

1.1广泛贝叶斯理论框架

当结构被已知荷载激励时，将进行相应的振动。然而对于结构来说，由于其质量较大，因此虽然已知激励是主导激励，但是环境的干扰一直存在，且带来的响应不容忽视，考虑到这个因素所获得的响应应该包括已知激励响应和环境振动响应两部分，并同时用预测误差εj表示测量和理论模型响应之间的差异，因此测量响应可以定义如下：

这里j＝1,...,n，n为抽样点数。该方法为频域方法，将时域信号转化为频域信号采用的是快速傅里叶变换，加速度对应的频域信号k的定义为：

这里δt表示采样间隔；i²＝-1.将快速傅里叶变换的实部和虚部构建一个新的向量：

n表示测量的自由度数目。应当注意的是，在之后的优化过程中主要采用的优化数据为所选取目标模态主导的频带内的fft数据所构建的新向量{zk}。

基于贝叶斯理论所开发方法的优点之一是将一个较难的基于已知响应数据获得重要振动模态参数的问题，通过巧妙的转化，保证在不影响精度的情况下，转化成一个较易实现的反问题。基于这种转化，使得在已知测量加速度数据和输入力{fj}的情况下，得以构建结构参数θ的后验分布p(θ|{zk})表达式：

p(θ|{zk})∝p(θ)p({zk}|θ)(4)

其中，p(θ)表示结构参数的先验分布，由于在实际测试中，其往往假定为常数，因此后验分布p(θ|{zk})直接正比于似然函数p({zk}|θ)。

下一步将推导似然函数的表达式。由于在理论模型中，环境振动作为显式表示，等于在这个表达式中，环境激励的响应和预测误差将一起考虑，这两个变量的统计学特性可以分别通过变量s(模态力功率谱密度)和se(预测误差功率谱密度)来表示。对于来说，当采集数据足够时，针对于给定的θ是一个零均值的变量，且其fft分量相互独立并满足高斯分布。对于{zk}，其分布在不同的频率成分也是独立的，并且满足联合高斯分布，对于分布中涉及的均值和协方差是需要在数据分析过程中识别的两组主要数据。zk的均值为μk(θ)∈r²ⁿ，主要由已知激励振动响应主导，并且其fft分量满足μk(θ)＝e[zk|θ]，这里e[·|θ]表示对于参数集θ的条件期望。对于协方差矩阵ck(θ)∈r^2n×2n，因为已知激励的响应可以完全由模态参数确定，因此协方差矩阵仅与环境激励和预测误差有关。基于以上分析，可以得到似然函数的表达式为

这里，det(.)表示行列式的值；ck(θ)是zk的协方差矩阵，表达式为

公式中φ＝[φ(1),φ(2),...,φ(m)]∈r^n×m表示振型矩阵；m为振型数目；se表示预测误差的功率谱密度(psd)；i2n∈r^2n×2n是单位矩阵；hk∈c^m×m是传递函数，其中第(i,j)个元素为:

式中

hjk同理。对应第i个元素，对应的频率为fk；βik＝fi/fk为频率比；fk是fft对应的频率；fi和ζi分别表示第i个频率和阻尼比；sij是环境激励模态力功率谱矩阵s的第(i,j)个元素。

1.2已知单激励动力学理论

在常用的激励形式中，由于采用多激振器存在同步等难题，大多数情况都会使用单激振器进行激励。下面我们将针对单激振器的情况推导似然函数的表达式。我们需要首先推导均值μk(θ)的表达式。针对于已知激励受迫振动的加速度响应表达式为

这里φ(i)∈rⁿ(i＝1,...,m)是第i阶振型，包含n个测试自由度；和ηfi(t)分别是受迫振动的模态加速度、速度和位移，满足以下公式：

其中ωi＝2πfi是第i阶圆频率(单位为弧度/秒)；pi(t)是第i阶模态力，表达式为

这里φ(i,i)表示在第i阶模态中激励位置i所对应的振型分量；f(t)是受迫振动的激振力；mi表示第i阶模态的模态质量。将式(11)中pi(t)的表达式带入式(10)，并进行相应的fft变换，可以得到的表达式：

这里f(fk)是激振力f(t)的fft值，对应频率成分为fk.在这个等式的推导过程中，我们假设以下两个等式成立:

这两个假设较容易满足，只需要在数据采集中保持初始和结束时的响应幅值较小即可。实际操作中，可以通过在采集受迫振动响应时，分别在开始之前以及结束之后采集一段环境激励的数据即可。基于式(12)，加速度响应的频域变量表达式可以转换为

φ(i)为第i阶模态的振型。对于式(1)取fft,从时域变换到频域，表达式变为

这里表示所采集的加速度响应的fft；表示所对应的理论fft值；表示环境激励响应的fft；εk∈cⁿ表示预测误差的fft，预测误差假设在所选取的频率段内满足平滑的功率谱密度谱，由于所选取的频率段通常比较窄，所以这个假设很容易满足。预测误差可以表示为

这里z1和z2∈rⁿ分别表示两组相互独立的实数高斯分布向量。

对于一个已知激振器激励的情况来说，激振力f(t)可以用下式进行表示：

f(t)＝-mss(t)(18)

其中，ms和s(t)分别表示激振器的质量和激励器质量上所采集的激励的加速度时程曲线。因为激振器质量块运动方向和结构被激励产生响应的加速度方向相反，需要添加负号。

对式(18)取fft，并将其带入式(15)，可以得到:

式中sk是激振器加速度s(t)的fft值；激振器质量与结构模态质量的比值用ri表示:

这里需要注意的是ri的值与振型的模有关，因为振型以隐式形式包含在mi中。

从式(19)的定义中，我们知道所需要识别的模态参数包含每一个模态(以第i阶模态为例)的自振频率fi，阻尼比ζi，质量比ri，以及振型，其中所测量自由度数目为n。另一组参数包括模态力功率谱s和预测误差功率谱se。如果所选取频率段包含模态数目为m，那模态力的功率谱矩阵s将包含m×m个实数。对于最后两组功率谱参数我们均假设其在所选取的频率段内满足常数幅值。对于这些识别参数，其最优值将满足似然函数最大的要求。通常实际操作中我们会采用负对数似然函数(nllf)进行识别计算，表达式如下：

nf为频率段内采样点的数目。在这个等式中，振型的模将在稍后进行确定。在识别过程中如果直接进行最小值优化是不合理且难以实现的，因为ck(θ)矩阵的欠秩现象，也就是说矩阵的维度2n通常比其秩2m大，这个问题将导致优化过程中的病态条件问题。直接优化所面临的另一个难题是同时识别的参数数目较大，需要同时识别的参数数目为：

np＝3m+mn+m²+1(22)

由于测量自由度的数目对识别总参数的数目影响很大，例如当测试自由度数目n＝24时,即使选取频率段内的模态数目为m＝2，需要识别的参数总数目也将高达np＝59，这些问题都造成了直接优化的困难和问题，也可能导致不收敛。

为了解决优化的问题下面将集中解决分离模态的情况，即所识别区间内包含单一模态的快速算法推导。

2分离模态快速算法推导

对于分离模态来说，即所选取的频率段内仅有一个模态，也就是m＝1。针对于这种情况，我们所需要识别的参数θ包括频率f，阻尼比ζ，质量比r，振型φ∈rⁿ，模态力功率谱s以及预测误差功率谱se。

2.1负对数似然函数nllf表达式重构

为了推导振型φ∈rⁿ，质量比r以及预测误差功率谱se的解析解，需要将式(21)中的nllf表达式进行重组变换，从而使得这些变量可以在求解过程中得到较清晰的表达式。由于环境激励的作用，ck(θ)∈r^2n×2n和μk(θ)∈r²ⁿ将受到振型模数的影响，求解过程中需要将隐式的形式转化为显式，这里采用的是特征空间分解的方式推导得到ck(θ)的行列式值和逆:

式中

其中是一个复传递函数，βk＝f/fk；in∈r^n×n是单位矩阵。将式(23)(24)带入(21),nllf的表达式可变换为

其中

式(26)对于自振频率和阻尼比是非线性的，使得直接求解非常困难，因此对于这两个参数将通过直接的优化迭代程序来进行最优值求取。对于振型来说，在激励位置所对应的振型元素φ(i)为四次幂，振型中其他元素为二次幂，这些特征均表明直接得到振型的解析解是很难实现的，而且振型模数的确定是在任何类型的算法中都需要考虑的一个因素。在本理论中，将采用单位模数，即

φ^tφ＝1(31)

由于质量比r和激励位置的振型元素φ(i)总是同时出现，因此定义一个新变量：

r'＝rφ(i)(32)

基于这个新变量，nllf表达式变为

式中

通过合并参数φ(i)和r，振型的模数约束被合理而巧妙的去除，从而得到一个新的无约束的变量r'，之后对于φ(i)和r，可以通过所识别的φ和r'值计算得到。通过这样的方式，nllf表达式对于振型φ来说只是二次幂，因此可以通过解析求解得到。

2.2振型φ的解析解求解过程

在最终的目标函数中加入振型模数φ^tφ＝1的约束，并引入λ作为拉格朗日乘数，得到：

j(θ)＝l(θ)+λ(φ^tφ-1)(35)

使得目标函数j对于振型φ的一阶导数为零，得到：

a'φ+q＝seλφ(36)

式中

将振型进行变换处理之后，约束φ^tφ＝1和式(36)组成了一个约束特征值问题，而且变为了一个常规的特征值求解问题。因此振型的解析解可以通过解特征值问题得到，对应于最大特征值的特征向量前半部分，所对应的特征值矩阵为

式中g∈r^2n×2n，得到的振型需要将模数归一化。

2.3质量比r'的解析解

质量比r'的最优值可以通过对于l(θ)的求导并让其导数为零得到：

求解得到的最优值为:

3低噪音假设下的近似解

低噪音假设，即假设se/sdk＜＜1，通过这个假设上述等式可以被大大的简化，从而可以得到下列重要参数的近似解，虽然这个假设会对于精度有些影响，但是对于仅提供合理的求解参数初始值的需要，精度是完全足够的：

ak～1-se/sdk(43)

a'～a0'(44)

其中

基于上面等式的简化，等式(39)可以简化为:

通过这个简化以及q的简化表达式，矩阵g在式(40)中是的，因此可以看做是对角矩阵：

基于上一个章节的分析，可以知道近似最优振型即为a0'的最大特征值所对应的特征向量

同时也将低噪音假设及相关参数的简化表达式(43)和(44)带入负对数似然函数表达式(26)，这使得se和s的近似最优值可以得到解析解，过程如下：

式中

基于数学知识可以知道，满足alnx+b/x这个形式的唯一最小值在x＝b/a处，因此得到se的最优值为：

其中

对于s的求解过程类似，其最优值为：

对于质量比r,初始值可以通过对l(θ)微分获得：

通过这个合理的假设，参数的等式可以得到合理的简化，而且对于振型的近似最优值仅与采集数据有关。这意味着一旦采集数据得到之后，可以直接计算得到振型的近似最优值。而对于其他参数，包括se,s和r,相对应的表达式仅涉及两个参数变量，即自振频率和阻尼比。

4、快速算法流程和基本步骤

为了加快优化过程，对于振型和质量比的最优值解析解可以用其他参数的表达式来计算，因此在直接优化过程中所涉及的参数数目被大大减少了，从而使得优化过程可以具有广泛应用性且十分快速，不受采集自由度数目的影响。快速算法适用于模态识别中较普遍的分离模态的情况，从而可以得到更多参数的解析解表达式。在所选取的频率区间内，所涉及的参数主要包括自振频率f，阻尼比ζ，振型φ，模态力功率谱密度s，预测误差功率谱密度se和质量比r。

如图1所示，整个优化过程的流程为：

步骤(1)，对于优化程序所需要输入的结构参数集合θ的初始值，

自振频率f可以通过计算测得加速度的功率谱，并针对频率变量进行绘制，通过拾取曲线的峰值所对应的频率值，作为自振频率的初始值输入优化算法；

阻尼比ζ的初始值，通常会取具有普遍性的1％；

对于质量比r'的初始值，通过式(54)

计算得到；

振型φ的初始值等于式(48)

模态力功率谱s初始值为式(53)

预测误差功率谱se，应用式(51)求取初始值。

步骤(2)，输入初始值之后，接下来需要通过优化算法得到参数{f,ζ,s,se}的最优值，这一步将通过最小化负对数似然函数式(33)

步骤(3)，通过对式(40)进行特征值分解,得到振型的最优值φ，等于最大特征值所对应的特征向量的前一半，之后对其进行归一化。

步骤(4)，对于质量比r，可以先通过式(42)计算得到r'，并基于式(32)r'＝rφ(i)进一步得到r的最优值，从而完成了所有需要参数的最优值求解过程。

本发明的关键点和保护点：

1)目标函数的推导过程

2)基于已知激励和环境激励共同考虑的模态识别方法和理念

3)基于模态参数最优值和其不确定性的目标函数

4)整个方法的实现过程

与现有技术相比，本发明方法主要有以下几个优点：

1)本发明的技术比传统方法更加便捷，可以用极短的时间完成参数及结果的识别，因此可以在工程场地直接判定测试数据的有效性，如果测试数据有问题，可以第一时间补测。

2)本发明除了传统方法可以得到的最优值，同时可以基于该技术对最优值对应的不确定性进行评估，从而用一个合理区间表达结构参数的真值，识别的模型参数比传统方法更加准确。

3)本发明与传统技术相比，从根本上考虑了环境激励对于已知受迫振动结果的影响，从而避免了环境激励对结果精度的影响，使得到的结果更加精确，更好的符合工程的实际情况。