基于环境激励的大跨度斜拉桥模态参数和索力识别

叶锡钧，颜全胜，李 健，王卫锋，朱添丰，刘明慧

(1.华南理工大学 土木与交通学院，广州 510641;2.伊利诺伊大学 香槟分校工学院，美国香槟 61801)

基于环境激励的大跨度斜拉桥模态参数和索力识别

叶锡钧1，颜全胜1，李 健2，王卫锋1，朱添丰1，刘明慧1

(1.华南理工大学 土木与交通学院，广州 510641;2.伊利诺伊大学 香槟分校工学院，美国香槟 61801)

模态参数和索力是评估斜拉桥健康状态的关键参数。建立了崖门大桥的有限元模型并对其进行了基于环境激励的模态测试和拉索振动测试;提出了基于ERA的多参考点稳定图算法，设置不同的参考点，利用自然激励技术结合特征系统实现算法(NExT－ERA)识别模态参数，通过阻尼比、基于输出矩阵的一致模态指标(CMI_O)和模态置信度(MAC)作为判别标准，识别出崖门大桥的竖向和横向模态参数，通过和增强频域分解法(EFDD)识别结果的比较，可知该算法具有良好的识别效果;分析了斜拉索与主梁的共振频率范围，通过二次拟合识别较长拉索的低阶频率，根据两种不同方法的索力识别结果可知，该桥的索力分布比较均匀对称。

环境激励;斜拉桥;模态参数识别;多参考点稳定图;索力识别

对于运营中的大跨度斜拉桥，基于环境激励的模态测试和频率法［1］分别是测定其动力特性和索力的最快捷方法，而且几乎都是唯一选择。

目前，模态参数识别的方法主要有峰值法(PP)、频域分解法(FDD)、增强频域分解法(EFDD)等频域方法和特征系统实现算法(ERA)、随机子空间算法(SSI)、时间序列分析法等时域方法。频域的识别方法简单、快速而且直观，但在识别阻尼比时需进行反傅里叶变换，在时域内通过指数衰减法来实现，受截断误差的影响，精度不高;而时域的识别方法，其共同的缺点是难以确定结构的系统阶次，难于区分真实模态与噪声模态［2－3］。近年来，稳定图方法［4－5］的应用越来越广，被认为是一种比较理想的确定系统阶次的方法，可以从诸多模态中判别真假模态，但如何设置判别标准是关键步骤。

对于索力识别，近年来国内外学者研究得很多。Zui等［6］提出了利用低阶频率估计索力的实用公式;王卫锋［7］以二次抛物线代替悬链线为基础，提出了考虑拉索垂度、抗弯刚度、边界条件和斜度的影响的拉索振动方程，通过迭代求解方程得到索力值。对于各种不同的方法，准确识别拉索的自振频率是关键。

本文建立了崖门大桥的有限元模型并对其进行了基于环境激励的模态测试和拉索振动测试;提出基于ERA的多参考点的稳定图算法，识别崖门大桥的模态参数;并利用二次拟合求出长索的低阶频率，根据Zui等［6］的经验公式和王卫锋［7］的方法识别索力，进行比较。

1 NExT－ERA

1.1 NExT

NExT(Natural Excitation Technique)方法是由James等［8］于1993年提出，其基本思想是利用白噪声环境激励下结构两点之间响应的互相关函数和脉冲响应函数具有相似表达式的特点，将其代替结构的自由振动响应或脉冲响应函数作为输入，运用时域模态识别方法进行模态参数识别。NExT法的计算有两种方法，一种是在时域里采用卷积算法直接计算得到，另一种是先计算出实测信号的互功率谱，再经傅里叶逆变换得到实测互相关函数［9］。本文采用后者，由于其采用了统计平均处理，即使会受到一定的谱泄漏的影响，但所得到的互相关函数的信噪比也有较大幅度的提高。

1.2 ERA

ERA(Eigensystem Realization Algorithm)是多输入多输出的时域模态参数识别方法，由 Juang等［10］于1984年提出，其利用实测的脉冲响应函数为基本模型构造广义Hankel矩阵，运用奇异值分解方法，确定系统的阶次和状态方程中的系统矩阵A、输入矩阵B和输出矩阵C，进而求解系统矩阵A的特征值问题，求得极点与留数，从而确定系统的模态参数，具体的算法公式详见文献［10］。

2 基于输出矩阵的一致模态指标

Pappa［11］提出了一致模态指标 CMI，其计算根据EMAC和MPC，公式如下：

