基于少量静态测量数据的结构损伤定位与程度识别方法

摘要

本发明公开了一种基于少量静态测量数据的结构损伤定位与程度识别方法，包括以下步骤：1）根据指标D

说明书

技术领域

本发明涉及结构健康检测损伤识别技术领域，更具体地，涉及一种利用少 量静态测量数据对结构进行损伤检测的方法。

背景技术

随着我国科学技术的发展，各种大型复杂结构随处可见。从地上的建筑、 桥梁、高铁，到水上的钻井平台，乃至太空中的航天结构等等，都属于所需要 考虑到的结构范畴。然而这些结构在其服役过程中，随着时间的推移，不能避 免地在某些构件中发生局部损伤。如果对这些损伤置之不理，任其发展扩大， 最终将导致结构的毁坏。结构的毁坏不但会对国家带来经济损失，严重的更会 对人民造成生命威胁。因此，为了保证各种结构在使用期间的安全性和可靠性， 结构损伤识别与健康监测技术引起了广泛的关注。通过对结构进行监测和诊断， 及时发现局部损伤，对结构进行安全评估，已成为必然要求。

结构损伤检测的方法多种多样。从最原始的用肉眼观察到现在用高科技仪 器，如超声波、红外线等，这些技术已广泛地应用于工程实际之中。在各种各 样的方法中，基于结构振动特性的损伤识别技术得到了最广泛的研究。基于结 构振动特性的损伤识别技术的基本理论是：局部损伤必然引起结构动力特征变 化。如果按照所用数据分类，目前这类方法还能可以细分为两类：静态识别技 术和动态识别技术。

对于动态识别技术，就是应用动力系统中的各种参数来进行损伤识别。目 前主要技术方向有频域和时域两种。频域方向的技术主要应用到动力系统中的 频率、振型、模态曲率或频响函数等；时域方向的技术主要应用到动力系统的 速度、加速度、位移响应。动态识别技术的优势在于能用于识别的参数非常多， 方法比较丰富。而其劣势在于测量数据受噪声影响较大，所应用数据量较大等。

对于静态识别技术，就是应用经历测量数据来进行损伤识别。由于静力平 衡方程只涉及结构刚度，因此识别目标比较明确。而且静力测量数据比较容易 获取，测量误差很小。但是就目前的静力识别技术，大多存在算法过程复杂、 计算量大等问题。

文献“基于静力响应的桥梁结构损伤识别（国外建材科技，2006，27（2）， 105～107）”介绍了一种应用静力测量数据的损伤识别方法。该方法通过一步求 解损伤控制方程来得到识别结果，识别偏差较大。而且由于没有进行损伤位置 预测，造成未知量较多，求解方程时计算量大的问题。而文献“基于振动响应 的耦合杆系统损伤识别（中山大学学报，2012，51（5），31～35）”提出了一种 基于振动响应的模型修正方法。该方法通过循环迭代来对损伤参数进行识别， 能得到较好的识别结果。然而由于应用到时域响应数据，要求测量一定时间内 的振动数据，所以所测量的数据量相对较多。同样的，由于没有对损伤位置先 做预测，会造成识别效率较低的问题。

发明内容

为了克服现有技术的不足，本发明提供了一种使用少量静态测量数据的结 构损伤识别方法，该方法检测所需数据量少，效率高，且具有较高的精度。

为了实现上述目的，本发明的技术方案为：

一种基于少量静态测量数据的结构损伤定位与程度识别方法，基于残余应 变能的概念，通过有限元模型缩聚方法减少需测量数据，从而进行可疑单元的 预测；然后应用基于静响应灵敏度的模型修正方法，对可疑单元的损伤参数进 行识别，以精确检测出损伤位置和损伤程度。

