0 引言

桥梁工程是一项关系国家民生发展的重要工程，伴随着中国桥梁建设的飞速发展，截止2023年末，全国公路桥梁达107.93万座、总里程达9 528.82万米^[1]。与此同时，服役长达10 a以上的桥梁数量也在逐年攀升，特别是对于量大面广的中小跨径公路桥梁，如不加以及时监测与检测，极有可能引发重大的安全事故^[2]。目前，针对中国公路桥梁检测主要分为无损检测、静载检测和动载检测，其主要形式包括外观损伤检测、内部缺陷检测、力学性能检测等^[3]。此种检测方法虽能达到对桥梁全面检测的目的，但所耗费的人力物力巨大，且检测的周期长，对交通的影响大，既不经济也不便捷。自Yang^[4]于2004年首次提出通过移动安装有传感器的检测车响应识别桥梁模态参数的间接测量法概念后，国内众多学者纷纷投入到桥梁间接测量技术的研究与应用中。采用这类方法，有望实现对中小跨径桥梁整体状态的快速评估，大幅提高检测的效率和极大降低成本^[5]。

如何通过获取的车体动力响应提取桥梁模态参数和损伤信息是研究的一大热点^[6]。2004年，Yang^[4]首次尝试将装有加速度传感器的检测车辆作为桥梁信息载体，提出利用获取的竖向加速度车体响应间接提取桥梁基频的方法，数值模拟验证了方法的可行性。此后，Yang^[7]指出，路面不平整度会对车桥耦合桥梁检测方法的频率识别产生不利影响，试图使用两辆连续相同的车辆来减少影响。Yin^[8]建立了斜拉桥-检测车的车桥耦合模型, 通过数值模拟的方式成功获取了车体的位移响应，并通过对位移响应进行POD分解，提取损伤特征分量，较为准确地判断了结构的损伤位置。Li^[9]提出一种基于广义模式搜索算法(GPSA)的桥梁模态参数间接识别方法。此方法通过采用GPSA算法对车桥耦合系统下的车体响应进行优化修正成功识别到了桥梁的刚度和一阶频率，数值模拟采用对比的方式验证了所提方法的优势。Zhang^[10-11]在总结前人研究后，提出直接采用短时傅里叶变换对获取的车体响应进行处理，成功获取了桥梁的模态振型，最后以模态振型平方(MOSS)作为损伤判定指标识别桥梁损伤，但同时也提出，检测车移动速度可能会对识别结果产生一定的影响。文献[12]以公路混凝土桥梁为研究对象，建立了双轴测试车和桥梁的车桥耦合有限元模型，采用EMD分解和Hilbert谱分析对车体信号处理后成功识别了结构的损伤。Yang^[13]指出桥面粗糙度影响测试车对桥梁频率的提取，提出利用随机子空间法对测试数据进行处理，数值模拟研究表明此方法可以有效抑制路面粗糙度对频率提取的影响。Tan^[14]利用单辆安装传感器的检测车提取桥梁阻尼，但是需要预先知道桥梁的振型。Zhang^[15-16]首次提出单轴测试车接触点响应概念，并从理论上推导了车桥接触点响应的理论解。通过数值模拟和现场试验验证表明：相比于车体响应，通过车桥接触点响应能够获取更多桥梁的模态信息。Nayek^[17]提出一种采用输入状态估计方法提取多自由度测试车接触点响应的方法，随后在考虑桥面不平整度的情况下识别了桥梁的振型。阳洋^[18-19]提出采用两辆参数一致的单轴测试车行驶过桥检测桥梁损伤的方法，该方法通过对前后测试车响应相减来消除路面粗糙的影响，数值模拟和实桥试验验证了方法的可行性。Locke^[20]考虑温度、车速、交通和表面粗糙度等因素引起的车体动力响应的影响。采用卷积神经网络对车辆加速度数据进行频域分析，识别了桥梁单元损伤，数值模拟和试验研究的结果证明了卷积神经网络在一系列现实噪声信号中检测桥梁损伤的能力。Yang^[21]采用单轴试验车对纵向两端简支、横向自由的板式桥梁进行损伤检测，提出采用均匀平移响应(UTR)来定位板式桥梁损伤。亓兴军^[22]利用有限元软件建立了双轴半车模型，采用希尔伯特变换对测试车响应进行处理，最后成功获取连续梁桥模态振型并以此识别桥梁挠度。数值研究表明：该方法对测试车车速不敏感，但路面粗糙度和环境噪声会对识别结果产生影响。