其中：i表示第i阶振型;j表示第j个输出点，k表示第k个输入点;EMACjik表示第i阶振型下，对应j－k(输出－输入对)的EMAC值;和分别是基于输出矩阵C和输入矩阵B的EMAC值。

Yun［12］对CMI进行改进，提出基于输出矩阵C的模态一致性指标(CMI_O)，并通过数值算例验证了在基于环境激励的模态分析过程中，CMI_O比CMI更具可靠性：

式中：i表示第i阶振型;j表示第j个输出点;q表示输出点的数量。从式(3)可以看出，CMI_O的计算基于，仅取决于输出响应。

3 基于ERA的多参考点稳定图算法

对于环境激励的情况，由于系统的输入激励未知，通过NExT/ERA方法识别模态参数时，选择某输出点作为参考点并计算其与其他输出点的互相关函数矩阵作为输入，参考点非常重要，选择不当则不能识别出某些较弱的模态。基于Pappa等［13］应用的真实模态判别方法，本文提出一种基于ERA的多参考点稳定图算法，以基于输出矩阵的CMI_O代替CMI作为判别准则，本文选取CMI_O阀值为75%。主要包括图1的三个过程：

图1 多参考点稳定图算法的三个过程Fig.1 3 parts of the multiple reference DOFs stabilization diagram algorithm

算法步骤如下：

(1)选定初始参考点r，设定不同的系统阶次(如，n+10，n+20，…，max_n)，运用NExT/ERA识别出不同阶次下的模态参数;每个阶次下识别出的模态参数，若能通过阀值“0% ＜阻尼比＜10%和CMI_O＞75%”，则认为是可信度较高的模态参数，进入下一步判别，共有N个模态参数;

(2)通过设置Δf，将第(1)步中的N个模态归结为M组(每个组里各个模态相互间的频率差值小于Δf，可以认为是同一阶的频率)，不符合该要求的模态将被记录下来加到当参考点为r+1时重复这一步判别;

(3)判断M组模态中，每一组里模态的数目T，当T＞10时，认为该模态是真实模态(否则将被记录下来加到当参考点为r+1时重复这一步判别);同时寻找每一组中CMI_O最大值对应的模态作为本阶段结果，得到P阶模态;

(4)通过步骤(3)可以在稳定图上清晰看到稳定轴，计算P阶模态中任意两阶之间的MAC值，小于70%的则认为都是真实模态;若大于70%，对应CMI_O值较大者则认为是真实模态;得到参考点为r时的模态参数组Set(r)，在稳定图上以红点标记;

(5)参考点位置为r+1，依步骤(1)至(4)得到Set(r+1)。比较Set(r)和Set(r+1)同一阶模态的CMI_O，若 CMI_O(modei(r))＞ CMI_O(modei(r+1))，则把Set(r+1)中的 modei(r+1)用 Set(r)中的modei(r)代替(CMI_O值较大者认为更准确);若Set(r)中有Set(r+1)没有出现的模态，则加到Set(r+1)中，得到新的参数组NewSet(r+1);

(6)改变参考点位置为r+2，重复步骤(1)至(4)得到Set(r+2)，再根据步骤(5)比较NewSet(r+1)和Set(r+2)得到NewSet(r+2)，直到最后一个测点作为参考点，得到最终的模态参数。

算法流程图如图2所示。

图2 多参考点稳定图算法流程图Fig.2 Flowchart of the multiple reference DOFs stabilization diagram algorithm

4 实例分析

4.1 工程概况

崖门大桥是广东省西部沿海高速公路建设项目的重点工程之一，大桥为双塔单索面预应力混凝土斜拉桥，采用塔梁固结柔性墩体系。每侧边跨设置一个辅墩，跨径组合为(50 m+115 m+338 m+115 m+50 m)，全桥总长668 m。主梁是单箱五室预应力混凝土梁，梁高3.48 m，桥面全宽26.8 m，设双向四车道。主墩为双壁柔性墩，横向宽12 m，双壁之间中心距6 m，墩高47.6 m。桥面以上塔柱高73.5 m，塔柱断面为单箱混凝土断面。全桥共有斜拉索200根，主梁标准索距6 m，如图3(a)。

利用ANSYS软件，按照塔梁固结形式，采用Beam4和Link8单元建立大桥的空间杆系有限元模型，斜拉索的初始应变值按照设计索力计算给定，除采用Ernst公式修正斜拉索的弹性模量外，其他参数按照设计参数给定，如图4。

4.2 模态测试与参数识别

