共扼梯度法对简支梁桥的损伤识别研究

朱泽文，方有亮

(河北大学建筑工程学院，河北保定071000)

国内外关于结构损伤识别工作主要开始于20世纪80年代初，经过30多年的发展，已经在结构损伤识别方面取得了一些成就，但离实际的应用还有一段距离。桥梁结构损伤识别主要有静力和动力损伤识别，以及两者结合的损伤识别。然而，因为动力识别的优越性，目前基于动力的结构损伤识别尤其受到关注。因为该方法具有非破坏性、方便、快速和廉价的优点且具有广阔的工程应用前景。从理论角度来讲，由于局部损伤的产生，该区域的结构刚度和承载能力将会有所下降, 而结构的质量特性一般不会受到影响，结构的模态参数(刚度、振型，模态阻尼)将随之而改变，因此通过研究结构的振动特性来识别结构是否发生损伤, 并确定损伤的位置和程度是一种可行的途径。

1 结构模型修正

关于结构振动的有限元模型修正中，常采用基于参数的模型修正和直接矩阵逼近的模型修正方法[1]、[2]。本文主要基于参数的模型修正，将正弦激励力的移动荷载作用在简支梁上，根据正问题所得的响应逆算结构的刚度参数。文献[3]中将共扼梯度法与遗传算法结合，应用在航天领域的结构损伤识别。文献[4]中，C·-H·HUANG将共扼梯度方法应用在非线性的有阻尼系统的外部力的识别上。文献[5]中，将共扼梯度和遗传算法相结合，对一桁架及悬臂梁进行损伤识别。基于上述文献中相关问题的研究，本文也将共扼梯度和遗传算法结合应用在简支梁桥的损伤识别上，不同的是将正弦激励移动荷载作用于简支梁桥上，通过所得的位移响应数据对简支桥梁的损伤程度进行评估。本文将根据简支梁中单元受损位置、损伤程度情况进行讨论分析损伤识别效果。

2 简支梁模型与动力学方程

图1 简支梁模型示意

(1)

本文中所用阻尼矩阵为Rayleigh阻尼：

C=αM+βK(α，β为比例常数)

(2)

3 逆问题分析

3.1 正问题

对于任意的K，C以及f(t)，利用数值计算方法(Newmak-β)对式(1)进行数值求解，将求得的响应逆算结构的损伤参数。

3.2 灵敏度问题

本文将弹性模量E作为刚度损伤因子，假定每个单元的弹性模量均不同。该刚度矩阵不考虑轴向力及轴向位移的作用。令K=f[E]，其中E=[E1,E2,…ENe]，Ne表示单元个数。为了推导每个未知刚度Ki的灵敏度问题，则假定每个未知的刚度参数同时被扰动，Ki的增量为ΔKiδij，δij是克罗内克函数，j=1,…Ne，同时，位移响应也发生微小的扰动Δxij(t)，在这Δxij表示当第j个单元的参数发生变化时，第i个节点位移变化量。将Ki用Ki+ΔKiδij，xi(t)用xi(t)+Δxij(t)替代，代入式(1)中，且与(1)式相减，忽略二阶微量式，则可推导出Ne个单元的灵敏度问题方程：

(3)

其中Δxj=[Δx1j,Δx2j，……，ΔxNj]，

ΔKj表示第j个单元的弹性模量发生改变后的单元刚度矩阵，因此当单元没有损伤时，此时的ΔKj为0

3.3 伴随问题及共扼梯度方程

逆问题的解答其实就是求解最优化的问题，定义目标函数为：实测数据与计算数据的差的平方和，即：

(4)

tf表示计算最终的时间，x(K,t)为计算响应值，xexp(t)为实际值。

为了求解式(4)的最小值，在这引入拉格朗日算子λ(t)，因此式(4)可写成：

(5)

将Ki用Ki+ΔKiδij，xi(t)用xi(t)+Δxij(t)替代代入式(5)并与之相减且忽略二阶微量式，可得：

(6)

对式(6)进行分部积分且将初值条件代入，即可推导出关于λ(t)的向量式子为：

(7)

(8)

求得伴随问题后，我们可以从式(6)和式(7)推导出：

(9)

由文献[8]可知，J[K]的方向导数可写成：

(10)

(11)

3.4 共扼梯度方法

共扼梯度方法的迭代程序中关于刚度变化参数用表示：

Kn+1=Kn-βnPnn=1,2,……

(12)

式中：β表示步长，Ρn表示共扼梯度方向(也称为最佳搜索方向)

Pn=J′n+γnPn-1(γ0=0)

(13)

γn为一线性组合因子，在这应注意，当n=0时，γ0=0，

(14)

在式(12)中出现的βn，可以通过式(4)，将Kn+1(t)换成Kn(t)

(15)

将式(15)中Kn+1=Kn-βn代入，利用泰勒展开即可求得：

(16)

3.5 终止准则

J[K]<ε*(10-15)

(17)

3.6 遗传算法

遗传算法是模拟达尔文的遗传选择和自然淘汰的生物进化过程，对群体反复进行基于遗传学的操作(遗传、交叉和变异)。根据预定的目标适应度函数对每个个体进行评价，依据适者生存、优胜劣汰的进化规则，不断得到更优的群体，同时以全局并行搜索方式来搜索优化群体中的最优个体，以求得满足要求的近似最优解。本文遗传算法的目标适应度函数与共扼梯度的目标函数相同：

(18)

式中：xexp(t)为实测响应数据；x(K,t)为更新后的刚度矩阵通过正算得到的响应数据。

4 移动荷载