奇异谱分析(SSA)方法是近年来兴起的一种研究非线性时间序列数据的强大方法。此方法通过对时间序列分解和重构达到提取信号中主要信息的目的^[23]。目前，SSA方法已应用于多个方面，如信号中噪声的分解和消除、数据预测和故障诊断等^[24]，但将该方法应用于结构的动力响应信号进行损伤诊断的相关研究并不多，本研究尝试将奇异谱分析方法应用于桥梁结构的损伤识别中，结合车桥耦合间接测量技术，达到快速高效的诊断目的。本研究首先从理论上推导了车桥耦合理论和奇异谱分析理论，初步说明方法的可行性，接着采用数值算例对本研究方法流程进行了详细地展示，最后进行了各种影响参数的分析，进一步说明本研究所提方法的可行性。

1 理论基础 1.1 车桥耦合理论

根据间接测量下的车桥接触点理论^[15]，桥梁被简化为欧拉-伯努利线弹性简支梁，梁长为L，单位长度质量为m，抗弯刚度为I_E。检测车被简化了弹簧-质量块体系，集中质量为m_v，质量块与简支梁用一根刚度为k_v的弹簧连接。检测车以速度v匀速行驶过桥梁，激励梁振动并产生响应，检测车车体竖向位移响应为y_v，接触点竖向位移响应为u_c，梁振动的竖向位移响应为u_b。车桥耦合简化模型如图 1所示。

不考虑检测车阻尼和桥梁阻尼比的情况下，车辆和桥梁的运动方程分别为：

$ m_{\mathrm{v}} \ddot{y}_{\mathrm{v}}+k_{\mathrm{v}}\left(y_{\mathrm{v}}-u_{\mathrm{c}}\right)=0, $

(1)

$ \bar{m} \ddot{u}_{\mathrm{b}}+E I u^{\prime \prime \prime \prime}=\left[k_{\mathrm{v}}\left(y_{\mathrm{v}}-u_{\mathrm{c}}\right)-m_{\mathrm{v}} g\right] \delta(x-v t), $

(2)

式中，(¨)表示对时间t求二阶导数；()″″表示对位置x对四阶导数；g为重力加速度；δ为狄克拉函数，其余参数同前。需要注意的是，车体的竖向位移响应是从车体平衡位置处算起。

根据模态叠加原理^[25]，梁的竖向位移响应可表示为：

$ u_{\mathrm{b}}=\sum\limits_n\left[q_{\mathrm{bn}}(t) \sin \left(\frac{n \mathtt{π} x}{L}\right)\right], $

(3)

式中，q_bn(t)为模态坐标，梁的振型ϕ_n(x)=$\sin \left(\frac{n \mathtt{π} x}{L}\right)$，即梁的竖向位移响应可由一系列的模态坐标和振型乘积叠加得到。

根据所设参数和结构动力学相关理论，车辆频率ω_v和桥梁频率ω_bn可表示为：

$ \omega_{\mathrm{v}}=\sqrt{\frac{k_{\mathrm{v}}}{m_{\mathrm{v}}}}, $

(4)

$ \omega_{\mathrm{bn}}=\frac{n^2 \mathtt{π}^2}{L^2} \sqrt{\frac{I_{\mathrm{E}}}{\bar{m}}}。$

(5)

在车辆运动过程中，接触点竖向位移响应存在如下关系：

$ u_{\mathrm{c}}=\left.u_{\mathrm{bn}}(x, t)\right|_{x=v t}。$

(6)

根据前述所列公式，推导简化可得接触点处的竖向位移响应的表达式：

$ u_{\mathrm{c}}=\sum\limits_{n=1}^M\left[\frac{\Delta_{stn}}{1-S_{\mathrm{n}}^2}\left(\sin \frac{n \mathtt{π} v t}{L}-S_{\mathrm{n}} \sin \omega_{\mathrm{bn}} t\right) \sin \frac{n \mathtt{π} v t}{L}\right]。$