斜拉索的布置为单索面，主梁中间有宽3 m的分隔带，由于交通限制，只在主梁中间分隔带布置竖向和横向测点，见图3(b)。模态测试的仪器采用中国地震局工程力学研究所的99－1型拾振器、941型8线放大器和杭州忆恒科技的AVANT动态数据采集系统。

4.2.1 竖向模态参数识别

全桥共布置45个竖向测点，以跨中的22、23测点为参考点，分为11组，每组6个测点(第8组除外，表1所示)，设置采样频率为160 Hz，每组采样时间20 mins。

表1 竖向模态测试方案Tab.1 Plan for vertical modal test

从有限元计算结果得知，斜拉桥的前几阶固有频率都在10 Hz以内，故选取较为平稳的数据段，把采样频率降到40 Hz，对数据进行去除趋势项和平滑去噪处理后采用增强频域分解法(EFDD)和本文算法进行识别。限于篇幅，本节只列出竖向测试的第4组数据的模态参数识别过程。典型的加速度时程见图5。

通过EFDD的识别结果可清晰看出，在0～3 Hz范围内有8阶主要的模态，其中第4和第8阶模态较弱，其为ERA的时域识别方法提供很好的基础，见图6。

图5 竖向第4组：基准测点22、23的加速度时程Fig.5 Acceleration responses of Channel 22＆23 in group 4

对于NExT－ERA，选择的参考点不同，识别出的模态数目可能不一样。以竖向测试的第4组为例，根据本文算法，选定参考点为1通道，通过第(1)～(3)步，由稳定轴观察到7阶模态(图7，T＞10的共有7阶模态)，再根据第(4)步比较7阶模态里相互间的MAC值确定是否真实模态，得到模态参数Set(1)(6阶，第7阶与第6阶的MAC值大于70%，排除)，在稳定轴上标记红点。

同理，改变参考点为2通道，经过(1)～(4)步，得到模态参数Set(2)，第2和4阶模态没被识别出来(图8);接下来第(5)步，除了比较两组模态 Set(1)和Set(2)的CMI_O值外，还把Set(1)中的第2和4阶模态加到Set(2)中，得到Newset(2)。当参考点为通道3时识别的参数组为Set(3)，比较Newset(2)和Set(3)，得到Newset(3)，直到参考点变为通道6得到第4组的最终模态参数。

同时，对于某些比较弱的模态，往往要通过设置较高的系统阶次才能识别出来。如图7，当通道1为参考点时，直到系统阶次为240左右，第1、2阶模态才开始出现;而参考点为通道2时，图8，第1阶模态在阶次为10时就出现了，但第2阶模始终没被识别出来。

图6 竖向第4组：EFDD识别结果Fig.6 Singular value by EFDD：group4，vertical

图7 竖向第4组：参考点为通道1时的稳定图Fig.7 Stabilization Diagram of vertical group 4：Ref.=Channel 1#

图8 竖向第4组：参考点为通道2时的稳定图Fig.8 Stabilization Diagram of vertical group 4：Ref.=Channel 2#

观察第1～11组的模态参数，在每组里面都出现的模态认为是真实模态。对于EFDD，求11组模态参数的平均值为最终结果;对于本文算法，对比各组里同一阶模态的CMI_O值，最大值为最终识别结果，结果列于表2。

从表2看出，本文算法与EFDD识别的结果基本符合，但识别结果与有限元计算结果存在一些差别。

表2 竖向模态参数识别结果Tab.2 Identified vertical modal parameters

按最大位移值归一化方法，得到桥面竖向振型。由于测试的分组数目较多，识别的振型误差会有累计，但可以看出，本文算法识别的振型更接近有限元结果。限于篇幅，本文只列出前三阶竖向振型，见图9～图11。

4.2.2 横向模态参数识别

横向测试共布置22个测点，只在图3(b)中边跨偶数和中跨基数测点，以19和21测点为参考点，分为7组，每组5个测点(第7组4个测点除外)。

通过EFDD的奇异值分布图，在0～0.7 Hz范围能识别出4 阶模态(0.332，0.469，0.566，0.664)，图 12。但通过本文算法最终确定3阶真实模态，图13和表3。识别的前3阶横向振型见图14～图16。

表3 横向模态参数识别结果Tab.3 Identified transverse modal parameters

图14 第1阶横向振型Fig.14 1st Transverse bending mode

图15 第2阶横向振型Fig.15 2nd Transverse bending mode

图16 第3阶横向振型Fig.16 3rd Transverse bending mode

4.3 斜拉索振动测试与索力识别

