1.本发明属于土木工程结构健康监测数据分析技术领域，涉及时变结构存在交叉模态的模态识别方法，具体为一种基于非线性调频模态分解的时变交叉模态识别方法。
背景技术：
2.结构模态参数包括模态频率、振型和阻尼比等，可以反应结构动力特性。在环境与时间的共同作用下，结构模态参数会随之发生变化，所以准确识别结构的时变模态参数十分关键。传统的模态参数识别方法假定结构是不随时间变化的，然而，运营状态下的土木工程结构动力特性往往会随时间发生改变，进而使模态参数表现出时变特性，传统的时不变模态识别方法不再适用。
3.时频分解方法能够捕捉信号中的时变频率，因此适用于时变模态参数的时变。经验模态分解方法和集合经验模态分解方法是经典的时频分解方法，但该方法缺乏数学理论依据，存在模态混叠与端点效应等问题。除此之外，变分模态分解算法通过求解变分约束问题的最优解来实现模态的分离，但无法直接获取时变结构参数。非线性调频模态分解方法能够用于频率变化范围较大时的信号分解，但对于模态识别并未有相关应用。
4.实际结构的模态参数将随环境的变化或结构损伤的变化而产生较大的变化，在变化的过程中将不同模态之间可能存在时间交叉的现象，导致不同模态的混淆和无法分离。因此，研究一种时变交叉模态参数识别方法具有重要的工程意义。
技术实现要素：
5.根据上述提出的技术问题，提供一种基于非线性调频模态分解的时变交叉模态识别方法。本发明采用的技术手段如下：
6.一种基于非线性调频模态分解的时变交叉模态识别方法，包括如下步骤：
7.步骤1、对结构上各加速度传感器采集到的加速度响应进行短时傅里叶变换，得到时频分布；
8.步骤2、计算时频分布对时间的偏导数，进一步计算i个传感器位置获得的时频脊线对应的频率fi(t,f)，i＝1,2,...,n，其中n为传感器个数，计算所有传感器位置对应的时频脊线频率的平均频率并将平均频率对应的时频点位置标记为
9.步骤3、取第一个传感器位置的时频分布x1(t,f)作为参考，计算第i＝1,2,...,n个位置在第k个频率fk处的时间振型系数φi(t,fk)；
10.步骤4、提取时频点位置对应的时间振型系数并将其按传感器位置的顺序构造成振型向量：将t＝0时的时间记为t0，提取t0时刻对应的脊线频率其中m为初始中心频率的个数。计算各脊线频率的时间振型相关系数m
l
，选取m
l
≥1-e时对应的时频点作为确定
交叉模态的初始中心频率的时频点，e为误差值，其中初始中心频率确定为时间的函数，即
11.步骤5、初始中心频率为计算模态分解过程中需要的调制和解调算子和计算两个解调信号和的初始值；
12.步骤6、给出最小化目标函数约束条件为采用优化算法对目标函数中的u
l
、v
l
、f
l
(t)进行更新，将最终得到的f
l
(t)，l＝1,2,...,m作为时变交叉模态的瞬时频率估计。
13.进一步地，步骤2中，采用如下方式计算第i个传感器位置获得的时频脊线对应的频率：
[0014][0015]
其中，xi(t,f)表示第i个传感器位置加速度响应的时频分布，表示xi(t,f)对时间t的偏导数，表示对t求偏导，j表示单位复数；
[0016]
步骤3中，采用如下方式计算第i个位置在第k个频率fk处的时间振型系数：
[0017][0018]
其中，φ表示振型系数，re[]表示提取信号的实部，xi(t,fk)表示时间t、频率fk对应的时频分布系数；
[0019]
步骤4中，采用如下方式计算各脊线频率的时间振型相关系数：
[0020][0021]
其中，m
l
表示脊线频率的时间振型相关系数，||表示取绝对值，|| ||表示计算2-范数；
[0022]
步骤5中，采用以下方式计算调制和解调算子及解调信号的初始值：
[0023][0024][0025][0026][0027]
其中，和分别表示初始调制和解调算子，diag[]表示取对角元素，其中α为给
定的惩罚因子，ξ为二阶差分因子：为初始中心频率。
[0028]
本发明的有益效果：
[0029]
本发明通过利用加速度响应的时频分布和时间振型相关系数提取具有交叉模态的时变模态初始中心频率，将初始中心频率作为非线性调频模态分解的输入参数，准确实现具有交叉模态的时变瞬时频率识别。
具体实施方式
[0030]
为使本发明实施例的目的、技术方案和优点更加清楚，下面将结合实施例对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。
[0031]
本发明实施例公开了一种基于非线性调频模态分解的时变交叉模态识别方法，包括以下步骤：
[0032]
步骤1、将结构n个传感器位置的加速度响应向量表示为x＝[x1,x2,...,xn]，采用短时傅里叶变换计算xi的时频分布，其表达式为xi(t,f)，其中i＝1,2,...,n表示传感器位置，t表示时间，f表示频率；
[0033]
步骤2、计算时频分布xi(t,f)对时间的偏导数，表示为其中表示对t求偏导，进一步计算时频脊线对应的频率：其中fi表示第i个传感器位置获得的时频脊线频率，j表示单位复数。计算所有传感器位置对应的时频脊线频率的平均频率并将平均频率对应的时频点位置标记为
[0034]
步骤3、取第一个传感器位置的时频分布x1(t,f)作为参考，计算第i＝1,2,...,n个位置在第k个频率fk处的时间振型系数：其中φ表示振型系数，re[]表示提取信号的实部，xi(t,fk)表示时间t、频率fk对应的时频分布系数；
[0035]
步骤4、提取时频点位置对应的时间振型系数并将其按传感器位置的顺序构造成振型向量：将t＝0时的时间记为t0，提取t0时刻对应的脊线频率其中m为初始中心频率的个数。计算各脊线频率的时间振型相关系数：其中m
l
表示脊线频率
的时间振型相关系数，| |表示取绝对值，|| ||表示计算2-范数。选取m
l
≥1-e时对应的时频点作为确定交叉模态的初始中心频率的时频点，e为误差值，其中初始中心频率确定为时间的函数，即
[0036]
步骤5、初始中心频率为计算模态分解过程中需要的参数和其中和分别表示初始调制和解调算子，diag[]表示取对角元素。计算两个解调信号和的初始值：其中α为事先给定的惩罚因子，ξ为二阶差分因子：
[0037]
步骤6、给出最小化目标函数约束条件为其中其中分别表示第l阶模态分量的调制和解调算子，t0到t
n-1
表示将连续的时间轴t离散化为n个时间点，f
l
(t)表示第l阶模态的瞬时频率，u
l
、v
l
表示解调信号，ε为容许误差。采用优化算法对目标函数中的u
l
、v
l
、f
l
(t)进行更新，将最终得到的f
l
(t)，l＝1,2,...,m作为时变交叉模态的瞬时频率估计。
[0038]
实施例
[0039]
采用一座简支梁桥的模拟数据分析进行说明。该简支梁桥的矩形截面尺寸为宽40cm、高85cm，长为10m。四个竖向加速度传感器分别固定在在2m、4m、6m、8m处梁顶端的位置。对2m到4m间的局部梁端的弹性模量进行时变损伤模拟，将弹性模量由3
×
104mpa经过1400s线性递减到0，再恢复到3
×
104mpa，继续持续600s，产生加速度的总时间长度为2000s。激励形式为白噪声，采样频率为200hz。该模拟结构共有五阶模态。在模拟损伤的过程中，第三阶模态频率随时间变化形式为0-1400s时从60hz下降到20hz，1401s-2000s时维持在60hz。第四阶模态频率随时间变化形式为0-1400s时从100hz下降到30hz，1401s-2000s时维持在100hz。因此，第三阶模态频率在1401s-2000s时变化到了第四阶模态频率的范围内，发生了交叉，从而不易区分60hz这一频率属于第几阶模态。
[0040]
将结构4个传感器位置的加速度响应向量表示为x＝[x1,x2,x3,x4]，采用短时傅里叶变换计算xi的时频分布，其表达式为xi(t,f)，其中i＝1,2,3,4表示传感器位置，t表示时间，f表示频率；
[0041]
计算时频分布xi(t,f)对时间的偏导数，表示为其中表示对t求偏导，进一步计算时频脊线对应的频率：其中fi表示第i个传感器位置获得的时频脊线频率，j表示单位复数。计算所有传感器位置对应的时频脊线频率的平均频率并将平均频率对应的时频点位置标记为
[0042]
取第一个传感器位置的时频分布x1(t,f)作为参考，计算第i＝1,2,3,4个位置在第k个频率fk处的时间振型系数：其中φ表示振型系数，re[]表示提取信号的实部，xi(t,fk)表示时间t、频率fk对应的时频分布系数；
[0043]
提取时频点位置对应的时间振型系数并将其按传感器位置的顺序构造成振型向量：将t＝0时的时间记为t0，提取t0时刻对应的脊线频率其中计算各脊线频率的时间振型相关系数：其中m
l
表示脊线频率的时间振型相关系数，||表示取绝对值，||||表示计算2-范数。选取m
l
≥1-e时对应的时频点作为确定交叉模态的初始中心频率的时频点，取e＝0.001，其中初始中心频率确定为时间的函数，即
[0044]
初始中心频率为计算模态分解过程中需要的参数和其中和分别表示初始调制和解调算子，diag[]表示取对角元素。计算两个解调信号和的初始值：其中α为事先给定的惩罚因子，ξ为二阶差分因子：
[0045]
给出最小化目标函数约束条件为其中
分别表示第l阶模态分量的调制和解调算子，t0到t
n-1
表示将连续的时间轴t离散化为n个时间点，n＝200
×
2000，共4
×
105个数据点，f
l
(t)表示第l阶模态的瞬时频率，u
l
、v
l
表示解调信号。采用优化算法对目标函数中的u
l
、v
l
、f
l
(t)进行更新，将最终得到的f
l
(t)，l＝1,2,3,4,5作为时变交叉模态的瞬时频率估计。