(7)

式中变量：

$ \Delta_{s t n}=-\frac{2 m_{\mathrm{v}} g L^3}{E \operatorname{In}^4 \mathtt{π}^4}, \quad S_{\mathrm{n}}=\frac{n \mathtt{π} v}{L \omega_{\mathrm{bn}}}。$

从式(7)中可以看出，车桥接触点响应仅含有与桥梁频率相关的信息，不含有车辆频率，利用车桥接触点响应对桥梁进行分析，可完全避免车桥频率对信息的干扰，大大提高了识别的精度^[26]。实际操作中，可利用放置于检测车上的传感器采集车体响应，然后采用常用的数值计算方法(newmark-β法、中心差分法等)进行处理便可求取接触点响应。

1.2 奇异谱分析

奇异谱分析^[27]在应用于含结构特征信息的时间序列时，能够有效地将含有结构特征分量和噪声及其他无用分量进行分离，适合于对结构的响应信号进行特征提取和损伤诊断。SSA分析包括分解和重构两个过程，相关理论如下。

1.2.1 分解

通过车桥耦合理论反演车体信号得到了车桥接触点信号，其一维时间序列为U_N=[u₁, u₂, …, u_N－1, u_N]，N为时间序列的采样点总数。定义嵌入窗口长度为L，L为1至N范围内的任意整数。将时间序列U_N嵌入到K=(N－L+1)个向量中，得到轨迹矩阵为U：

$ \boldsymbol{U}=\left[\boldsymbol{U}_1, \boldsymbol{U}_2, \cdots, \boldsymbol{U}_{N-1}, \boldsymbol{U}_N, \right]=\left(\begin{array}{ccc} u_1 & \cdots & u_K \\ \vdots & \ddots & \vdots \\ u_L & \cdots & u_N \end{array}\right)。$

(8)

由式(8)可得，轨迹矩阵U实质上是一个Hankel矩阵。

定义矩阵C=U^TU，矩阵C即为轨迹矩阵U的协方差矩阵，对协方差矩阵C进行奇异值分解，分解形式为：

$ \boldsymbol{C}=\boldsymbol{P} {\varSigma} \boldsymbol{V}^{\mathrm{T}}, $

(9)

式中，Σ由L个降序排列的非负特征值λ₁, λ₂, …, λ_L－1, λ_L组成的对角矩阵。P为由特征向量P₁, P₂, …, P_L－1, P_L组成的矩阵，设轨迹矩阵的秩(非零特征值的个数)为d，轨迹矩阵U可以表示为：

$ \boldsymbol{U}=\boldsymbol{U}_1+\boldsymbol{U}_2+\cdots+\boldsymbol{U}_d, $

(10)

其中，$\boldsymbol{U}_i=\sqrt{\lambda_i} \boldsymbol{P}_i \boldsymbol{V}_i^{\mathrm{T}}, \quad \boldsymbol{V}_i=\boldsymbol{U}^{\mathrm{T}} \boldsymbol{P}_i / \sqrt{\lambda_i}$。

为了确定各个初等矩阵U_i在轨迹矩阵U中占据的比重，定义贡献率系数为：

$ a_i=\lambda_i / \sum\limits_{i=1}^d \lambda_i。$

(11)

由式(11)可得，特征值λ_i越大，系数a_i就越大，其对应的初等矩阵U_i所包含原始信号中结构特征信息就越多。

1.2.2 重构

按照式(11)所计算的贡献率大小将上述矩阵U_i进行重新分组，然后将每个组内的矩阵相加得到新的矩阵S^[28]。由对角平均化公式可以将上述时间序列成分重构成时间长度为N的各个分量，令S是一个维数为L×K的矩阵，其中对角平均化的计算公式如下：

$ s_k= \begin{cases}\frac{1}{k} \sum\limits_{m=1}^k s_{m, k-m+1}^* & 1 \leqslant k \leqslant L^* \\ \frac{1}{L^*} \sum\limits_{m=1}^{L^*} s_{m, k-m+1}^* & L^* \leqslant k \leqslant K^* \\ \frac{1}{N-K+1} \sum\limits_{m=k-K^*+1}^{N-K^*+1} s_{m, k-m+1}^* & K^* \leqslant k \leqslant N\end{cases}, $

