中、小跨径的简支梁桥数量大、分布广，但简支梁桥在使用荷载、结构缺陷、环境因素、材料衰变等影响下，极易出现混凝土开裂病害，降低桥梁的承载能力。但简支梁桥结构形式简单，数量众多，不利于构建桥梁结构健康监测系统。故在桥梁使用过程中，对混凝土表观裂缝进行简单高效的识别是保障桥梁运营安全的重要手段^[1]。实现该类桥梁损伤识别的前提是掌握带裂缝简支梁桥在运营荷载作用下的工作机理。不少学者认为裂缝在桥梁振动过程中一直处于张开状态，即开口裂缝模型，并针对其动力特性和损伤识别做了大量研究^[2-4]。1983年，Gudmunson^[5]通过试验发现裂缝在振动过程中将出现张开闭合的周期性交替过程，为准确模拟这一机制，各国学者提出了呼吸裂缝模型^[5-12]。Cheng等^[9]以带呼吸裂缝的悬臂梁为研究对象，发现其频率降低量远小于开口裂缝模型。Chondros等^[10]、Vigneshwaran等^[13]针对呼吸裂缝的研究均得到类似规律。因此，以开口裂缝模型计算的频率、振幅将会低估结构的损伤程度，从而造成安全隐患^[12]。

当前针对呼吸裂缝的研究多集中在数值模拟模型及动力特性方面，而采用呼吸裂缝模型进行损伤识别的研究较少，常采用简化的移动力模型、移动质量模型^[14]或单自由度弹簧质量模型^[15]，不能客观地模拟车桥之间的耦合作用，在一定程度上降低了识别结果的准确性。文献[16]采用车桥耦合振动系统数值模拟手段，详细对比了呼吸裂缝模型与开口裂缝模型对桥梁位移响应幅值、瞬时频率等动力响应的影响，讨论了采用小波分析进行裂缝识别的可行性。

文中选取曲率曲线型呼吸裂缝模型，将其融入简支梁车桥耦合振动系统中，从而在运营荷载下建立考虑裂缝非线性振动的精细化数值计算模型。基于改进的HHT算法，从简支梁某点的加速度响应中提取损伤识别指标，建立裂缝的非线性损伤识别方法，并分析车重、车速等不同参数对文中方法识别效果的影响，以期采用单个传感器实现混凝土简支梁桥开裂的定位和定量识别。

1 运营环境下呼吸裂缝的模拟 1.1 车桥耦合振动系统

为准确模拟运营环境下车辆与简支梁的相互作用，将简支梁和以速度v行驶的车辆简化为图 1所示的平面车桥耦合振动系统^[14]。其中，简支梁离散为长度相等的N个单元；车辆简化为半车模型，其车体采用转动惯量为I₀、质量为m₀的刚体表示。车辆前后轮组的质量分别为m₁、m₂，车体质心距车体前后距离分别为b₁、b₂，前后轮组通过弹簧单元k_i和阻尼单元c_i(i=1, 2, 3, 4)与车体相连并与简支梁接触。

但车辆作用于桥面时，将造成桥面发生变形，故设车辆前轮、后轮组所对应的桥面接触点变形分别为y₁、y₂，并假设车轮不会脱离桥面。

根据车、桥之间位移与作用力的协调关系，可建立车桥耦合共振的动力学方程^[17]：

$ \left[ {\begin{array}{*{20}{c}} {{\mathit{\boldsymbol{M}}_{\rm{v}}}}&{}\\ {}&{{\mathit{\boldsymbol{M}}_{\rm{b}}}} \end{array}} \right]\left[ {\begin{array}{*{20}{c}} {\ddot D}\\ {\ddot Y} \end{array}} \right] + \left[ {\begin{array}{*{20}{c}} {{C_{\rm{v}}}}&{}\\ {}&{{C_{\rm{b}}}} \end{array}} \right]\left[ {\begin{array}{*{20}{c}} {\dot D}\\ {\dot Y} \end{array}} \right] + \left[ {\begin{array}{*{20}{c}} {{\mathit{\boldsymbol{K}}_{\rm{v}}}}&{}\\ {}&{{\mathit{\boldsymbol{K}}_{\rm{b}}}} \end{array}} \right]\left[ {\begin{array}{*{20}{c}} D\\ Y \end{array}} \right] = \left[ {\begin{array}{*{20}{c}} {{F_{\rm{v}}}}\\ {{F_{\rm{b}}}} \end{array}} \right], $