上述方法可以简单概括为包括以下步骤：1）根据指标D_i选取加载方案；2） 计算各单元的残余应变能（RSE）值；3）比较各RSE的值，选择RES值较大 的单元定义为有可能出现损伤的可疑单元；4）基于可疑单元的损伤参数{α}， 计算灵敏度矩阵[S]以及待修正模型和真实模型的响应差值{δU_m}；5）求解方 程[S]{δα}={δU_m}；6）更新可疑单元的损伤参数α^it+1=α^it+δα^it；7）若未达到 精度要求则回到4）循环迭代，否则输出损伤参数作为识别结果。该方法定义 了残余应变能的概念，能对损伤位置进行预测，减少需要识别的参数数量。然 后应用静态响应灵敏度的方法进行模型修正。所需数据量少，识别效率和精度 较高。

上述方法具体可以包括以下步骤：

1）用有限元方法将待检测结构进行简化建模，并把结构划分为nel个单元；

2）基于应变能的概念，选择静载施加方案；所使用的应变能准则为：

$D_i=\frac{1}{\frac{1}{\mathrm{nel}}\underset{j=1}{\overset{\mathrm{nel}}{\Sigma }}{\mathrm{SE}}_{\mathrm{ij}}}\sqrt{\underset{j=1}{\overset{\mathrm{nel}}{\Sigma }}({\mathrm{SE}}_{\mathrm{ij}}-\frac{1}{\mathrm{nel}}\underset{j=1}{\overset{\mathrm{nel}}{\Sigma }}{\mathrm{SE}}_{\mathrm{ij}})}$

其中，是第j个单元在第i组预设外载作用下的单元应变 能；其中下标i表示预设加载方案的编号，U_i为完好模型在第i个加载方案下的 位移向量，而则是完好状态下第j个单元的单元刚度矩阵；其中D_i越小，说 明该加载方案越适合用于该结构以进行识别。

3）选择测点，其中力的作用点必须为测点；基于选择的测点和加载方案， 应用有限元模型缩聚方法，进行有限元模型缩聚；

4）对损伤状态下的结构加载静力，测量各测点静态位移数据，得到测量的 位移向量U_m；同时计算完好模型受到同样静载，在测点位置所产生的位移向量 U，计算各个单元的残余应变能RSE的值，残余应变能RSE的计算公式为：

${\mathrm{RSE}}_j^R=\delta U_m^T{K}_j^{\mathrm{eR}}\delta U_m,j=\mathrm{1,2},···,\mathrm{nel}$

其中上标T和R分别表示矩阵的转置和有限元模型的缩聚状态，下标j代 表单元编号，δU_m=U-U_m；

然后，通过比较各RSE的值，选择RES值较大的单元定义为有可能出现 损伤的可疑单元；

5）对步骤4）中得到的可疑单元的损伤参数α进行进一步的精确识别；计 算得到静态响应灵敏度为：

$\frac{\partial U_d}{\partial {\alpha }_j}=-K_d^{-1}{K}_j^e{U}_d$

其中，U_d和K_d分别为迭代识别过程中，根据损伤参数的值计算得到的位 移向量和结构总体刚度矩阵；

6）建立模型修正方程[S]{δα}={δU}，其中[S]是由静态响应灵敏度组成的 灵敏度矩阵，δU=U_m-U_d；利用Tikhonov正则化方法和L-curve方法求解该模 型修正方程，得到δα的值；

7）可疑单元的损伤参数更新为α=α+δα，重新计算待修正模型的刚度矩 阵K_d；

8）若损伤参数改变量δα没达到预设精度要求，则回到步骤5）继续迭代； 否则，以步骤7）中得到的损伤参数α为最终识别结果。

上述步骤6）中的模型修正方程为病态方程，最小二乘法不能满足对此类 方程的求解，为了求解该修正方程，还可采用以下方法：迭代改善发、加权迭 代改善法、共轭梯度法，条预优法和遗传算法及其他进化算法等来解得损伤参 数增量δα的值。

上述所述α=1表示没有损伤，α=0表示彻底损坏，通过查看损伤参数的 编号，对应有限元模型单元编号，得到结构损伤的精确位置以及损伤程度。此 处的有限元建模的时候就会有单元编号，每个单元对应有一个损伤参数，所以 单元编号和损伤参数编号是一一对应的。