(12)

式中，L^*=min(L, K), K^*=max(L, K)，运用式(12)即可将矩阵S还原为对应的时间序列[s₁, s₂, …, s_N]。由此，应用对角平均公式对式(10)中子矩阵U_i还原可得到d个重构序列RC，按照对应的特征值由大到小降序排列对应的重构序列RC₁, RC₂, …, RC_d。一般认为取前几项特征值较大的对应重构序列即可反应结构的特征信息，本研究取前两项RC₁, RC₂进行说明。

1.3 移动平均滤波

1.2节对奇异谱分析进行了理论推导，但得到的前几项重构序列往往具有一定的趋势和含有一定的周期成分，而周期成分可能导致识别效果不佳，因此采用文献[29]建议的移动平均滤波器(MAF)^[29]对重构序列进行滤波和平滑处理。其理论公式为：

$ \overline{\boldsymbol{RC}}_n(i)=\frac{1}{T_k f_s}\sum\limits_{j=i-k / 2}^{j=i+k / 2} \boldsymbol{R} \boldsymbol{C}_n(i+j) \quad(n=1, \quad 2), $

(13)

式中，T_k为需要去除的周期成分；f_s为采样频率。

1.4 桥梁损伤识别流程

本研究对车辆接触点响应进行奇异谱分析的桥梁损伤快速识别流程如图 2所示。

2 数值算例 2.1 数值模型

按照实际桥梁参数范围建立车桥耦合模型，如图 1所示，具体参数设置：桥梁被简化为跨长L=30 m的简支梁，弹性模量为E=3.45×10¹⁰ N/m²，单位长度质量m=12 000 kg/m，梁的截面惯性矩I=3.00 m⁴，桥梁阻尼比取ξ=0.01。检测车为简化的单自由度质量弹簧模型，质量m_v=5 000 kg，刚度k_v=500 kN/m，不考虑检测车阻尼，检测车从左端支座处以速度v=1 m/s匀速行驶过桥梁。采用MATLAB软件建立车桥耦合VBI单元进行数值分析^[4]，采用文献[18-19]所建议的采样频率100 Hz。为便于设置桥梁损伤，将桥梁等间距划分为15个单元并依次编号，对目标单元刚度进行折减以模拟桥梁损伤，桥梁单元划分示意图如图 3所示。

首先对检测车采集的车体位移响应进行反算得到车桥接触点位移响应，然后以车桥接触点位移响应作为分析信号进行奇异谱分析，得到贡献率最大的前两项重构时间序列，最后通过序列曲线进行桥梁损伤快速识别。下面将以数值模型参数为例，具体展示数值分析的步骤。

2.2 基于奇异谱分析的桥梁损伤识别

为便于对比说明，设置无损工况和有损工况，无损工况即为桥梁单元截面抗弯刚度无折减，有损工况为：对单元⑥设置30%的损伤，即该单元的抗弯刚度I_E6=0.7×I_E，其余单元抗弯刚度无折减。按照图 2识别流程和2.1节单元划分进行数值分析，得到两种工况下的车体位移响应和反算的车桥接触点位移响应如图 4，5所示。

从图 4和图 5可以看到，桥梁在发生损伤后，其响应曲线相比于无损情况下存在偏差，但这样的对比并不能直观准确地定位桥梁的损伤位置。因此，利用奇异谱分析分别对两种工况下的接触点位移响应进行处理，得到前两项重构序列RC₁, RC₂曲线，如图 6所示。

图 6中得到有损工况下的RC₁曲线较无损工况有明显的偏差，但通过图 6(a)并不能直观地反映其损伤位置，这是由于重构序列RC₁所对应的是奇异值分解下的最大特征值，它对原始序列的贡献率是最大的，因此该序列包含了原始序列的整体趋势，也就无法反映其局部的特征^[30]。但在图 6(b)中，重构序列RC₂凸显了原始序列的局部特征，相较于无损工况下的曲线，可以明显地看到有损工况RC₂曲线在单元⑥处发生了突变，准确地定位了损伤。

对于重构序列RC₂后续序列的分析，本节利用上述损伤工况，采用数值模拟的方式对RC₃, RC₄做进一步对比分析，结果如图 7所示。