索力测试采用Imote2无线传感器节点和SHM_A传感器板［14］，将其固定在离减震器2 m的位置，测量斜拉索在竖向平面的振动，图17。设置采样频率100 Hz，滤波频率30 Hz，采样时间为3 mins，对东塔和西塔中跨各5对不同长度的斜拉索(1、7、13、19、25，每对上、下游各1根，共20根)进行振动信号采集(见图3(a))。分析采用的降采样频率为50 Hz，NFFT=1 024，频率分辨率为0.005 Hz。

4.3.1 频率识别

在测试过程中发现，在无风或风荷载很小的时候，有些斜拉索也会产生大振幅的剧烈振动，同时伴有桥面的振动。将各索的自功率谱密度函数依次缩小一定的倍数和主梁的竖向振动频谱图画在同一张图上(图19)，可以发现，较长的斜拉索(13、19、25)在0 ～3 Hz范围内的频率与桥面的振动频率相同;可以认定上述拉索的大振幅振动是由桥面振动激发引起的受迫振动，其实测的频率是受迫振动的频率，而不是其自振频率，从频谱图上难以观察出其基频的位置，只能识别出其第3或4阶以后的频率(图中箭头);而在7号索以前的斜拉索，由于其基频较大(大于2 Hz)，故能从频谱图上清晰地识别出来。图19中，能通过峰值识别出1号索的自振基频为4.023 Hz，7号索的自振基频为2.031 Hz。

图19 东塔上游斜拉索索梁耦合振动的频谱图Fig.19 Power spectrum of Deck and Cables：East pylon

对于不能识别前几阶频率的较长的13、19、25号斜拉索，由于使用频率间距法确定拉索的低阶频率误差较大，本文采用二次拟合的方法计算其前几阶的频率。

图20(a)是东塔上游25号索的各阶频率拟合结果，从频谱图识别第4到20阶频率，以二次曲线进行拟合。实测的和拟合的各阶频率差值见图20(b)。可见，实测各阶频率的差值最小为0.703 Hz，最大为0.79 Hz，且各阶频率差值的分布无序;而拟合的各阶频率的差值变化是平缓的，从0.731 Hz到0.767 Hz，与文献［15］的结论较为吻合。

图20 东塔上游25号索的频率识别结果Fig.20 Identification result of Cable 25#(upstream)of East Pylon

4.3.2 索力识别

由于安装减震器后拉索的长度变短，频率变大，如按原索长计算索力必然造成很大的误差。故认为2002年竣工验收［16］时安装减震器前后的恒载索力不变，按照安装减震器前后的频率对计算索长进行修正。利用MATLAB根据Zui等［6］的实用公式编制索力计算程序，以100 m为标准将索初步判定为长索或短索;再以Γ和ξ作为依据进行迭代判断，该使用7个实用公式中的哪一个;长索用第4阶频率为计算，短索则用第1阶或第2阶频率为计算。王卫锋［7］的索力计算方法采用跟Zui的方法同样的频率阶数和频率值。这里应当指出，对于长索且满足Γ和ξ条件的，Zui的公式只用到第4阶频率，故上一步中拟合的低阶频率误差对该索力的识别没有影响。

图21 东、西塔实测索力结果Fig.21 Result of cable tension estimation of East and West Pylon

从图21中看出，索力的整体分布比较对称均匀，两种方法计算的索力值很接近，误差在5%之内。除了西塔1号索外，其余斜拉索的实测索力都比设计索力要大，这主要是由于施工过程中为保证主梁线型而作了索力调整，而且经过对比修正后的计算索长也不能保证十分精确。

5 结论

从模态识别结果看出，对于环境激励的模态分析，使用CMI_O代替CMI作为阀值，通过设置不同的参考点计算互相关函数作为输入，基于ERA的多参考点稳定图算法能很好地识别出大桥的低阶模态参数，有效地分辨虚假模态和避免真实模态的遗漏，准确地识别出主梁竖向和横向的前3阶振型。

但在分析过程中发现，有三个地方的取值不能是一个定值：第(3)步中判断模态出现的次数T时(本文取值T=5)，对于某些模态测试环境理想、激励充分的情况下T值能取大些;另在判别准则的取值上还是取决于经验，如第(1)步中设置阀值时CMI_O的取75%;第(4)步中MAC的取70%，都是根据经验取值的，而且在一定程度上取决于模态试验的效果。在这些问题上需要工作者能进行一些经验判断。

