卢海林，吴 玏，黄民水，徐文胜

(1．武汉工程大学环境与城市建设学院，湖北 武汉 430074； 2．华中科技大学土木工程与力学学院，湖北 武汉430074)

0 引 言

传统桥梁安全性评估主要依据表观损伤调查的结果.随着检测技术和计算机技术的发展，对于大型桥梁结构进行损伤识别，都通过建立有限元模型进行仿真分析[1].但是，按照桥梁设计图建立的有限元模型往往不能完全反映桥梁结构的真实状态，这就需要对有限元模型进行修正.本文采用基于优化理论的模型修正法修正大型桥梁结构动力分析模型.该方法在不影响桥梁正常使用的情况下，通过环境激励下桥梁的振动测试，采用功率谱峰值法来获得大型桥梁主要模态特征值，然后挑选对于模态数据有足够灵敏度的参数进行修正.参数型模型修正计算复杂，具有非线性特点，修正后的有限元动力分析模型可以更高精度的反映桥梁结构的真实状态.

通过该方法修正的动力分析模型用于桥梁结构损伤识别可以提高损伤量的判断和损伤的定位的精度.本文应用该方法对黄河某特大桥的有限元动力分析模型进行了修正.

1 模态识别

传统模态识别方法，为了建立频响函数必须同时知道系统的输入-输出，然而对大型桥梁结构，人工激励复杂昂贵，且效果并不理想.基于环境激励的模态识别方法克服了这些缺点，直接利用行人、风、地面微动等环境作用下，桥梁所受的激励是平稳的各态历经宽带随机力，其响应的主谐量是在其固有频率附近的振动，然后获得结构各阶频率.

环境激励模态参数识别方法主要有峰值法(PP)、频域分解法(FDD)、NExT-ERA法、随机减量法(RD-ITD)、时间序列分析法(AR/ARMA)、随机子空间法(SSI)[2]等.功率谱峰值法在工程上广泛应用.其原理为假设输入为平稳白噪音激励时，环境振动响应的自功率谱函数和频响函数具有相似关系[3]，如式(1).

Syy(ω)=H(ω)RuuHH(ω)

(1)

Syy(ω)为自功率谱函数，H(ω)为频响函数，Ruu为常数矩阵.如果结构阻尼比较小，模态频率近似平均分布，则在固有频率ωi处，可以将功率谱模态近似分解式为：

(2)

理论上，系统的各阶频率都应该可以在功率谱图中找到对应的峰值点.选取峰值点，就可以确定结构的频率.根据动力响应信号的功率谱密度曲线进行参数识别，快速准确，是桥梁检测中最常用的方法之一.但也存在峰值选取具有主观性的特点，因此需要检测人员具有相应的工程实际经验[4].

2 参数化模型修正

2.1 参数化有限元模型

考虑后期要通过ANSYS的优化功能进行参数化选代修正，前期有限元模型(FE)建立须考虑以下几个问题：

a.参数化模型建立.对于大型桥梁建模，利用ANSYS的GUI界面繁琐复杂，且不利于进行后期参数化迭代计算，所以使用APDL语言(ANSYS Parametric Design Language,参量化设计语言)建立模型几乎成了唯一的选择.

b.采用BEAM188单元.BEAM188单元适合分析从细长到中等粗短的梁结构，该单元基于铁木辛哥梁结构理论，并考虑了剪切变形.SECTYPE,SECDATA, SECWRITE,SECREAD, SECOFFSET等命令可以将自定义的梁截面写入，并调整截面中原始参考点，生成变截面梁；利用 *do,i,,命令可以进行大量循环创建常截面和变截面梁单元.同时，相对于Solid45单元建模后迭代计算费时太长，BEAM188单元计算耗时相对较少，适合大型工程实际应用.

c.建立较为完善的边界条件，选用适当的物理参数，对不同强度混凝土采用不同的弹性模量和密度.

2.2 参数型模型修正

模型修正是伴随有限元理论发展而发展起来的，它利用振动测试得到的动力特性来修正理论模型，使得修正后理论模型的模态参数与实测参数接近.类型大致可以分为：参数型模型修正，矩阵型模型修正[4].