(1)

其中：M_v、C_v和K_v分别表示车辆的质量、阻尼和刚度矩阵；D=[d₁  d₂  d₃  d₄]^T，为半车模型的位移向量；M_b、C_b、K_b分别为简支梁的质量、阻尼和刚度矩阵；Y为简支梁各节点竖向平动和转动组成的位移向量；F_v为半车模型的外激励向量^[17]：

$ {F_{\rm{v}}} = {\left[ {\begin{array}{*{20}{c}} 0&0&{{k_3}{y_1} + {c_3}{{\dot y}_1}}&{{k_4}{y_2} + {c_4}{{\dot y}_2}} \end{array}} \right]^{\rm{T}}}, $

(2)

其中，F_b为作用在简支梁上的节点力向量，F_b=N₁f₁+N₂f₂，其中N₁、N₂均为有限元插值函数，取较为常用的三次Hermite多项式；f₁、f₂分别为前后车轮施加在桥上的作用力^[16]：

$ \left\{ {\begin{array}{*{20}{l}} {{f_1} = - {k_3}({y_1} - {d_3}) - {c_3}({{\dot y}_1} - {{\dot d}_3}) - {m_1}g - \frac{{{b_2}}}{{{b_1} + {b_2}}}mg,}\\ {{f_2} = - {k_4}({y_2} - {d_4}) - {c_4}({{\dot y}_2} - {{\dot d}_4}) - {m_2}g - \frac{{{b_1}}}{{{b_1} + {b_2}}}mg}。\end{array}} \right. $

(3)

采用Newmark-β法求解式(1)，即可得到车桥耦合振动系统的动力响应时程。

1.2 曲率型呼吸裂缝模型

采用曲率型呼吸裂缝模型，该模型通过曲率表征裂缝所引起的刚度变化，设模型表征的裂缝为切口裂缝，且在深度范围内贯通整个横截面，若结构第i单元有裂缝，则在第t时刻该单元的单元刚度k_i^b(t)可表述为^[13]：

$ k_i^b(t) = k_i^c + \frac{1}{2}(k_i^o - k_i^c)\left[ {1 + \frac{{y_i^{\prime \prime }(t)}}{{y_{i{\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\rm{max}}}^{\prime \prime } }}} \right], $

(4)

其中，k_i^c、k_i^o分别为裂缝闭合、完全张开时第i单元的单元刚度。k_i^o可采用附加柔度法^[18]计算：

$ {k_{\rm{o}}} = T{{\mathit{\boldsymbol{\tilde C}}}^{ - 1}}{\mathit{\boldsymbol{T}}^{\rm{T}}}, $

(5)

$ \mathit{\boldsymbol{T}} = {\left[ {\begin{array}{*{20}{c}} { - 1}&{ - L}&1&0\\ 0&{ - 1}&0&1 \end{array}} \right]^{\rm{T}}}, $

(6)

其中，$ \mathit{\boldsymbol{\tilde C}} $为结构混凝土开裂后的附加柔度矩阵，通过求解剪力弯矩关于附加应变能的二阶偏导进行计算，具体计算公式详见文献[18]; y″_i(t)为第i单元当前时刻的挠度曲率；y″_{i max}为整个振动过程中该单元的最大挠度曲率值，往往考虑k_i^b(t)=k_i^o的情况下，通过车辆过桥全过程计算并取最大的y″_i(t)作为初始值，再经多次迭代获得^[14]。在每时间步更新第i单元的单元刚度，则可实现呼吸裂缝的模拟。

2 裂缝非线性损伤识别方法 2.1 HHT算法基本原理

Hilbert-Huang变换(HHT)是处理非线性、非平稳信号的一种时频分析方法^[19]。该方法的基本思路是将原始信号x(t)进行经验模态分解(EMD)：

$ x(t) = \sum\limits_{i = 1}^n {{c_i}} (t) + {r_n}(t), $

(7)

其中，c_i(t)为第i个本征模态函数(IMF)分量；n为IMF分量的个数，r_n(t)为残差。

将IMF分量进行Hilbert变换：

$ \mathop {{c_i}}\limits^ \wedge (t) = \frac{1}{\pi }\int_{ - \infty }^\infty {\frac{{{c_i}(\tau )}}{{t - \tau }}} {\rm{d}}\tau , $

(8)

以c_i(t)为实部，$ \widehat {{c_i}} $(t)为虚部构造解析函数：

$ {z_i}(t) = {c_i}(t) + j\mathop {{c_i}}\limits^ \wedge (t) = {a_i}(t){{\rm{e}}^{i\varphi i(t)}}, $

(9)

其幅值函数a_i(t)和相位函数φ_i(t)分别为

$ {a_i}(t) = \sqrt {{c_i}(t) + \mathop {{c_i}}\limits^ \wedge (t)} , $

(10)

$ {\varphi _i}(t) = {\rm{arctan}} \frac{{\mathop {{c_i}}\limits^ \wedge (t)}}{{{c_i}(t)}}, $

(11)

瞬时频率表示为

$ {f_i}(t) = \frac{1}{{2\pi }}{\omega _i}(t) = \frac{1}{{2\pi }}\frac{{{\rm{d}}{\varphi _i}(t)}}{{{\rm{d}}t}}, $

(12)

RE表示取实部, 则原信号的Hilbert谱为

$ H(\omega ,t) = {\rm{RE}} \sum\limits_{i = 1}^n {{a_i}} (t){{\rm{e}}^{j\int {{\omega _i}} (t){\rm{d}}t}}。$

(13)

应指出的是，为了消除HHT变换的端点效应和模态混叠现象，采用了镜像延拓方法和EEMD方法进行改进，这2种方法的详细原理见文献[20-21]，引入这2种方法后的HHT方法为改进的HHT方法。

2.2 基于能量指标的损伤识别方法

对于运营状态下的混凝土简支梁，当车辆通过裂缝截面时，裂缝处刚度的变化将带来冲击效应，使加速度响应信号增加额外的能量。但该能量的变化较小，采用改进的HHT算法提取能量指标，实现裂缝的损伤识别，具体步骤为：

1) 能量时程曲线的计算。根据式(13)将所得Hilbert谱向时间轴投影，得到能量随时间的变化曲线h(t)：