通过频率识别和索力识别结果，可知主梁和索在0～3 Hz范围内有共振现象;二次拟合法比频率间距法能更准确地识别出较长拉索的低阶频率;崖门大桥的中跨斜拉索索力分布比较对称均匀;但在如何精确地修正安装减震器后的索长问题上还得进一步研究。

［1］王卫锋，韩大建.斜拉桥的索力测试及其参数识别［J］.华南理工大学学报(自然科学版)，2001，29(1)：18－21.

［2］续秀忠，华宏星，陈兆能.基于环境激励的模态参数辨识方法综述［J］.振动与冲击，2002，21(3)：2－5.

［3］王 彤，张令弥.运行模态分析的频域空间域分解法及其应用［J］. 航空学报，2006，27(1)：62－66.

［4］Peeters B， De Roeck G.. Reference based stochastic subspace identification for output－ only modal analysis［J］.Mechanical Systems and Signal Processing，1999，13(6)：855－878.

［5］常 军，孙利民，张启伟.基于两阶段稳定图的随机子空间识别结构模态参数［J］.地震工程与工程振动，2008，28(3)：47－51.

［6］ Zui H，Shinke T，Namita Y.Practical formulas for estimation of cable tension by vibration method［J］.Journal of Structural Engineering，ASCE，1996，122(6)：651 －656.

［7］王卫锋.PC斜拉桥的施工监测［D］.广州：华南理工大学建筑学院，2000.

［8］James G H，Carne T G，Lauffer J P.The natural excitation technique(NExT)for modal parameter extraction from operating wind turbines［R］.Albuquerque：Sandia National Laboratory report，SAND92 －1666.UC －261，1993.

［9］王 济，胡 晓.MATLAB在振动信号处理中的应用［M］.北京：中国水利水电出版社，2006.

［10］ Juang J N，Pappa R S.An eigensystem realization algorithm for modal parameter identification and model reduction［J］.Journal of Guidance，Control and Dynamics，1985，8(5)：620－627.

［11］ Pappa R S，Elliott K B.Consistent－mode indicator for the eigensystem realization algorithm［J］.Journal of Guidance Control and Dynamics ，1993，16(5)：852 －858.

［12］ Yun G J.Modal identification and damage detection for structural health monitoring under ambient vibration environment［C］//Structures Congress 2009.Austin，USA：American Society of Civil Engineers，2009.

［13］ Pappa R S，James G H，Zimmerman D C.Autonomous modal identification of the space shuttle tail rudder［J］.Journal of Spacecraft and Rockets ，1998，35(2)：163 －169.

［14］ Illinois SHM Project：http：//shm.cs.uiuc.edu.

［15］侯立群，欧进萍，赵雪峰，等.哈尔滨四方台斜拉桥模态参数和索力识别［J］. 振动与冲击，2009，28(5)：106－110.

［16］胡利平.崖门大桥主桥斜拉桥－工程检测试验报告［R］.广东省交通建设工程质量检测中心，2002.

Modal identification and cable tension estimation of
long span cable-stayed bridge based on ambient excitation

YE Xi-jun1，YAN Quan-sheng1，LI Jian2，WANG Wei-feng1，ZHU Tian-feng1，LIU Ming-hui1
(1.School of Engineering and Transportation，South China University of Technology，Guangzhou 510641，China;2.College of Environmental Engineering，University of Illinois at Urbana－Champaign，Champaign 61801，USA)

For health monitoring of long span cable－stayed bridge，modal parameters and cable tension are the key parameters to assess the condition of the bridge.A finite element model of Yamen bridge was built.Modal test and cable vibration test of the bridge were performed under ambient excitation.An improved multiple reference DOFs stabilization diagram algorithm based on ERA(Eignsystem Realization Algorithm)was presented.By setting different reference DOFs in each group of data，NExT(Natural Excitation Technique)－ERA was used to identify modal parameters.Damping ratio，consistent mode indicator from observability(CMI－O)and modal assurance criterion(MAC)were used as thresholds to identify the most accurate modal parameters.Lower order frequencies were estimated by quadratic fit method，and the cable tension was estimated by two different methods.Based on the analysis of deck and cable vibration，it is evident that the vertical vibration of the bridge deck is tightly coupled with the cable vibrations within the frequency range of 0～3Hz.

ambient excitation;cable－stayed bridge;modal parameter identification;multiple reference DOFs stabilization diagram;cable tension estimation

TU311;O327

A

广东省交通运输厅科技项目(2010－02－015)

2011－03－30 修改稿收到日期：2011－09－09

叶锡钧 男，博士生，1984年10月生