参数型修正法则是对矩阵元素或结构的设计参数(如材料的E，ρ，A，I等)进行修正.这种方法物理意义明确，实用性强，能保持系统矩阵的对称带状特征，便于实际大型桥梁分析计算，并与其它结构优化设计过程密切相关.其主要缺点是计算复杂，对非线性目标函数的优化采用迭代算法，其修正效果也依赖于计算方法和测试精度，而且难以引入振型修正[5].

参数型修正方法又可以分为模态法与频响函数法，本文采用模态法，其思路与结构优化理论类似，基本原理如下[6]：

无阻尼结构动力方程为:

MA{x″}+KA{x}={0}

(3)

MA、KA为系统质量和刚度矩阵.理论频率和振型分别为ωA、ΦA，它们与实际结构的动力学特性，即通过实验模态分析分别获得的频率ωt和振型Φt，有一定的误差.两者之间差值为：

ΔΦ=Φt-ΦA

(4)

Δω=ωτ-ωA

(5)

引入参数ε，则实测值和理论值的关系可以表达为：

ωt=ωA+εΔω

(6)

Φr=ΦA+εΔΦ

(7)

设Mt、Kt具有ε的幂次级数形式：

Mt=MA+εΔM1+ε2ΔM2+L

(8)

Kt=KA+εΔK1+ε2ΔK2+L

(9)

由振型归一化正交性条件有：

(10)

KtΦt=MtΦtωt

(11)

对质量矩阵进行修正，将(5)、(6)代入(8)，可得

(ΦA+εΔΦ)T(MA+εM1+…)(ΦA+εΦ)=I

(12)

展开式(10)，收集同次幂的项，得

(13)

(14)

……

由式(12)得到

(15)

式中ΔΦ、ΦA和MA为已知，可以得到M的修正量ΔM.然后对K进行修正，现将式 (5)、(6)、(7)代入式(9)得：

(KA+εΔK1+…)(ΦA+εΔΦ)=

(MA+εΔM1+…)(ΦA+εΔΦ)(εA+εΔω)

(16)

展开式(14)，收集同次幂的项，得

KAΦA=MAΦAωA

(17)

式(15)是原方程的特征方程.

KAΔΦ+ΔK1ΦA=MAΦAΔω+ΔM1ΦAωA+MAΔΦωA

(18)

由式(16)得到:

ΔK1ΦA=MAΦAΔω+ΔM1ΦAωA+MAΔΦωA-KAΔΦ

(19)

ΔΦ、Δω都是已知的，若ΔM=0，可以得到K的修正量ΔK[7].

3 工程实例

3.1 有限元模型建立与模态分析

黄河某特大桥(图1)起点桩号K598+882.600，终点桩号K599+702.750，桥梁全长为820.15 m.孔跨布置为40+88+2×160+88(连续刚构)+7×40(m)(先简支后连续T梁)，桥共分4联.两侧各设50 cm防撞护栏，预制梁顶设6 cm混凝土调平层，FYT-1型防水材料、10 cm沥青混凝土桥面铺装，主桥下部桥墩为钢筋混凝土双肢薄壁变截面矩形空心墩，墩高86～100 m，墩顶双壁中距10 m.

图1 黄河某特大桥Fig.1 Grand bridge in the Huang River

a.主桥结构布置：88+2×160+88 m(连续刚构)预应力混凝土箱梁.

b.桥面宽度：10.0+2×0.5 m(防撞护栏).

c.道路等级：公路-Ⅰ级.

d.设计车速：60 km/h.

本文采用BEAM188，link2(模拟预应力钢筋)单元，采用参数化方法建立仿真模型；主梁混凝土采用C50，桥墩混凝土采用C40，预应力钢筋采用fpk=1 860 MPa钢绞线.对有限元模型(图2)进行模态分析，结果如表1所示.

图2 有限元模型Fig.2 Finite element model

3.2 振动测试