由图 7可以看出，重构序列RC₃, RC₄曲线甚至不能得到类似RC₁, RC₂光滑曲线，无法反映桥梁损伤。这是由于前两项重构序列RC₁, RC₂占据绝大部分桥梁特征信息，且文献[31-32]指出靠后的重构序列所含噪声成分更多，因此，为限于篇幅和便于对比，后续分析将仅展示重构RC₂曲线。

但图 6(b)中RC₂曲线并不平滑，若曲线周期性过大将给识别带来一定的困难，因此，采用1.3节建议的移动平均滤波器对曲线进行平滑处理，得到图 8。从图 8可以看到，平滑后的曲线更加准确地定位了桥梁的损伤，突出了本研究方法的可行性。

3 影响参数分析

在实际桥梁损伤识别过程中，许多影响因素直接或间接地影响到动力信号的采集，进而影响损伤识别的精度。为了进一步验证本研究方法的可行性，将对车桥耦合分析下常考虑的影响因素(桥梁阻尼比、桥梁单元划分、行车速度、环境噪声、不同损伤程度和位置)等进行模拟分析。

3.1 桥梁阻尼比的影响

对桥梁阻尼比设置为0.000(无阻尼)，0.005，0.010，0.020，这4组阻尼比工况，检测车和桥梁其余参数同2.1节，按照图 2和2.2节流程识别桥梁损伤，损伤工况同2.2节，得到桥梁损伤下的重构序列RC₂曲线，如图 9所示。

图 9所示结果表明，不同桥梁阻尼比下重构序列RC₂曲线基本一致，说明本研究所提方法对桥梁阻尼比不敏感。

3.2 桥梁单元划分的影响

桥梁单元划分的数目会直接影响损伤识别的精度，因此桥梁单元的划分对于最终的损伤判定具有重要意义。保持桥梁跨长30 m不变，桥梁单元划分数目和损伤位置相关信息列于表 1中。检测车和桥梁其余参数同2.1节，按照图 2和2.2节流程识别桥梁损伤，得到不同工况桥梁损伤的重构序列RC₂曲线，如图 10所示。

从图 10可得出，在不同桥梁单元划分数的情况下，本研究所提方法均能准确识别损伤单元位置。但需要注意，实际桥梁单元划分数目不宜过少也不宜过多。如果单元数过少，对桥梁损伤的判断会十分粗糙；如果单元数过多，不仅对测试分析时间略有增加，而且过于精细的划分会导致测试车在每个单元采集的信号点减少，进而导致获得桥梁信息减少并且受噪声干扰增大，增加识别的误差，因此对于实际工程应用是不适用的。需要注意的是，本研究所提方法主要针对中小跨径的简支梁桥进行了分析，因此综合上述图表分析结果，针对中小跨径桥梁(30 m左右)建议按照长跨比(单元长度与跨长比值)为0.07左右划分单元数。

3.3 车辆行驶速度的影响

对检测车车速设置为1，3，5，8，10 m/s，5组车速工况，检测车和桥梁其余参数同2.1节，按照图 2和2.2节流程识别桥梁损伤，损伤工况同2.2节，得到不同车速工况桥梁损伤的重构序列RC₂曲线，如图 11所示。

从图 11可以得到，当车速8 m/s及以下，重构序列RC₂曲线均在损伤单元处发生突变，成功识别了桥梁损伤，表明本研究所提方法能较好适应较大范围的检测车速度。但当检测车车速达到10 m/s时，其损伤位置处曲线突变不明显，这是因为当采集仪器采样频率固定时，车速过快会导致采样点数减少，获取的桥梁信息也相对稀少，因此不建议检测车车速过快。

3.4 环境噪声的影响

实际测量并不能消除环境噪声的干扰。因此，需要考虑噪声的影响，本节采用信噪比(R_SN)表示噪声水平，R_SN表示为：

$ R_{\mathrm{SN}}=10 \times \lg \left(\frac{\frac{1}{N} \sum\limits_{i=1}^N y_i^2}{\frac{1}{N} \sum\limits_{i=1}^N \delta_i^2}\right), $

(14)

