桥梁结构在外界荷载的长期作用下会出现构件损伤或材料性能退化等问题, 累积损伤最终可能会导致桥梁坍塌.及时准确地识别桥梁结构损伤或材料性能退化一直是桥梁工程领域研究的热点问题.国内外已有许多国家和地区的桥梁布设了结构健康监测系统, 用来开展桥梁结构状况评估, 张启伟等^[1-5]基于结构健康监测数据作了大量的研究工作.结构损伤识别作为结构健康监测的重要任务, 得到广泛关注^[6-7].

中性轴为受弯梁中性层与梁横截面的交线, 是所有梁式结构共有的结构参数.理论上, 在未损伤的状况下受弯梁横截面的中性轴位置不会发生变化.当横截面发生损伤(开裂、分层、材料性能退化等)时, 横截面质心随之改变, 中性轴位置会发生移动.中性轴位置的改变可以反映梁截面刚度的变化, 也可以反映其他影响梁截面质心的结构行为, 如预应力损失、截面材料性能退化、截面应力重分布等.中性轴位置可以用于梁式结构的结构损伤识别、性能预测与状态评估.

目前, Sigurdardottir等^[8]研究基于中性轴位置的结构损伤识别.Plude等^[9]对多主梁组合钢箱梁桥采用有限元模拟, 开展损伤识别研究, 发现在固有频率、应变分布与中性轴位置3种损伤识别参数中, 中性轴位置具有较高的敏感性.O’Brien等^[10]指出桥面板上荷载分布位置会引起中性轴位置改变.Gangone等^[11]在一座即将摧毁的桥梁上开展控制试验, 发现中性轴位置对桥梁结构的不均匀沉降具有敏感性.

本文采用九堡大桥结构健康监测系统获得的应变数据, 基于中性轴位置参数进行损伤识别分析.基于小波多分辨率分析对长期监测数据进行处理, 去除应变监测数据中的低频(温度)组分.采用卡尔曼滤波方法对高频(交通荷载)组分进行滤波处理.采用滤波后的应变数据得到梁截面中性轴的变化情况, 对梁的损伤状况进行识别.

1 受弯梁中性轴位置确定 1.1 受弯梁变形分析

基于材料力学的基本理论可知, 受弯梁在纯弯矩荷载下发生弯曲变形.假设梁是由无数层纵向纤维组成, 中性层是梁上拉伸区与压缩区的分界面, 中性层与横截面的交线为中性轴, 如图 1所示.o-o、o′-o′为变形前、后的中性轴, c-c、c′-c′为变形前、后的纵向纤维；Δx、Δx′为变形前、后的纤维长度；h为c-c纤维距中性轴的距离；ρ为受弯梁的曲率半径；dθ为夹角.c-c纤维的应变可以通过下式计算：