文献[6]中，将移动荷载模拟成一质量块与结构通过弹簧连接，这样移动荷载与桥梁就存在了相互的作用，然后通过不同的速度来研究结构的动力特性。文献[7]中，利用动态规划的方法识别结构中的移动荷载。基于以上文献的考虑，本文将移动荷载的思路运用到桥梁识别上，移动速度为18 km/h，从梁的左端移动到梁的右端。利用Newmak-β数值方法计算结构动态响应，计算时间8 s，时间步长0.01 s，γ=0.6，β=(0.5+γ)2/4

移动荷载：

f=104×[1+0.5×sin(2πt)+0.5sin(5πt)]

(19)

当荷载作用在第i单元非节点上时，根据结构力学矩阵位移法中关于等效结点荷载处理可知，将集中荷载等效到该单元的两个端节点上。

5 损伤识别

本节仅用共扼梯度法在有噪音的情况下进行损伤识别，按照下式给损伤后的结构响应数据加入噪音：

xext(t)=x(t)×(1+σR)

(20)

其中σ为标准误差，R为随机的数([0,1])。在这里σ取1%，仅为检验该方法对噪音的抗干扰能力，实际的噪音可能更高。

5.1 单一损伤

图2 5号15%损伤

图3 1号15%损伤

图2、图3表示当跨中的5号和支座处1单元受损伤15%时，在有噪音干扰情况下，计算与实际损伤值的对比。可以看出，识别效果虽有误差，但相差不大。 这说明，对于单个单元受损，不管位置在哪，损伤程度如何，仅用共扼梯度方法，任给初始值，都能识别得到较好的结果。

5.2 多单元受损

以下主要针对三种情况进行损伤识别：(1)损伤程度相同，损伤位置对称(靠近支座处)。( 2)损伤程度不同，损伤位置对称(靠近支座处) 。(3)损伤程度不同，损伤位置不对称。

图4 1、2、9、10号损伤20%

图5 1、2号损伤20%，9、10号损伤50%

图6 1、2号损伤20% ，5、6号损伤50%

由图4可知，当结构有四个单元受损伤时，损伤程度一样时，对称的布置于两支座处时，损伤识别效果并不理想。由图5可知，当损伤的位置对称布置于两支座处，但是损伤程度不一样时，加入噪音后，识别效果基本没有变化。由图6可知，当结构受损伤的程度不一样，损伤的位置也不对称时，虽跨中的两个单元计算损伤值与实际的损伤值相差不大，但这主要与初始的估计值有关。整体的识别效果不理想。

6 用混合方法对单元损伤进行识别

本节中引入遗传算法计算初始值，主要解决上节中仅用共扼梯度法，任给初始值不能很好识别结构损伤问题。因为共扼梯度法计算需要一个初始值，然而当多单元损伤，多程度损伤时，这时的初始值的选择对于识别的结果至关重要。所以在本节中，主要是利用遗传算法先估算出一个初始值，再用该结果作为共扼梯度算法的初始值。这种混合方法的运用可以比较好的识别结构的损伤。

图7 1、2号损伤 20%，5、6号损伤 50%

图8 1、2号损伤 20%，9、10号损伤 50%

从图7、图8可知，利用混合方法多单元进行损伤识别(加入了噪音)，虽然结果和真实值不是很吻合，但整体的效果相比仅用共扼梯度法进行识别要好的多。

7 结论

本文将正弦激励移动荷载作用在简支梁桥损伤识别上，可以更好地激励桥梁结构的响应，从而可以更加有效地利用速度响应数据对结构进行损伤识别。

通过以上的识别效果表明，当仅用共扼梯度法对桥梁结构从理论上进行识别时，对于单个单元受损，识别效果很好，对于多单元损伤，当损伤的为跨中附近单元且损伤程度一样时，可以较好识别，然而当损伤程度不一样时，识别效果较差。这主要与初始值的选择有关，当用共扼梯度和遗传算法联合运用对移动荷载作用下桥梁结构的损伤识别从理论上能得到较好的识别效果。加入噪音后，噪音仅对靠近支座处损伤程度较低的单元识别有影响，对于其它位置的损伤识别影响不大，这说明该方法对噪音的抗干扰能力较强。

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

[2] 唐晓峰，王皓，唐国安.有限元动力学模型的综合修正方法[J].上海航天，2008(2)

[3] Leonardo D, Chiwiacowsky，Paolo Gasbarri，Haroldo F. de Campos Velho，Damage assessment of large space structures through the variational Approach[J]. Acta Astronautica 2008(62):592-604

[4] C·-H·HUANG. A INVERSE NON-LINEAR VIBRATION PROBLEM OF ESTIMATING THE EXTERNAL FORCES FOR A SYSTEM WITH DISPLACEMENT-DEPENDENT PARAMETERS，Journal of Sound and Vibration 2001, 242(5):749-765

[5] LEONARDO D.CHIWIACOWSKY，HAROLDO F.DE CAMPPS VELHO and PAOLO GASBARRI A variational approach for solving an inverse vibration problem[J]. Inverse Problems in Science and Engineering 2006,14(5):557-577

[6] Y·B·Yang，C·W·Lin.Vehicle-bridge interaction dynamics and potential applications，Journal of Sound and Vibration 2005(284)：205-226

[7] S·S·Law Y·L·Fang. MOVING FORCE IDENTIFICATION: OPTIMAL STATE ESTIMATION APPROACH[J]. Journal of Sound and Vibration 2001, 239(2):233-254

[8] Y·Jarny，M·N·OZISIK，J·P·Bardon.A general optimization method using adjoint equation for solving multidimensional inverse heat conduction[J].International Journal of Heat and Mass Transfer 1991,34(11):2911-2919