$ h(t) = \int_0^{\omega {\rm{p}}} H (\omega ,t){\rm{d}}\omega , $

(14)

式中，ω_p为截止频率。

2) 有无损伤的判定。因裂缝处桥梁响应的能量有突变，故通过h(t)时程曲线有无明显尖峰来判定结构是否开裂。

3) 损伤定位的判别。提取损伤尖峰的最大值作为特征点，记该特征点出现的时刻为t₂，车辆入、出桥的时间分别记为t₁、t₃，则裂缝和车辆入桥端的距离l为

$ l = \frac{{{t_2} - {t_1}}}{{{t_3} - {t_1}}} \times L, $

(15)

其中，L为简支梁的长度。

因在车桥耦合振动系统中，裂缝添加于某个单元中，计算添加裂缝所在的单元号m，以便于对比：

$ m = {\rm{int}} (l/{l_e}) + 1, $

(16)

其中，int(x)为取整函数；l_e为单元长度。

4) 损伤程度识别。记能量时程曲线h(t)在t₂时刻的幅值为w_A，在t₁时刻的幅值为w_U，定义这2个时刻幅值比f为

$ f = \frac{{{w_A}}}{{{w_{{U_1}}}}} \times 100\% 。$

(17)

在同一车辆作用下，该比值越大，裂缝处响应能量值越大，即裂缝深度越深，所以用于开裂程度的判定。