测试系统采用DH3817N，传感器有4个(如图3)，分别为DH610V(1个)，地震所941B型传感器(1个，V08358)，DH610H(2个)，3个测点，每个测点1-2号传感器测竖向加速度，3号测横桥向加速度，4号测纵桥向加速度.在第3、5跨选择3/4截面、1/2截面和1/4截面作为控制截面，每截面选择一个测点，该测点沿竖向布置两个拾振器，桥纵向和桥横向各布置一个拾振器.将拾振器布置在各测点，采用脉动法测试桥梁结构的脉动位移曲线，如图4.进行功率谱分析后得到桥梁结构的自振特性.图4中峰值都代表可能的固有频率，但还必须排除由于车辆随机制动力和竖向>冲击力造成的功率谱波段杂乱，所以要与ANSYS模型分析的频率进行比对和排查.理论值与实测值结果见表2.

图4 功率谱峰值法Fig.4 Method of power spectrum peak

表2 理论值与实测值Table 2 Theoretical values and measured values

3.3 修正参数

ANSYS可以利用PSD命令从主要参数(包括主梁高度h1,顶板厚度T1,底板厚度T2,主梁弹性模量E1、高墩弹性模量E2、箱梁混凝土密度ρ1)中选取模态参数灵敏度高的值.研究中发现，相对于各阶频率，主梁的混凝土弹性模量变化最为明显.假设主梁初始弹性模量设计值为3.55×1010，当混凝土弹性模量提高15%，竖向基频增加约28.8%，横向基频增加约48.5%，纵向基频增加约5.4%，即竖向频率与横向频率对弹性模量变化比较敏感，而纵向频率对主梁弹性模量的变化反应不明显，并且在后续损伤识别过程中，通过表征刚度的参数来判别损伤是常用手段之一.根据结构动力学原理，刚度元素中包含表达式EI或EA，都与弹性模量有直接的关系，所以，本文选择弹性模量EI为主要修正参数.

利用ANSYS的APDL语言，将第三跨1号单元到22号单元划分为四段，1到2号单元弹性模量设置为EX1，3到9号单元弹性模量设置为EX2，10到16号单元弹性模量设置为EX3，17到22号单元弹性模量设置为EX4.同理，根据对称性，第六跨可以布置为四段，分别为EX22，EX23，EX24，EX25.三个0号箱分别设置为EX5，EX13，EX21，第四跨分别为EX6，EX7，EX8，EX9，EX10，EX11，EX12，第五跨分别为EX14，EX15，EX16，EX17，EX18，EX19，EX20.划分原则是将跨中部位箱梁的顶板和底板尺寸变化不大者，划分为同一弹性模量段，将梁高、顶底板厚度按照二次曲线变化划分为三段.这样，就将桥梁弹性模量划分为25个区段，保证修正后能初步判断箱梁刚度下降段位置.

初次计算时，赋予EX1～EX25=3.55×1010，运用： opvar，EXi，DV，0，4.0e10，1e4将EX1至EX25定为设计变量，其上限为4.0×1010，下限为0，步长为1.0×104；也可以将下限值适当提高，适当缩小搜索范围，降低收敛运算次数.

将理论频率和实测频率之间相对误差的平方和定义为目标函数，分别提取表2中的两组值，得到公式(20)和(21).

(20)

(21)

再以桥梁结构在环境激励下的动力响应作为约束条件和目标函数，在一定范围内搜索设计变量的取值.寻求f(x)最小值的过程，就是不断调整目标函数，从而使主梁弹性模量不断向实际值逼近的过程.利用ANSYS命令：*get，freqi，mode，i，freq 来提取需要模态的频率，然后使用ANSYS命令：opvar，f(x) ，obj，，，1e-4 来实现目标函数的最小化.其中步长为1×10-4[8].

4.4 修正结果

ANSYS通过循环30次，其中第8次迭代为最优解，最终结果如表3.

表3 迭代后弹性模量值Table 3 The elastic modulus value after iteration

通过迭代修正后，频率与修正前频率对比如表4所示.

从修正后的模型可以快速判断桥梁健康状况，也可为后续的损伤识别提供依据.

a.总体动力特性良好；修正后除个别梁段外，弹性模量均高于规范值，反映出桥梁总体满足承载能力的设计要求.