$ \varepsilon = \mathop {\lim }\limits_{\Delta x \to 0} \frac{{\Delta x' - \Delta x}}{{\Delta x}} = \mathop {\lim }\limits_{\Delta x \to 0} \frac{{\left( {\rho + h} \right){\rm{d}}\theta - \rho {\rm{d}}\theta }}{{\rho {\rm{d}}\theta }} = \frac{h}{\rho }. $

(1)

式中：ε为c-c纤维的应变.

纵向纤维的应变与其距中性轴的距离成正比, 即截面应变沿梁厚度方向呈线性分布.另外, 与中性轴距离为h的应变可由下式求解：

$ \varepsilon = \frac{{Mh}}{{EI}}. $

(2)

式中：M为弯矩, E为材料弹性模量, I为截面惯性矩, h方向取中性轴向下为正.

1.2 中性轴位置

基于欧拉-伯努利梁理论可知, 受弯梁横截面应变可以根据几何关系确定.受弯梁横截面的应变分布见图 2.图中, ε₁、ε₂为应变, 均取应变的绝对值；h₁、h₂为应变距中心轴的距离.根据几何关系, 可知：

$ \frac{{{\varepsilon _1}}}{{{\varepsilon _2}}} = \frac{{{h_1}}}{{{h_2}}}. $

(3)

式中：ε₁、ε₂可以通过传感器测量.

取r表示受弯梁中性轴位置, 中性轴位置可以通过下式^[12]求解：

$ r = \frac{{{\varepsilon _1}}}{{{\varepsilon _1} + {\varepsilon _2}}} = \frac{{{h_2}}}{h}. $

(4)

2 应变监测数据处理方法 2.1 小波多分辨率分析

一般情况下, 根据结构应变监测数据推算结构应力不容易, 因为应变监测数据通常为多组分数据(包含荷载效应、温度效应与噪声等).一方面, 各组分所引起的结构应力不同.比如, 温度效应引起的桥面应变虽然很大, 但伸缩缝位移与温度应变呈线性关系^[13], 通过伸缩缝位移不一定会产生结构应力.解决该问题的有效方法是将多组分数据进行单组分分解.小波多分辨率分析可以很好地实现该目标^[13].

信号变换域处理为使用一组基函数将信号从信号域转换到转换域.小波变换作为一种变换域信号处理方法, 选择某一特定函数(母小波)表示信号, 通过改变母小波的尺度参数, 即伸缩平移运算, 将信号分解并重构成多种不同尺度的信号.基于小波变换的多分辨率分析过程是对信号进行多尺度分解与重构的过程.

使用函数F(t)表示母小波, 小波变换过程可以表示为

$ {\left( {{W_F}f} \right)_{\left( {a,b} \right)}} = \frac{1}{{{{\left| a \right|}^{1/2}}}} \times \int_{ - \infty }^\infty {f\left( t \right)\bar F\left( {\frac{{t - b}}{a}} \right){\rm{d}}t} ;a > 0. $

(5)

式中：(W_Ff)_{(a, b)}为小波变换过程, a为尺度参数, b为时间参数, F(t)中上横线表示共轭, f(t)为输入信号函数.f(t)可以通过对(W_Ff)_{(a, b)}双重积分得到：

$ f\left( t \right) = \frac{1}{{{C_F}}} \times \int_{ - \infty }^\infty {\int_{ - \infty }^\infty {{W_F}f\left( {a,b} \right)\bar F\left( {\frac{{t - b}}{a}} \right)\frac{1}{{{a^2}}}{\rm{d}}a{\rm{d}}b} } . $

(6)

式中：C_F为容许条件参数.

一般情况下, 由于a与b的离散性, 通常采用离散小波变换进行实际信号的处理.为了提高小波变换过程的有效性, 尺度参数与时间参数均采用二进制表示.选择F(t)作为母小波表示原始信号, 相应的离散小波{F_{j, k}}可以用尺度参数与时间参数表示为

$ {F_{j,k}}\left( t \right) = {2^{j/2}} \times F\left( {{2^j}t - k} \right). $

(7)

式中：{F_{j, k}}为属于L²(R)空间的标准正交基.f(t)可以表示为

$ f\left( t \right) = \sum\limits_j {\sum\limits_k {{\alpha _{j,k}}{F_{j,k}}\left( t \right)} } , $

(8)

$ {\alpha _{j,k}} = \int_{ - \infty }^\infty {f\left( t \right){F_{j,k}}{\rm{d}}t} . $

(9)

式中：α_{j, k}为变换参数.如果用D_j、A_j表示j级的细节部分与近似部分, 则

$ {D_j} = \sum\limits_{k \in Z} {{\alpha _{j,k}}{F_{j,k}}\left( t \right)} , $

(10)

$ {A_J} = \sum\limits_{j > J} {{D_j}} , $

(11)

$ f\left( t \right) = {A_J} + \sum\limits_{j \le J} {{D_j}} . $

(12)

式中：J为级数.

式(7)~(12) 给出小波多分辨率分析的树形结构分解和重建过程.通过小波多分辨率分析可知, 原始信号可以分解为多尺度信号.其中, 细节部分具有原始信号的高频信号特征, 近似部分具有原始信号的低频信号特征.

选取的小波种类不同会引起小波变换结果的差异, 小波函数的正交性、支撑集、正则性、消失矩阶数等特性决定了在数据处理应用中如何进行选择^[14].Daubechies小波简称db小波, 是90年代初由法国学者Daubechies提出.db族小波具有正交性并且是紧支撑集, 小波分解系数可以精确重构并降低分解信号时的能量散失.正则性和消失矩阶数特性可以保证小波更好地刻画信号特征, 具有较好的滤噪能力^[14], 缩短因支撑集太大造成的计算时间冗长.

2.2 卡尔曼滤波

1960年, 匈牙利裔鲁道夫·卡尔曼提出采用递归方法进行离散数据线性滤波的方法, 即卡尔曼滤波方法.卡尔曼滤波过程是一个信号处理与更新过程, 实质为根据实时观测的数据(包含噪声及其他误差), 以最小均方误差估计作为最佳估计准则, 获得更趋于真实的观测数据^[15-18].卡尔曼滤波方法的优势是可以适用于平稳与非平稳数据的滤波, 具有很好的抗噪声干扰能力^[12].

卡尔曼滤波理论认为线性系统可以描述为由状态方程与观测方程构成.

$ {\mathit{\boldsymbol{X}}_k} = {\mathit{\boldsymbol{A}}_k} \times {\mathit{\boldsymbol{X}}_{k - 1}} + {\mathit{\boldsymbol{B}}_k} \times {\mathit{\boldsymbol{u}}_k} + {\mathit{\boldsymbol{w}}_k}, $

(13)

$ {\mathit{\boldsymbol{Z}}_k} = {\mathit{\boldsymbol{C}}_k} \times {\mathit{\boldsymbol{X}}_{k - 1}} + {\mathit{\boldsymbol{v}}_k}. $

(14)

式中：X_k、Z_k为系统k时刻的状态变量与观测变量, X_k-1为系统k-1时刻的状态变量；u_k为系统k时刻的控制变量；w_k、v_k为系统k时刻过程噪声变量与观测噪声变量, 通常假设两者为均值为0、相互独立的白噪声；A_k、B_k、C_k分别为状态转移矩阵、状态控制矩阵和观测矩阵.式(13) 为状态方程, 式(14) 为观测方程.

卡尔曼滤波过程可以通过下列迭代方程实现：

$ {{\mathit{\boldsymbol{\bar X}}}_k} = {\mathit{\boldsymbol{A}}_k} \times {\mathit{\boldsymbol{X}}_{k - 1}} + {\mathit{\boldsymbol{B}}_k} \times {\mathit{\boldsymbol{u}}_{k - 1}}, $

(15)

$ {\mathit{\boldsymbol{P}}_k} = {\mathit{\boldsymbol{A}}_k} \times {{\mathit{\boldsymbol{\bar P}}}_{k - 1}} \times \mathit{\boldsymbol{A}}_k^{\rm{T}} + \mathit{\boldsymbol{R}}, $

(16)

$ {\mathit{\boldsymbol{X}}_k} = \frac{{{\mathit{\boldsymbol{P}}_k} \times \mathit{\boldsymbol{C}}_k^{\rm{T}}}}{{{\mathit{\boldsymbol{C}}_k} \times {\mathit{\boldsymbol{P}}_k} \times \mathit{\boldsymbol{C}}_k^{\rm{T}} + \mathit{\boldsymbol{Q}}}}, $

(17)

$ {\mathit{\boldsymbol{X}}_k} = {{\mathit{\boldsymbol{\bar X}}}_k} + {\mathit{\boldsymbol{K}}_k} \times \left( {{\mathit{\boldsymbol{Z}}_k} - {\mathit{\boldsymbol{C}}_k} \times {{\mathit{\boldsymbol{\bar X}}}_k}} \right), $

(18)

$ \begin{array}{l} {{\mathit{\boldsymbol{\bar P}}}_k} = \left( {\mathit{\boldsymbol{G}} - {\mathit{\boldsymbol{K}}_k} \times {\mathit{\boldsymbol{C}}_k}} \right) \times {\mathit{\boldsymbol{P}}_k} \times {\left( {\mathit{\boldsymbol{G}} - {\mathit{\boldsymbol{K}}_k} \times {\mathit{\boldsymbol{C}}_k}} \right)^{\rm{T}}} + \\ \;\;\;\;\;\;\;{\mathit{\boldsymbol{K}}_k} \times \mathit{\boldsymbol{R}} \times \mathit{\boldsymbol{K}}_k^{\rm{T}}. \end{array} $

(19)

式中：G为单位矩阵, R、Q为过程噪声协方差矩阵与观测噪声协方差矩阵, u_k-1为系统k-1时刻的控制变量, X_k为k时刻的估计值, P_k为k时刻的误差协方差, P_k-1、P_k为k-1、k时刻的更新估计误差协方差, K_k为卡尔曼增益.式(15) 为状态估计方程, 若没有确定的控制变量, 则右边第2项可以略去；式(16) 为状态估计下的误差协方差方程；式(17) 为卡尔曼增益方程；式(18) 为更新状态估计方程；式(19) 为更新估计误差协方差方程.

3 基于中性轴位置损伤识别应用 3.1 应变监测系统

一般情况下, 大跨径钢箱梁桥在交通荷载作用下发生弯曲时, 可以将钢箱梁视为受弯梁.九堡大桥是一座横跨钱塘江的钢混组合拱桥, 全长1 855 m, 如图 3所示.桥梁建成后布设了长期结构健康监测系统, 通过对多种荷载作用下的桥梁结构响应进行实时采集与分析, 评估大桥在运营期间的健康状况和工作性能.结构应变监测是健康监测系统的组成部分, 九堡大桥在拱肋、钢箱梁内部共布设166个光纤光栅应变传感器, 对应变状况进行在线监测, 钢箱梁光纤光栅应变传感器布设见图 4.图中, STR_G11_1表示传感器的编号, 其中STR表示应变传感器, G11表示梁截面编号, 1表示该截面上传感器的序号.如图 5所示为实际安装在桥面板上的光纤光栅传感器实物图.

3.2 数据处理

选择钢箱梁下游截面(G11) 中的STR_G11_1与STR_G11_4两个传感器的监测数据进行截面损伤识别.如图 6所示为钢箱梁顶板与底板一周时间内的应变时程.可以发现, 温度的日循环周期变化导致应变时程存在周期性变化的趋势；温度的周期性变化引起钢箱梁发生顺桥向周期性的收缩或膨胀.

采用小波多分辨率分析方法, 对顶板应变时程进行多尺度分析.采用db4小波对应变时程进行多层分解, 基于变换参数进行细节应变(交通荷载组分)与近似应变(温度组分)的重构.应变时程多分辨率分解如图 7所示.图 7表明, 温度引起的应变比荷载引起的大, 但因伸缩缝位移会将温度引起的应变释放, 对钢箱梁应力的影响较小, 在进行中性轴位置损伤识别时要去除温度的影响.使用小波多分辨率分析方法分离信号中高、低频组分, 具有很好的效果.

3.3 确定中性轴位置

选择一周内某一天3个时间段(0:00~1:00, 7:00~8:00, 8:00~9:00) 的应变数据, 采用卡尔曼滤波方法滤波, 确定中性轴位置的变化情况, 结果如图 8所示.图 8表明, 在交通荷载的作用下, 滤波后的中性轴位置的波动幅度较滤波前的小很多.在不同时间段, 卡尔曼滤波方法均具有很好的稳定性.在不同荷载作用下, 经过卡尔曼滤波后的中性轴位置没有发生较大的移动, 均在0.5附近.

卡尔曼滤波通过对系统状态的初始值和过程噪声与观测噪声协方差矩阵参数赋值, 基于前一时刻对当前时刻的估计值与当前时刻的观测值, 根据迭代公式实现对当前时刻的估计值的更新, 从而完成观测数据的滤波工作, 有效地降低了噪声对应变数据的影响.滤波后的结果表明, 卡尔曼滤波具有很好的稳定性和较强的抗噪声干扰能力.

3.4 损伤识别

标准差能够反映一组数据的离散程度, 采用标准差对卡尔曼滤波效果进行评定.滤波前、后的标准差对比见表 1.可以发现, 经过卡尔曼滤波方法处理的信号标准差较小, 卡尔曼滤波方法具有很好的抗噪声能力.在不同时间段内, 滤波后的中性轴位置与均值的变化波动为-2.24%与2.87%, 滤波前的中性轴位置与均值的变化波动为10.6%与-14.7%, 表明在交通荷载作用下, 基于应变监测数据计算的中性轴位置会发生波动, 但均围绕理论均值位置变动, 即发生可回复性变动, 钢箱梁截面未产生损伤.

若采用未进行滤波的数据对结构进行损伤识别, 因数据存在较大的离散程度, 不同时间段内的均值具有较大的误差, 基于该均值进行损伤识别容易造成误判.

4 结论

(1) 基于中性轴位置方法, 对九堡大桥桥面板进行损伤识别, 应变数据来自于九堡大桥结构健康监测系统.采用小波多分辨率分析方法去除监测数据中的温度组分, 采用卡尔曼方法对监测数据进行滤波.通过研究发现, 基于交通荷载作用下的应变数据计算的中性轴位置会发生变动, 若截面未发生损伤, 则中性轴会围绕理论均值变动, 不会发生不可回复的变动.

(2) 采用小波多分辨率分析与卡尔曼滤波方法, 获得了较好的效果, 但小波多分辨率分析具有高频段分辨率低的缺点, 卡尔曼滤波方法需要已知监测数据中噪声的先验知识, 初始值的选择对滤波结果有较大影响.

(3) 在结构未发生损伤的情况下, 基于应变监测数据分析获得桥面板中性轴位置发生的可回复移动特性.对于如何实现基于中性轴位置的结构损伤定量评估以及如何确定中性轴变动阈值, 需要进一步的研究.