3 数值算例 3.1 损伤工况

根据参考文献[16]，定义图 1所示简支梁和车辆荷载模型的具体参数，如表 1所示。

设裂缝深度系数S=a/h，其中，a为裂缝深度，h为简支梁截面高度。在不同部位添加裂缝及不同的S设置不同的损伤工况，如表 2所示。设车辆过桥的速度为2 m/s，选取简支梁跨中处计算所得加速度响应为输出信号。

3.2 损伤识别结果

按照第2节方法，计算工况1~4时的能量时程曲线如图 2所示。由图 2可知，在混凝土开裂处，能量时程曲线的幅值均出现增大现象，说明通过裂缝处幅值变化能正确识别简支梁出现的损伤。应注意的是，在车辆入、出桥时，能量幅值均较大，但在完好状态和不同损伤工况下，能量时程曲线均表现出该现象，并不影响对混凝土是否开裂的判定。

取图 2中各曲线尖峰点为特征点，计算各工况下特征点的距离l和裂缝所在单元，如表 3所示。在表 3中，计算得到的裂缝单元信息与在有限元模型中添加裂缝的单元信息完全吻合，表明文中方法能准确识别混凝土开裂位置。

计算损伤工况3、5、6、7、8下幅值比f，如图 3所示。图 3表明，随着裂缝深度系数的增大，f值呈单调上升趋势，即通过f值的大小能够有效识别该裂缝的相对损伤程度。

3.3 车辆条件影响分析

对于运营条件下的简支梁桥，过桥车辆的重量和速度具有较大的随机性，所以需分析在不同车重和车速条件下文中方法的识别效果。

3.3.1 车重影响

在裂缝位于简支梁左端l=15.5 m处，裂缝深度系数S=0.5，车速v=2 m/s的条件下，改变车重参数设置损伤工况9~11，如表 4所示。

根据文中方法，计算以上4种工况下裂缝所在单元和幅值比，如表 4所示。由表可知，在不同车重情况下，计算得到裂缝所在单元号m均为16，与裂缝预设位置吻合；计算的能量幅值比相差很少，表明车重即使变化超过100%，也不影响对裂缝损伤的定位和程度的识别。

3.3.2 车速影响

为考查车速对开裂识别的影响，选定工况2的参数，分别取车速为2、4、6、8、10、15、20 m/s，计算得到的损伤位置l和幅值比f，如图 4所示。

从图 4(a)中可知，除在车速v=15 m/s条件下识别得到的损伤位置与预设损伤单元有少许偏差外，其余各车速条件下损伤定位结果均在预设损伤单元内，表明车速对损伤定位的影响较小。

从图 4(b)中可知，随着车辆行驶速度的增大，裂缝损伤的幅值比f也相应增大，说明车速对损伤程度识别的影响较大。这主要是由于过大的车速会造成能量尖峰(也包含特征点)之间的重合，从而无法找到能量时程曲线上的尖峰特征点，例如，选取车速v为4、8、15 m/s，过裂缝桥段能量时程曲线如图 5所示。

应说明的是，虽然车速对损伤程度识别有较大影响，但从图 2可知，每辆车入桥、出桥均有明显的特征，可利用能量时程曲线计算得到车辆过桥的车速v。在长期监测中可获得大量车辆过桥的数据，选用车速相差较少数据进行损伤程度的识别。

4 结论

建立了运营荷载下考虑裂缝呼吸效应的简支梁精细化数值计算模型，基于改进的HHT方法计算车辆过桥的能量时程曲线。提出了基于能量指标的混凝土简支梁桥裂缝位置和程度的识别方法。模拟算例表明：

1) 文中损伤识别方法采用单个传感器即可准确定位裂缝的损伤位置以及区分损伤程度的相对大小；

2) 不同车重对裂缝部位及损伤程度识别的影响较小；

3) 车速大小不影响裂缝部位的识别精度，但对损伤程度识别有较大影响，这一缺点可通过选用相同或相近车速条件下的信号在一定程度上予以克服。