式中，N为检测车行驶过桥所采集的数据点总数；y_i为测试车在第i时刻采集的含有噪声的位移响应数据；δ_i为测试车在第i个时刻对应的噪声值。信噪比R_SN的单位为dB，其值越小，表示受噪音的影响越大。

对检测车的位移响应添加白噪声来模拟含噪信号，其公式为：

$ x_{\text {noise }}=x+E_{\mathrm{p}} \sigma(x) N_{\text {noise }}, $

(15)

式中，x_noise为含有噪声的检测车位移响应；x为检测车原始位移响应；E_p为加入噪声的程度；σ(x)为检测车原始位移响应的标准差；N_noise为高斯白噪声。其中噪声程度E_p=10%，E_p=3.16%，E_p=1%分别对应信噪比R_SN=20 dB，R_SN=30 dB，R_SN=40 dB。在数值模拟部分，具体加噪的设置方法是采用Matlab软件中的“awgn”函数实现。

对检测车位移响应设置无噪音，20，30，40 dB，4组噪声工况，检测车和桥梁其余参数同2.1节，按照图 2和2.2节流程识别桥梁损伤，损伤工况同2.2节，得到不同噪声工况桥梁损伤的重构序列RC₂曲线，如图 12所示。

从图 12可以得出，不同的信噪比R_SN情况下，本研究所提方法均有效识别了桥梁损伤单元位置，且在信噪比40 dB工况下识别的重构序列RC₂曲线基本与无噪工况保持一致。可以看到，在信噪比20 dB和30 dB工况下，虽然重构序列RC₂曲线在其他无损位置出现了不平滑的现象，但在损伤位置仍有较大的突变，进一步突出了本研究所提方法具有较强的抗噪能力，后期可针对环境噪声影响下的曲线平滑处理进行更加深入的研究。

3.5 不同损伤程度和位置分析

前述已针对桥梁单损伤的情况，验证了所提方法在桥梁阻尼比、桥梁单元划分、行车速度、环境噪声等影响下良好的损伤识别效果。为进一步评估所提方法在复杂情况下对桥梁损伤的检测能力，对桥梁不同位置处的单元设置不同损伤程度，具体设置如表 2所示。

根据表 2设置情况，工况1~4为单损伤工况，工况5~7为双损伤工况。检测车和桥梁其余参数同2.1节，按照图 2和2.2节流程识别桥梁损伤，得到不同损伤工况下桥梁损伤的重构序列RC₂曲线，如图 13，14所示。

从图 13和图 14可以看到，不论是单损伤工况还是多损伤工况，本研究所提方法均能准确地识别桥梁单元损伤位置。同时也可以看到，当损伤程度增大时，曲线在损伤位置的突变更加明显，突出了本研究所提方法的适用性。

4 结论

将奇异谱分析理论应用于车桥耦合间接测量技术中，提出了基于奇异谱分析的车桥耦合桥梁损伤识别方法。通过数值算例对所提出方法的可行性和适用性进行了验证，详细分析了多种参数对方法适用性的影响，得到以下结论。

(1) 通过理论推导、数值模拟研究分析，成功验证了本研究方法能很好地识别桥梁损伤位置。由于本研究所提方法仅需对检测车时域响应进行处理，因此该方法操作简便，识别效率高，可为间接测量技术的桥梁结构健康监测和检测提供新的思路和方法。

(2) 针对各种影响参数的数值分析可以得到：本研究所提方法对桥梁阻尼比、桥梁单元划分、行车速度不敏感，即使在较强的环境噪音下仍能很好地识别桥梁单元的损伤，体现了本研究方法在各种因素影响的适用性和较强的鲁棒性。但为了满足实际需要和获取更多桥梁特征信息，针对中小跨径桥梁(30 m左右)建议单元长度与跨长比值在0.07为宜，同时车速不宜超过8 m/s。

本研究方法目前只针对单跨的简支梁桥进行了理论推导和数值分析，同时本研究在奇异谱分析研究过程中，主要以重构序列RC₂曲线进行了不同参数下的研究，并证明了其在对车桥耦合作用下桥梁损伤识别的能力，同时后续会对进一步提取RC成分信息中所包含局部特征信息和不同桥型做进一步深入的研究。