本发明的有益效果在于：本发明应用静态测量数据作为输入数据，较于文 献“基于振动响应的耦合杆系统损伤识别（中山大学学报，2012，51（5），31～35）” 使用动力测量数据，静力试验在现实工程中比较容易实现，并且静态数据的测 量误差相对较小。此外，与一步求解损伤控制方程的文献“基于静力响应的桥 梁结构损伤识别（国外建材科技，2006，27（2），105～107）”相比较，运用模 型修正方法，提高识别精度。最后，本发明通过两步法，即先预测可疑单元再 做精确识别，提高了识别效率与识别精度。

附图说明

图1为本发明实施例1的简支梁结构示意图。

图2为本发明实施例1的两种待选外载方案示意图。

图3为本发明实施例1预测可疑单元的结果示意图：各个单元的RSE值。

图4为本发明实施例1的最终识别结果示意图：各单元的损伤程度。

图5为本发明实施例2的连续梁结构示意图。

图6为本发明实施例2预测可疑单元的结果示意图：各个单元的RSE值。

图7为本发明实施例2的最终识别结果示意图：各单元的损伤程度。

图8为本发明的流程图。

具体实施方式

下面结合附图对本发明做进一步的描述，但本发明的实施方式并不限于此。

本发明的基于少量静力测量数据的损伤识别方法，主要特征为两步。首先 定义残余应变能的概念，通过有限元模型缩聚方法减少需测量数据，从而进行 可疑单元的预测；然后应用基于静响应灵敏度的模型修正方法，对可疑单元的 损伤参数进行识别，以精确检测出损伤位置和损伤程度。

1）静载的选取

假设待检测结构被离散为nel个单元，则完好状态下，在第i组预设外载作 用下，其静力平衡方程可以表达为：

KU_i=F_i

若表示第j个单元在总体自由度下的单元刚度矩阵，则第 j个单元在第i组预设外载作用下的单元应变能为：

${\mathrm{SE}}_{\mathrm{ij}}=\frac{1}{2}{U}_i^T{K}_j^e{U}_i$

而选取力的指标定义为：

$D_i=\frac{1}{\frac{1}{\mathrm{nel}}\underset{j=1}{\overset{\mathrm{nel}}{\Sigma }}{\mathrm{SE}}_{\mathrm{ij}}}\sqrt{\underset{j=1}{\overset{\mathrm{nel}}{\Sigma }}({\mathrm{SE}}_{\mathrm{ij}}-\frac{1}{\mathrm{nel}}\underset{j=1}{\overset{\mathrm{nel}}{\Sigma }}{\mathrm{SE}}_{\mathrm{ij}})}$

D_i越小，说明该加载方案越适合用于该结构以进行识别。

2）残余应变能预测可疑单元

如果原来的有限元模型的自由度为ndof，测量数据量为m，则未测量的自 由度数量为s=ndof-m。其中力的作用点必须为测点，重新排列自由度的顺序， 则静力平衡方程可以重新表达为：

$\left(\begin{array}{cc}{K}_{\mathrm{mm}}& K_{\mathrm{ms}}\\ K_{\mathrm{sm}}& K_{\mathrm{ss}}\end{array}\right)\left(\begin{array}{c}{U}_m\\ U_s\end{array}\right)=\left(\begin{array}{c}F\\ 0\end{array}\right)$

引入转换矩阵：

$T=\left(\begin{array}{c}I\\ -K_{\mathrm{ss}}^{-1}{K}_{\mathrm{sm}}\end{array}\right)$

其中I为单位对角矩阵。

则缩聚的有限元模型为：

K^RU_m=F

其中K^R=T^TKT。如此类推，缩聚下的单元刚度 $K_j^{\mathrm{eR}}=T^T{K}_j^{e}T(j=\mathrm{1,2},···,\mathrm{nel})$

在缩聚模型下，计算各个单元的残余应变能值。残余应变能的计算公式定 义为：

${\mathrm{RSE}}_j^R=\delta U_m^T{K}_j^{\mathrm{eR}}\delta U_m(j=\mathrm{1,2},···,\mathrm{nel})$