b.预应力混凝土构件由于全截面(无裂缝)或接近全截面(裂缝控制)参与工作，所以其结构刚度本应该比普通混凝土构件大[9]，但从修正结果可以看出，局部梁段混凝土弹性模量比规范值略低，分别集中在：第三跨悬臂梁端部、第四跨左侧根部、第五跨右侧根部，如图5所示.图中颜色深度表示弹性模量下降的程度.分析其原因主要有以下几点.

表4 迭代后模型频率与实测值对比Table 4 After iteration model frequency compared with the measured values

图5 梁段弹性模量Fig.5 Elasticity modulus of beam sections

①现场观察，大型运煤车队通过该桥频繁来往于附近煤矿，车队超载较为严重.主梁根部负弯矩和剪力较大.同时超载也可能使主梁跨中挠度超过限定值，进一步增大梁根部弯矩[10].

②混凝土浇筑质量问题.箱梁内部表面可见蜂窝和裸露钢筋，可能是混凝土浇模时，保护层垫块移位或混凝土漏振、振捣不均匀造成的.

③第三跨悬臂梁处弹性模量下降可能由于引桥部分桥墩不均匀沉降、支座侧斜或转动不自如等造成的，原因需要进一步调查.

④混凝土的弹性模量本身具有一定的离散性，这可能也是修正后的弹性模量与规范值有偏差的原因之一.

4 结 语

a.有限元模型的修正使模型的动力特性与实际情况更加接近，各阶频率修正后的误差基本在2%以下，精度较高.

b.本文所采用的基于优化理论的有限元动力模型修正方法，利用现有的有限元软件ANSYS，不用再次输出繁琐的M、K和C矩阵即可进行有限元模型修正，快速有效，适合于大型桥梁工程实际应用.

c.经过有限元模型修正，理论模态值与实测模态值基本吻合，可以作为桥梁健康监测和损伤识别的依据，也可为后续桥梁预防性养护等工作提供理论指导.

致谢

本研究得到湖北省教育厅和武汉工程大学的经费资助，在此一并表示衷心地感谢！

参考文献：

[1] 王元清，姚南，张天申，等.基于最优化理论的多阶段模型修正及其在桥梁安全评估中的应用[J].工程力学，2010,27(1):91-97.

WANG Yuan-qing,YAO Nan,ZHANG Tian-shen,et al. An application of multistage model updating based on optimization theory to the safety appraisal of bridge[J], Engineering Mechanics, 2010,27(1):91-97.(in Chinese)

[2] 陈宇．有限元模型修正技术在桥梁工程中的应用研究[D]. 成都：西南交通大学,2009.

CHEN Yu.Research on the Applic-ation of Finite Element Model Updating Technique in Bridge Engineering[D].Chengdu:Southwest Jiao-tong University,2009.(in Chinese)

[3] 刘继承，周传荣.一个基于优化的有限元模型修正方法[J].振动与冲击，2010，22(2):330-340.

Liu Jicheng,Zhou Chuanrong. FE model updating method based on optimization[J].Journal of Vibration and Shock, 2010，22(2):330-340. (in Chinese)

[4] 荣见华，郑建龙，许飞鸿.结构动力修改及优化设计[M].北京：人民交通出版社，2011.

[5] Salama M, Rose T, Garba J.Optimal placement of exciters and sensors for verification of large dyna-mical systems[C]//Kedward K-T. Proceedings of SDM conference,AIAA-87-0782,1987:1024-1031.

[6] Imamovic N.Model validation of large finite element model using test data[D].London: Imperial College,1998.

[7] Berman A, Nagy E J. Improvement of a large analy-tical model using test data[J].AIAA J,1983,21(8):1168-1173.

[8] 朱宏平，徐斌，黄玉盈.结构动力模型修正方法的比较研究及评估[J].力学进展，2002,32(4)：513-525.

Zhu Hongping, Xu Bin, Huang Yuying. Comparison and evaluation of analytical approaches to structural dynamic model correction[J]. Advances in Mecha-nics，2002,32(4)：513-525.(in Chinese)

[9] 白宝玉，王丽荣.桥梁工程[M].北京：高等教育出版社，2002.

[10] 邵旭东，程翔云，李立峰.桥梁设计与计算[M].北京:人民交通出版社，2006.