其中δU_m=U-U_m，U_m为由测量数据组成的位移向量，U为完好模型受 到同样静载，对应地在测点位置所产生的位移向量。在柱状图中，通过比较各 RSE值，选择RES值明显大于其临近的RES值所对应的单元为可疑单元。这 些单元的刚度参数α将应用于模型修正。

3）基于静响应灵敏度的模型修正

假设结构系统的损伤模型为：

$K_d=\underset{j=1}{\overset{\mathrm{nel}}{\Sigma }}{\alpha }_j{K}_j^e$

其中，α=1表示没有损伤，α=0表示彻底损坏。则计算得到关于损伤参 数的静响应灵敏度为：

$\frac{\partial U_d}{\partial {\alpha }_j}=-K_d^{-1}{K}_j^e{U}_d$

利用静响应灵敏度组成灵敏度矩阵：

其中k为测量数据的数量，n为可疑单元的数量。

然后求解损伤识别方程：

[S]{δα}={δU}

其中δU=U_m-U_d。一般来说，此方程为病态方程，所以应用Tikhonov正则 化方法和L-curve方法求解。求得每一步迭代（第it步）中的δα，则下一次迭 代（第it+1步）中的损伤参数变为：

α^it+1=α^it+δα^it

然后重新计算待修正模型的总体刚度K_d，灵敏度矩阵S和位移差δU。如 此循环迭代，直到满足预设精度值tolerance为止：

$\frac{\left|\right|{\alpha }^{\mathrm{it}+1}-{\alpha }^{\mathrm{it}}\left|\right|}{\left|\right|{\alpha }^{\mathrm{it}}\left|\right|}\le \mathrm{tolerance}$

具体识别流程如图8所示。

实施例1：对一简支梁进行损伤识别

如图1所示矩形界面简支梁，几何参数如图中所示，结构参数分别为：杨 氏模量E=3.0×10¹⁰N/m²，材料密度ρ=2800kg/m³。将该简支梁分解为如图 1所示的20个欧拉-伯努利梁单元。假设第5和第14个单元发生损伤，其刚度 分别折减了15%和20%。

1）第一步：选取外载方案

图2为预设的两种待选方案，分别计算其选取指标可以得到：D₁=3.9750， D₂=2.9814。通过比较得知D₂较小。说明图2（b）中的外载施加方案较优。

2）第二步：预测可疑单元

根据选定的外载方案，测量节点4，6，8，11，14，16和18的挠度值。应 用测量数据，计算各单元的RSE值。得到结果如图3所示。从图中可以看到， 单元5，6，12，13，14和15的RSE值明显大于与他们临近的单元。所以这6 个单元被选定为可疑的损伤单元。

3）第三步：模型修正

利用上述的基于静响应灵敏度的模型修正方法，选取单元5，6，12，13， 14和15的损伤参数进行识别。经过33次迭代，得到如图4所示的识别结果。 可以从图中看出，损伤位置和损伤程度均能精确地识别出来。

实施例2：对一三跨连续梁进行损伤识别

如图5所示的矩形截面三跨连续梁，其几何参数如图中所示。结构参数为： 杨氏模量E=3.0×10¹⁰N/m²，材料密度ρ=2800kg/m²。把该梁如图5所示离 散为24个单元。假设单元4和21存在损伤，损伤程度均为20%。加载方案已 经预先选定，见图5。

1）预测可以单元

根据力的作用位置，选取测点。测量节点3，5，7，11，13，15，19，21 和23的挠度值。利用这9个测量数据，计算24个单元的各自的RSE值。得到 结果如图6所示。可以看到，单元4，5，20和22的RSE值明显突出。则选取 这4个单元为可疑单元进行进一步识别。

2）模型修正

对单元4，5，20和22的损伤参数进行识别。应用响应灵敏度的方法，经 过23次迭代，得到如图7所示的识别结果。识别结果再次说明，本发明的方法 可以很好地识别损伤位置和损伤程度。

以上所述的本发明的实施方式，并不构成对本发明保护范围的限定。任何 在本发明的精神原则之内所作出的修改、等同替换和改进等，均应包含在本发 明的权利要求保护范围之内。