Movatterモバイル変換

[0]ホーム
URL:

CN114692465A - Nondestructive identification method of bridge damage position, storage medium and equipment - Google Patents

Nondestructive identification method of bridge damage position, storage medium and equipmentDownload PDF

Info

Publication number
CN114692465A
CN114692465ACN202210398866.7ACN202210398866ACN114692465ACN 114692465 ACN114692465 ACN 114692465ACN 202210398866 ACN202210398866 ACN 202210398866ACN 114692465 ACN114692465 ACN 114692465A
Authority
CN
China
Prior art keywords
covariance
time step
state quantity
bridge
unscented kalman
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Granted
Application number
CN202210398866.7A
Other languages
Chinese (zh)
Other versions
CN114692465B (en
Inventor
张延哲
丁勇
张立平
卜建清
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Shijiazhuang Tiedao University
Harbin Institute of Technology Shenzhen
Original Assignee
Shijiazhuang Tiedao University
Harbin Institute of Technology Shenzhen
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Shijiazhuang Tiedao University, Harbin Institute of Technology ShenzhenfiledCriticalShijiazhuang Tiedao University
Priority to CN202210398866.7ApriorityCriticalpatent/CN114692465B/en
Publication of CN114692465ApublicationCriticalpatent/CN114692465A/en
Application grantedgrantedCritical
Publication of CN114692465BpublicationCriticalpatent/CN114692465B/en
Activelegal-statusCriticalCurrent
Anticipated expirationlegal-statusCritical

Links

Images

Classifications

Landscapes

Abstract

A nondestructive identification method, a storage medium and equipment for a bridge damage position belong to the technical field of bridge engineering health monitoring. In order to solve the problem that the existing detection technology cannot simply and conveniently identify the damage position of the bridge structure, the invention determines the damage position based on the unscented Kalman filter algorithm and the self-adaptive unscented Kalman filter algorithm, calculates the sensitive parameter corresponding to each step in the filtering process, and calculates the sensitive parameter according to the sensitive parameter etakSelecting the next time course curveThe filter algorithm of time step, need form a new diagonal square matrix according to all main diagonal elements of the covariance of the state quantity while carrying out the algorithm of the adaptive unscented Kalman filter, enlarge the covariance value that every elastic modulus parameter of the square matrix corresponds to sequentially, and carry out the one-step complete unscented Kalman filtering operation and obtain the sensitive parameter etakAnd finding out the damage position; and expanding the covariance corresponding to the position and continuously executing the (k +1) th time step until all the damage positions are obtained at the end of the whole circulation.

Description

Translated fromChinese
桥梁损伤位置的无损识别方法、存储介质及设备Nondestructive identification method, storage medium and equipment for bridge damage location

技术领域：Technical field:

本发明属于桥梁工程健康监测技术领域，具体涉及桥梁损伤位置的无损识别方法、存储介质及设备。The invention belongs to the technical field of bridge engineering health monitoring, and particularly relates to a nondestructive identification method, storage medium and equipment of bridge damage positions.

背景技术：Background technique:

桥梁作为交通运输生命线工程中的节点和咽喉，在促进区域链接、保障物资运输、推动经济发展等方面有着不可替代的作用。随着桥梁服役年限的增长，重载交通作用、环境侵蚀、材料老化、疲劳效应等对桥梁结构的性能影响越来越大，致使桥梁局部抗力下降、损伤累积，诱发破坏甚至坍塌。如何尽早获取桥梁损伤位置，对于桥梁及时养护维修、预防桥梁大范围破坏、保障人员安全及节能环保方面意义重大。As the nodes and throats of transportation lifeline projects, bridges play an irreplaceable role in promoting regional links, ensuring material transportation, and promoting economic development. As the service life of bridges increases, heavy-duty traffic, environmental erosion, material aging, and fatigue effects have an increasing impact on the performance of bridge structures, resulting in a decrease in the local resistance of bridges, accumulation of damage, and damage or even collapse. How to obtain the location of bridge damage as soon as possible is of great significance for timely maintenance and repair of bridges, preventing large-scale bridge damage, ensuring personnel safety, and energy conservation and environmental protection.

目前，按照是否破坏原结构，桥梁损伤位置检测方法主要分为有损检测法和无损检测法两类，其中有损检测法主要是借助工具获取桥梁某位置的实体样本，如通过打孔取样等，并进一步观察或试验分析，这无疑会对原结构产生影响，且检测速度慢，效率低；无损检测法包括超声波探测技术、声发射探测技术、探地雷达检测技术等，但这都需要专门的设备，检测成本比较高，且需要专业的技术人员操作。因此，亟需发展适用于桥梁结构且易于操作的无损检测方法。At present, according to whether the original structure is damaged, bridge damage location detection methods are mainly divided into two categories: destructive testing method and non-destructive testing method. Among them, the non-destructive testing method mainly uses tools to obtain a physical sample of a certain position of the bridge, such as sampling by punching, etc. , and further observation or experimental analysis, which will undoubtedly have an impact on the original structure, and the detection speed is slow and the efficiency is low; non-destructive testing methods include ultrasonic detection technology, acoustic emission detection technology, ground penetrating radar detection technology, etc., but these all require special The detection cost is relatively high and requires professional technicians to operate. Therefore, there is an urgent need to develop a non-destructive testing method suitable for bridge structures and easy to operate.

发明内容：Invention content:

本发明为了解决现有的检测技术不能简便识别桥梁结构损伤位置的问题，进而提出一种桥梁损伤位置的无损识别方法。In order to solve the problem that the existing detection technology cannot easily identify the damaged position of the bridge structure, the invention further proposes a non-destructive identification method of the damaged position of the bridge.

一种桥梁损伤位置的无损识别方法，包括以下步骤：A method for non-destructive identification of bridge damage positions, comprising the following steps:

针对桥梁结构，确定桥梁结构对应状态的初始值，并组成初始状态量χ0，并依据卡尔曼滤波原理确定初始状态量的协方差矩阵，简称初始状态量协方差P0；其中χ0和P0分别称作第0时间步的状态量和状态量协方差；For the bridge structure, the initial value of the corresponding state of the bridge structure is determined, and the initial state quantity χ0 is formed, and the covariance matrix of the initial state quantity is determined according to the Kalman filter principle, referred to as the initial state quantity covariance P0 ; where χ0 and P0 is called the state quantity and state quantity covariance of the 0th time step, respectively;

基于无迹卡尔曼滤波器算法进行初步识别，基于无迹卡尔曼滤波器算法进行初步识别的过程中，需要基于无迹卡尔曼滤波器算法的量测更新步计算第k时间步的观测误差εk和第k时间步的量测预测协方差Pyy,k，并基于εk和Pyy,k计算输出每步对应的灵敏参数

Figure BDA0003598743550000011
Preliminary identification based on the unscented Kalman filter algorithm, in the process of preliminary identification based on the unscented Kalman filter algorithm, it is necessary to calculate the observation error ε of the kth time step based on the measurement update step of the unscented Kalman filter algorithm The measured and predicted covariance Pyy,k ofk and the kth time step, and the sensitivity parameters corresponding to each step are calculated and output based on εk and Pyy,k
Figure BDA0003598743550000011

然后绘制输出灵敏参数ηk的时程曲线，如果ηk时程曲线有峰值脉冲出现，则基于自适应无迹卡尔曼滤波器算法识别损伤位置，在基于自适应无迹卡尔曼滤波器算法进行识别的过程中，需要判断计算的ηk与灵敏参数阈值η0的大小，如果ηk＜η0，则基于无迹卡尔曼滤波器算法继续识别；如果ηk≥η0，则继续执行以下步骤：Then draw the time-history curve of the output sensitivity parameterηk . If the peak pulse appears in theηk time-history curve, the damage location is identified based on the adaptive unscented Kalman filter algorithm. During the identification process, it is necessary to judge the size of the calculated ηk and the sensitive parameter threshold η0 . If ηk <η0 , continue to identify based on the unscented Kalman filter algorithm; if ηk ≥ η0 , continue to perform the following step:

设置一个初始值为0且维度等于初始状态量χ0的向量，用字母L表示；取状态量协方差

Figure BDA0003598743550000021
的全部主对角元素组成一个新的对角方阵
Figure BDA0003598743550000022
且保持原对角元素的位置不变；依次扩大
Figure BDA0003598743550000023
中每一个弹性模量参数对应的协方差值，且每次只对一个协方差值进行扩大，并执行一步完整的无迹卡尔曼滤波运算，针对协方差值对应的行号或列号为z的扩大，执行一步完整的无迹卡尔曼滤波运算得到灵敏参数ηk，并令L(z)＝ηk；找到最小的L对应的位置号zmin，该位置就是损伤位置；然后令
Figure BDA0003598743550000024
扩大20倍，并继续执行第(k+1)时间步，并再次判断ηk+1与灵敏参数阈值η0的大小，直到全部循环结束得出所有的损伤位置。Set a vector with an initial value of 0 and a dimension equal to the initial state quantity χ0 , represented by the letter L; take the state quantity covariance
Figure BDA0003598743550000021
All the main diagonal elements of form a new diagonal square matrix
Figure BDA0003598743550000022
And keep the position of the original diagonal element unchanged; expand in turn
Figure BDA0003598743550000023
The covariance value corresponding to each elastic modulus parameter in the , and only one covariance value is expanded at a time, and a complete unscented Kalman filtering operation is performed, aiming at the row number or column corresponding to the covariance value. For the expansion of the number z, perform a complete unscented Kalman filter operation to obtain the sensitive parameter ηk , and let L(z)=ηk ; find the position number zmin corresponding to the smallest L, which is the damage position; then make
Figure BDA0003598743550000024
Expand 20 times, and continue to perform the (k+1)th time step, and judge the size of ηk+1 and the sensitive parameter threshold η0 again, until the end of all cycles to obtain all damage locations.

进一步地，利用无迹卡尔曼滤波器算法进行识别的过程中和利用自适应无迹卡尔曼滤波器算法进行识别的过程中，选择各梁单元节点的竖向位移响应作为观测值y。Further, during the identification process using the unscented Kalman filter algorithm and the identification process using the adaptive unscented Kalman filter algorithm, the vertical displacement response of each beam element node is selected as the observation value y.

进一步地，所述的桥梁结构对应状态包括桥梁各梁单元的弹性模量、桥梁结构的位移和速度。Further, the corresponding state of the bridge structure includes the elastic modulus of each beam element of the bridge, the displacement and the speed of the bridge structure.

进一步地，在基于自适应无迹卡尔曼滤波器算法进行识别的过程中，计算灵敏参数ηk的过程包括以下步骤：Further, in the process of identifying based on the adaptive unscented Kalman filter algorithm, the process of calculating the sensitive parameter ηk includes the following steps:

步骤7.1、基于无迹卡尔曼滤波器算法的UT变换原理，使用第(k-1)时间步的状态量χk-1和状态量协方差Pk-1生成(2n+1)个sigma点，并通过状态方程求解出每个sigma点对应的状态量

Figure BDA0003598743550000025
其中k从1开始，且k∈[1,N]，N为总的时间步数，n为状态量的维数，i为第i个sigma点，且i∈[1,2n+1]；Step 7.1. Based on the UT transformation principle of the unscented Kalman filter algorithm, use the state quantity χk-1 and the state quantity covariance Pk-1 of the (k-1)th time step to generate (2n+1) sigma points , and solve the state quantity corresponding to each sigma point through the state equation
Figure BDA0003598743550000025
Where k starts from 1, and k∈[1,N], N is the total number of time steps, n is the dimension of the state quantity, i is the i-th sigma point, and i∈[1,2n+1];

步骤7.2、基于无迹卡尔曼滤波器算法的时间更新步完成从第(k-1)时间步到第k时间步的状态量和状态量协方差的更新，分别记作

Figure BDA0003598743550000026
和
Figure BDA0003598743550000027
公式如下所述：Step 7.2. The time update step based on the unscented Kalman filter algorithm completes the update of the state quantity and the state quantity covariance from the (k-1)th time step to the kth time step, which are respectively recorded as
Figure BDA0003598743550000026
and
Figure BDA0003598743550000027
The formula is as follows:

Figure BDA0003598743550000028
Figure BDA0003598743550000028

Figure BDA0003598743550000029
Figure BDA0003598743550000029

式中，

Figure BDA00035987435500000210
和
Figure BDA00035987435500000211
分别为第k时间步第i个sigma点的权重值，
Figure BDA00035987435500000212
为第k时间步第i个sigma点对应的状态量估计值，Qk为第k时间步的噪声；In the formula,
Figure BDA00035987435500000210
and
Figure BDA00035987435500000211
are the weight values of the i-th sigma point at the k-th time step, respectively,
Figure BDA00035987435500000212
is the estimated state quantity corresponding to the i-th sigma point at the k-th time step, and Qk is the noise at the k-th time step;

步骤7.3、基于无迹卡尔曼滤波器算法的UT变换原理，使用步骤7.2中更新的

Figure BDA00035987435500000213
和
Figure BDA00035987435500000214
生成(2n+1)个sigma点，并通过观测方程求解出每个sigma点对应的观测估计值
Figure BDA0003598743550000031
Step 7.3, based on the UT transformation principle of the unscented Kalman filter algorithm, use the updated in step 7.2
Figure BDA00035987435500000213
and
Figure BDA00035987435500000214
Generate (2n+1) sigma points, and solve the observation estimate corresponding to each sigma point through the observation equation
Figure BDA0003598743550000031

步骤7.4、基于无迹卡尔曼滤波器算法的量测更新步计算输出第k时间步的量测预测值

Figure BDA0003598743550000032
且
Figure BDA0003598743550000033
Step 7.4. Calculate and output the measurement prediction value of the kth time step based on the measurement update step of the unscented Kalman filter algorithm
Figure BDA0003598743550000032
and
Figure BDA0003598743550000033

式中，

Figure BDA0003598743550000034
为第k时间步第i个sigma点的权重值，
Figure BDA0003598743550000035
为第k时间步第i个sigma点对应的观测估计值；In the formula,
Figure BDA0003598743550000034
is the weight value of the i-th sigma point at the k-th time step,
Figure BDA0003598743550000035
is the observed estimated value corresponding to the i-th sigma point at the k-th time step;

步骤7.5、基于无迹卡尔曼滤波器算法的量测更新步计算第k时间步的观测误差εk，且

Figure BDA0003598743550000036
Step 7.5. Calculate the observation error εk of the kth time step based on the measurement update step of the unscented Kalman filter algorithm, and
Figure BDA0003598743550000036

式中，yk为第k时间步的观测值，

Figure BDA0003598743550000037
为第k时间步的量测预测值；where yk is the observed value at the kth time step,
Figure BDA0003598743550000037
is the measurement prediction value of the kth time step;

步骤7.6、基于无迹卡尔曼滤波器算法的量测更新步计算第k时间步的量测预测协方差Pyy,k，且

Figure BDA0003598743550000038
Step 7.6. Calculate the measurement prediction covariance Pyy,k at the kth time step based on the measurement update step of the unscented Kalman filter algorithm, and
Figure BDA0003598743550000038

式中，

Figure BDA0003598743550000039
为第k时间步第i个sigma点的权重值，
Figure BDA00035987435500000310
为第k时间步第i个sigma点对应的观测估计值，
Figure BDA00035987435500000311
为第k时间步的量测预测值，Rk为第k时间步的噪声；In the formula,
Figure BDA0003598743550000039
is the weight value of the i-th sigma point at the k-th time step,
Figure BDA00035987435500000310
is the estimated observation value corresponding to the i-th sigma point at the k-th time step,
Figure BDA00035987435500000311
is the measurement prediction value at the kth time step, and Rk is the noise at the kth time step;

步骤7.7、基于步骤7.5及步骤7.6计算的εk和Pyy,k构造灵敏参数ηk，且

Figure BDA00035987435500000312
并计算输出每步的ηk值。Step 7.7. Construct the sensitive parameter ηk based on εk and Pyy,k calculated in step 7.5 and step 7.6, and
Figure BDA00035987435500000312
And calculate the value of ηk for each step of the output.

进一步地，在基于自适应无迹卡尔曼滤波器算法进行识别的过程中，在得到计算灵敏参数ηk后还要计算第k时间步

Figure BDA00035987435500000313
与
Figure BDA00035987435500000314
的互协方差Pxy,k并进行数据更新，具体过程包括以下步骤：Further, in the process of identification based on the adaptive unscented Kalman filter algorithm, after the calculation sensitivity parameter ηk is obtained, the kth time step needs to be calculated.
Figure BDA00035987435500000313
and
Figure BDA00035987435500000314
The cross-covariance Pxy,k of , and the data is updated. The specific process includes the following steps:

步骤7.8、基于无迹卡尔曼滤波器算法的量测更新步计算第k时间步

Figure BDA00035987435500000315
与
Figure BDA00035987435500000316
的互协方差Pxy,k，
Figure BDA00035987435500000317
Step 7.8. Calculate the kth time step based on the measurement update step of the unscented Kalman filter algorithm
Figure BDA00035987435500000315
and
Figure BDA00035987435500000316
The cross-covariance Pxy,k of ,
Figure BDA00035987435500000317

步骤7.9、更新第k时间步的卡尔曼增益矩阵：

Figure BDA00035987435500000318
Step 7.9. Update the Kalman gain matrix at the kth time step:
Figure BDA00035987435500000318

步骤7.10、更新并输出第k时间步的状态量：

Figure BDA00035987435500000319
Step 7.10, update and output the state quantity of the kth time step:
Figure BDA00035987435500000319

步骤7.11、更新并输出第k时间步的状态量协方差：

Figure BDA00035987435500000320
Step 7.11, update and output the state quantity covariance of the kth time step:
Figure BDA00035987435500000320

进一步地，找到最小的L对应的位置号zmin的具体过程包括以下步骤：Further, the specific process of finding the position number zmin corresponding to the smallest L includes the following steps:

步骤7.14、将第k时间步状态量协方差

Figure BDA00035987435500000321
的主对角线元素中第一个弹性模量参数对应的协方差值的行号或列号记为m，将状态量协方差
Figure BDA00035987435500000322
的主对角线元素中最后一个弹性模量参数对应的协方差值的行号或列号记为l，则总的弹性模量参数个数为(l-m+1)；Step 7.14. Calculate the state quantity covariance of the kth time step
Figure BDA00035987435500000321
The row number or column number of the covariance value corresponding to the first elastic modulus parameter in the main diagonal element of
Figure BDA00035987435500000322
The row number or column number of the covariance value corresponding to the last elastic modulus parameter in the main diagonal element of , is denoted as l, and the total number of elastic modulus parameters is (l-m+1);

取当前时间步中

Figure BDA0003598743550000041
的所有主对角元素组成一个新的对角方阵
Figure BDA0003598743550000042
且保持原对角元素的位置不变；依次扩大
Figure BDA0003598743550000043
中每一个弹性模量参数对应的协方差值，令
Figure BDA0003598743550000044
且每次只对一个协方差值进行扩大，
Figure BDA0003598743550000045
的其余元素值保持不变，其中
Figure BDA0003598743550000046
和
Figure BDA0003598743550000047
分别代表行数列数均为z位置处的协方差值，为一个标量；λ为扩大倍数；Take the current time step
Figure BDA0003598743550000041
All the main diagonal elements of form a new diagonal square matrix
Figure BDA0003598743550000042
And keep the position of the original diagonal element unchanged; expand in turn
Figure BDA0003598743550000043
The covariance value corresponding to each elastic modulus parameter in , let
Figure BDA0003598743550000044
And only one covariance value is expanded at a time,
Figure BDA0003598743550000045
The remaining element values of , remain unchanged, where
Figure BDA0003598743550000046
and
Figure BDA0003598743550000047
Respectively represent the covariance value at the z position with the number of rows and columns, which is a scalar; λ is the expansion multiple;

z起始值为m，且z∈[m,l-m+1]；依据z＝m时计算的

Figure BDA0003598743550000048
和第k时间步的状态量
Figure BDA0003598743550000049
执行一步完整的无迹卡尔曼滤波运算，即执行步骤7.1～步骤7.11，并输出ηk，并令L(z)＝ηk；The initial value of z is m, and z∈[m,l-m+1]; calculated according to z=m
Figure BDA0003598743550000048
and the state quantities at the kth time step
Figure BDA0003598743550000049
Perform a complete unscented Kalman filter operation, that is, perform steps 7.1 to 7.11, and output ηk , and let L(z)=ηk ;

步骤7.15、令z＝m+1，继续执行步骤7.14，直到z＝l-m+1时结束；Step 7.15, let z=m+1, continue to perform step 7.14, until z=l-m+1;

步骤7.16、不计零值，找到L中最小值对应的位置号，记为zmin，即zmin为损伤位置。Step 7.16, ignoring the zero value, find the position number corresponding to the minimum value in L, and record it as zmin , that is, zmin is the damage position.

进一步地，所述扩大倍数λ＝1×10ω，10ω等于初始状态量协方差P0中最小数量级的倒数。Further, the enlargement factor λ=1×10ω , and 10ω is equal to the reciprocal of the smallest order of magnitude in the initial state quantity covariance P0 .

进一步地，所述桥梁结构对应的运动控制微分方程或有限元模型是基于欧拉-伯努利梁单元建立的。Further, the motion control differential equation or finite element model corresponding to the bridge structure is established based on Euler-Bernoulli beam elements.

一种存储介质，所述存储介质中存储有至少一条指令，所述至少一条指令由处理器加载并执行以实现所述的一种桥梁损伤位置的无损识别方法。A storage medium, wherein at least one instruction is stored in the storage medium, the at least one instruction is loaded and executed by a processor to implement the method for non-destructive identification of bridge damage locations.

一种桥梁损伤位置的无损识别设备，所述设备包括处理器和存储器，所述存储器中存储有至少一条指令，所述至少一条指令由处理器加载并执行以实现所述的一种桥梁损伤位置的无损识别方法。A non-destructive identification device for a bridge damage location, the device includes a processor and a memory, the memory stores at least one instruction, and the at least one instruction is loaded and executed by the processor to realize the bridge damage location. non-destructive identification method.

有益效果：Beneficial effects:

本发明仅基于桥梁局部竖向位移响应作为观测值即可识别出桥梁的损伤位置，操作简单，且不破坏原有结构。同时，利用本发明对桥梁损伤位置识别，有益于桥梁结构的及时养护维修，可以有效防止大面积破坏，符合绿色、环保、节能、低碳的发展理念。The invention can identify the damaged position of the bridge only based on the local vertical displacement response of the bridge as the observation value, the operation is simple, and the original structure is not damaged. At the same time, using the invention to identify the damaged position of the bridge is beneficial to the timely maintenance and repair of the bridge structure, can effectively prevent large-scale damage, and conforms to the development concept of green, environmental protection, energy saving and low carbon.

附图说明：Description of drawings:

为了易于说明，本发明由下述的具体实施及附图作以详细描述。For ease of description, the present invention is described in detail by the following specific implementations and accompanying drawings.

图1为灵敏参数时程曲线图。Figure 1 is a graph of the sensitive parameter time history.

图2为灵敏参数时程曲线有峰值脉冲的效果图。Fig. 2 is the effect diagram of the peak pulse in the time-history curve of the sensitive parameter.

图3为图2的局部视图。FIG. 3 is a partial view of FIG. 2 .

图4为实施例的车-桥耦合系统图，其中车以一定速度从桥的一端行驶到另一端，图中，1-四分之一车辆模型的车体质量，用字母m1表示；2-四分之一车辆模型的车体与轮胎之间的悬挂刚度，用字母k1表示；3-四分之一车辆模型的车体与轮胎之间的悬挂阻尼，用字母c1表示；4-四分之一车辆模型的轮胎质量，用字母m2表示；5-四分之一车辆模型的轮胎与桥梁之间的接触刚度，用字母k2表示；6-轮胎与桥梁接触点；7-简支边界约束的固定端；8-梁单元；9-梁单元节点；10-简支边界约束的滑动端。Fig. 4 is the vehicle-axle coupling system diagram of the embodiment, wherein the vehicle travels from one end of the bridge to the other end at a certain speed, in the figure, 1-a quarter of the vehicle body mass of the vehicle model is represented by the letter m1 ; 2 - Suspension stiffness between the body and tires of the quarter vehicle model, denoted by the letter k1 ; 3 - Suspension damping between the body and tires of the quarter vehicle model, denoted by the letter c1 ; 4 - tire mass of the quarter vehicle model, denoted by the letter m2 ; 5 - contact stiffness between the tire and the bridge of the quarter vehicle model, denoted by the letter k2 ; 6 - tire-bridge contact point; 7 - Fixed end of simply supported boundary constraints; 8 - Beam elements; 9 - Beam element nodes; 10 - Sliding ends of simply supported boundary constraints.

具体实施方式：Detailed ways:

为使本发明的目的、技术方案和优点更加清楚明了，下面通过附图中示出的具体实施例来描述本发明。但是应该理解，这些描述只是示例性的，而并非要限制本发明的范围。此外，在以下说明中，省略了对公知结构和技术的描述，以避免不必要地混淆本发明的概念。In order to make the objectives, technical solutions and advantages of the present invention clearer, the present invention will be described below through specific embodiments shown in the accompanying drawings. It should be understood, however, that these descriptions are exemplary only, and are not intended to limit the scope of the present invention. Also, in the following description, descriptions of well-known structures and techniques are omitted to avoid unnecessarily obscuring the concepts of the present invention.

在此，还需要说明的是，为了避免因不必要的细节而模糊了本发明，在附图中仅仅示出了与根据本发明的方案密切相关的结构和/或处理步骤，而省略了与本发明关系不大的其它细节。Here, it should also be noted that, in order to avoid obscuring the present invention due to unnecessary details, only the structures and/or processing steps closely related to the solution according to the present invention are shown in the drawings, and the related structures and/or processing steps are omitted. Other details not relevant to the invention.

损伤位置识别时需要先制定损伤指标，并基于此判断桥梁结构是否发生损伤。一般认为桥梁结构的质量恒定，阻尼为质量和刚度的函数，且为了方便分析，一般都假定为瑞利阻尼，由此，桥梁损伤指标一般通过刚度参数反映。进一步，刚度由弹性模量和截面惯性矩组成，且一般忽略桥梁横截面的变化。因此，更深层次的损伤指标由弹性模量这一微观参数表示。此外，考虑到损伤位置识别的时效性及桥梁的尺寸特征，一般基于欧拉-伯努利梁单元构建桥梁的有限元模型，或者基于欧拉-伯努利梁单元建立桥梁结构的运动控制微分方程；并将此作为真实桥梁结构的合理简化。桥梁结构的位移、速度等响应可通过梁单元节点自由度的对应值表征，且每个梁单元的弹性模量属性都可独立赋予。由此，可根据不同梁单元的属性判断桥梁具体的损伤位置。When identifying the damage location, it is necessary to first formulate the damage index, and based on this, it is necessary to judge whether the bridge structure is damaged. It is generally believed that the mass of the bridge structure is constant, and the damping is a function of mass and stiffness. For the convenience of analysis, Rayleigh damping is generally assumed. Therefore, the bridge damage index is generally reflected by the stiffness parameter. Further, stiffness consists of elastic modulus and cross-section moment of inertia, and changes in bridge cross-section are generally ignored. Therefore, a deeper damage index is represented by the microscopic parameter of elastic modulus. In addition, considering the timeliness of damage location identification and the size characteristics of bridges, the finite element model of bridges is generally constructed based on Euler-Bernoulli beam elements, or the motion control differential of bridge structures is established based on Euler-Bernoulli beam elements. equation; and use this as a reasonable simplification of the real bridge structure. The response of the bridge structure such as displacement and velocity can be characterized by the corresponding values of the degrees of freedom of the beam element nodes, and the elastic modulus properties of each beam element can be assigned independently. Therefore, the specific damage location of the bridge can be determined according to the properties of different beam elements.

具体实施方式一：Specific implementation one:

本实施方式为一种桥梁损伤位置的无损识别方法，由于本发明是基于桥梁局部竖向位移响应进行识别的，所以本发明是一种基于桥梁局部竖向位移响应识别桥梁损伤位置的无损识别方法，本质上是一种基于自适应无迹卡尔曼滤波器算法的识别方法。This embodiment is a non-destructive identification method for bridge damage location. Since the present invention is based on the local vertical displacement response of the bridge, the present invention is a non-destructive identification method for identifying the bridge damage location based on the local vertical displacement response of the bridge. , which is essentially a recognition method based on the adaptive unscented Kalman filter algorithm.

本实施方式所述的一种桥梁损伤位置的无损识别方法，包括以下步骤：A method for nondestructive identification of bridge damage locations described in this embodiment includes the following steps:

步骤1、根据达朗贝尔原理、结构动力学及有限元理论，推导结构的运动控制微分方程或建立结构有限元模型，所述结构为桥梁结构，为方便表述，简称结构；Step 1. According to D'Alembert's principle, structural dynamics and finite element theory, the motion control differential equation of the structure is deduced or a finite element model of the structure is established. The structure is a bridge structure, which is abbreviated as structure for convenience of expression;

在此过程中，基于欧拉-伯努利梁单元建立桥梁结构的运动控制微分方程或有限元模型，并将各个梁单元的弹性模量作为损伤指标；In this process, the motion control differential equation or finite element model of the bridge structure is established based on Euler-Bernoulli beam elements, and the elastic modulus of each beam element is used as the damage index;

步骤2、将结构位移X、结构速度

Figure BDA0003598743550000062
以及弹性模量E组成的向量称作状态量，用符号
Figure BDA0003598743550000061
表示；Step 2. The structure displacement X, the structure speed
Figure BDA0003598743550000062
And the vector composed of the elastic modulus E is called the state quantity, with the symbol
Figure BDA0003598743550000061
express;

在一些实施例中，根据状态量χ，依据步骤1中结构的运动控制微分方程，并基于线性代数矩阵运算及数值分析的数值微分和积分运算，或，基于线性代数矩阵运算及结构动力学的Newmark-β法等数值积分方法推导出求解状态量的方程关系；In some embodiments, according to the state quantity x, according to the motion control differential equation of the structure instep 1, and based on linear algebraic matrix operations and numerical differential and integral operations of numerical analysis, or, based on linear algebraic matrix operations and structural dynamics Numerical integration methods such as Newmark-β method are used to derive the equation relationship for solving state quantities;

或者，or,

在另一些实施例中，根据状态量χ，依据步骤1中结构有限元模型进行状态量的输出设置，并将其作为状态方程；该过程中考虑噪声影响；In other embodiments, according to the state quantity χ, the output setting of the state quantity is performed according to the structural finite element model instep 1, and it is used as the state equation; the influence of noise is considered in this process;

步骤3、将传感器采集的各梁单元节点的竖向位移响应作为观测值y。Step 3. Take the vertical displacement response of each beam element node collected by the sensor as the observation value y.

观测值一般是指包含位移、速度、加速度、应变、应力、力等一系列能通过传感器测量的物理量，对于桥梁结构，本实施方式中使用位移测量值作为观测值。The observed value generally refers to a series of physical quantities that can be measured by sensors, including displacement, velocity, acceleration, strain, stress, force, etc. For the bridge structure, the displacement measured value is used as the observed value in this embodiment.

本发明中，状态方程是状态量的微分运算，状态量写成位移和速度的形式，是方便状态方程的推导，因为对位移求导是速度，对速度求导是加速度，而根据结构的运动控制微分方程很容易就能推出状态方程。另外，当观测值为位移时，因为状态量中的位移和观测值中的位移含义相同，那么观测方程中的位移关系也很容易得到。概括说，状态量中的位移和速度是为状态方程推导以及观测方程中的位移关系服务的，观测值中的位移是为修正观测方程计算的位移服务的。由于此过程的具体内容为公知常识，所以本发明中不再赘述。需要说明的是，本发明中的速度相当于中间量，不需要输出，它的存在只是为状态方程的推导服务的，它的计算由算法递推迭代完成。In the present invention, the state equation is the differential operation of the state quantity, and the state quantity is written in the form of displacement and velocity, which is convenient for the derivation of the state equation, because the derivation of the displacement is the velocity, and the derivation of the velocity is the acceleration, and according to the motion control of the structure Differential equations can easily deduce the equation of state. In addition, when the observed value is displacement, because the displacement in the state quantity and the displacement in the observed value have the same meaning, the displacement relationship in the observation equation is also easy to obtain. In general, the displacement and velocity in the state quantity serve for the derivation of the state equation and the displacement relationship in the observation equation, and the displacement in the observation value serves for the displacement calculated by the correction of the observation equation. Since the specific content of this process is common knowledge, it will not be repeated in the present invention. It should be noted that the speed in the present invention is equivalent to an intermediate quantity and does not require an output. Its existence is only for the derivation of the state equation, and its calculation is completed by the algorithm recursive iteration.

根据观测值类型，依据步骤1中结构的运动控制微分方程并基于线性代数矩阵运算、数学移项、数学合并同类项等知识推导出求解观测值的方程关系；According to the type of observation value, according to the motion control differential equation of the structure instep 1, and based on the knowledge of linear algebra matrix operation, mathematical shifting, mathematical merging of similar terms, etc., the equation relationship for solving the observation value is deduced;

或者，or,

根据观测值类型，依据步骤1中结构有限元模型进行相应观测值的输出设置，并将其作为观测方程；该过程中考虑噪声影响；According to the type of observation value, according to the structural finite element model instep 1, the output setting of the corresponding observation value is carried out, and it is used as the observation equation; the influence of noise is considered in this process;

步骤4、将结构的初始位移、结构的初始速度以及各弹性模量的初始值组成的向量称作初始状态量，用符号χ0表示，并根据卡尔曼滤波器算法原理得出初始状态量的协方差矩阵，简称初始状态量协方差，用符号P0表示，其中χ0和P0分别称作第0时间步(启始步)的状态量和状态量协方差；Step 4. The vector composed of the initial displacement of the structure, the initial velocity of the structure and the initial value of each elastic modulus is called the initial state quantity, which is represented by the symbol χ0 , and the initial state quantity is obtained according to the principle of the Kalman filter algorithm. The covariance matrix, referred to as the initial state quantity covariance, is represented by the symbol P0 , wherein χ0 and P0 are respectively called the state quantity and state quantity covariance of the 0th time step (starting step);

步骤5、基于无迹卡尔曼滤波器算法进行初步识别，过程如下：Step 5. Perform preliminary identification based on the unscented Kalman filter algorithm. The process is as follows:

步骤5.1、基于无迹卡尔曼滤波器算法的UT变换原理，使用第(k-1)时间步的状态量χk-1和状态量协方差Pk-1生成(2n+1)个sigma点，并通过状态方程求解出每个sigma点对应的状态量

Figure BDA0003598743550000071
其中k从1开始，且k∈[1,N]，N为总的时间步数，n为状态量的维数，i为第i个sigma点，且i∈[1,2n+1]；Step 5.1. Based on the UT transformation principle of the unscented Kalman filter algorithm, use the state quantity χk-1 and the state quantity covariance Pk-1 of the (k-1)th time step to generate (2n+1) sigma points , and solve the state quantity corresponding to each sigma point through the state equation
Figure BDA0003598743550000071
Where k starts from 1, and k∈[1,N], N is the total number of time steps, n is the dimension of the state quantity, i is the i-th sigma point, and i∈[1,2n+1];

步骤5.2、基于无迹卡尔曼滤波器算法的时间更新步完成从第(k-1)时间步到第k时间步的状态量和状态量协方差的更新，分别记作

Figure BDA0003598743550000072
和
Figure BDA0003598743550000073
公式如下所述：Step 5.2. The time update step based on the unscented Kalman filter algorithm completes the update of the state quantity and the state quantity covariance from the (k-1)th time step to the kth time step, which are respectively recorded as
Figure BDA0003598743550000072
and
Figure BDA0003598743550000073
The formula is as follows:

Figure BDA0003598743550000074
Figure BDA0003598743550000074

Figure BDA0003598743550000075
Figure BDA0003598743550000075

式中，

Figure BDA0003598743550000076
和
Figure BDA0003598743550000077
分别为第k时间步第i个sigma点的权重值，
Figure BDA0003598743550000078
为第k时间步第i个sigma点对应的状态量估计值，Qk为第k时间步的噪声；In the formula,
Figure BDA0003598743550000076
and
Figure BDA0003598743550000077
are the weight values of the i-th sigma point at the k-th time step, respectively,
Figure BDA0003598743550000078
is the estimated state quantity corresponding to the i-th sigma point at the k-th time step, and Qk is the noise at the k-th time step;

步骤5.3、基于无迹卡尔曼滤波器算法的UT变换原理，使用步骤5.2中更新的

Figure BDA0003598743550000079
和
Figure BDA00035987435500000710
生成(2n+1)个sigma点，并通过观测方程求解出每个sigma点对应的观测估计值
Figure BDA00035987435500000711
步骤5.4、基于无迹卡尔曼滤波器算法的量测更新步计算输出第k时间步的量测预测值
Figure BDA00035987435500000712
且
Figure BDA00035987435500000713
Step 5.3, based on the UT transformation principle of the unscented Kalman filter algorithm, use the updated in step 5.2
Figure BDA0003598743550000079
and
Figure BDA00035987435500000710
Generate (2n+1) sigma points, and solve the observation estimate corresponding to each sigma point through the observation equation
Figure BDA00035987435500000711
Step 5.4. Calculate and output the measurement prediction value of the kth time step based on the measurement update step of the unscented Kalman filter algorithm
Figure BDA00035987435500000712
and
Figure BDA00035987435500000713

式中，

Figure BDA00035987435500000714
为第k时间步第i个sigma点的权重值，
Figure BDA00035987435500000715
为第k时间步第i个sigma点对应的观测估计值；In the formula,
Figure BDA00035987435500000714
is the weight value of the i-th sigma point at the k-th time step,
Figure BDA00035987435500000715
is the observed estimated value corresponding to the i-th sigma point at the k-th time step;

步骤5.5、基于无迹卡尔曼滤波器算法的量测更新步计算第k时间步的观测误差εk，且

Figure BDA00035987435500000716
Step 5.5. Calculate the observation error εk of the kth time step based on the measurement update step of the unscented Kalman filter algorithm, and
Figure BDA00035987435500000716

式中，yk为第k时间步的观测值，

Figure BDA00035987435500000717
为第k时间步的量测预测值；where yk is the observed value at the kth time step,
Figure BDA00035987435500000717
is the measurement prediction value of the kth time step;

步骤5.6、基于无迹卡尔曼滤波器算法的量测更新步计算第k时间步的量测预测协方差Pyy,k，且

Figure BDA00035987435500000718
Step 5.6. Calculate the measurement prediction covariance Pyy,k at the kth time step based on the measurement update step of the unscented Kalman filter algorithm, and
Figure BDA00035987435500000718

式中，

Figure BDA00035987435500000719
为第k时间步第i个sigma点的权重值，
Figure BDA00035987435500000720
为第k时间步第i个sigma点对应的观测估计值，
Figure BDA00035987435500000721
为第k时间步的量测预测值，Rk为第k时间步的噪声；In the formula,
Figure BDA00035987435500000719
is the weight value of the i-th sigma point at the k-th time step,
Figure BDA00035987435500000720
is the estimated observation value corresponding to the i-th sigma point at the k-th time step,
Figure BDA00035987435500000721
is the measurement prediction value at the kth time step, and Rk is the noise at the kth time step;

步骤5.7、基于步骤5.5及步骤5.6计算的εk和Pyy,k构造灵敏参数ηk，且

Figure BDA00035987435500000722
并计算输出每步的ηk值；Step 5.7, construct the sensitive parameter ηk based on the εk and Pyy,k calculated in step 5.5 and step 5.6, and
Figure BDA00035987435500000722
And calculate and output the ηk value of each step;

步骤5.8、基于无迹卡尔曼滤波器算法的量测更新步计算第k时间步

Figure BDA0003598743550000081
与
Figure BDA0003598743550000082
的互协方差Pxy,k，
Figure BDA0003598743550000083
Step 5.8. Calculate the kth time step based on the measurement update step of the unscented Kalman filter algorithm
Figure BDA0003598743550000081
and
Figure BDA0003598743550000082
The cross-covariance Pxy,k of ,
Figure BDA0003598743550000083

步骤5.9、更新第k时间步的卡尔曼增益矩阵：

Figure BDA0003598743550000084
Step 5.9, update the Kalman gain matrix at the kth time step:
Figure BDA0003598743550000084

步骤5.10、更新并输出第k时间步的状态量：

Figure BDA0003598743550000085
Step 5.10, update and output the state quantity of the kth time step:
Figure BDA0003598743550000085

步骤5.11、更新并输出第k时间步的状态量协方差：

Figure BDA0003598743550000086
Step 5.11, update and output the state quantity covariance of the kth time step:
Figure BDA0003598743550000086

步骤5.12、时间步变为(k+1)，重复步骤5.1～步骤5.11直到最大时间步N完成，即直到循环结束。Step 5.12, the time step becomes (k+1), and steps 5.1 to 5.11 are repeated until the maximum time step N is completed, that is, until the cycle ends.

步骤6、绘制步骤5.7中输出的灵敏参数时程曲线，如果曲线整体平稳，无脉冲响应出现(参见图1)，则无需调用自适应无迹卡尔曼滤波器算法，按常规无迹卡尔曼滤波器算法，即步骤5(但忽略步骤5.5和步骤5.7)识别即可；如果ηk时程曲线有峰值脉冲出现(参见图2)，则需调用自适应无迹卡尔曼滤波器算法(步骤7)识别损伤位置，并令峰值脉冲前出现的最大灵敏参数值等于灵敏参数阈值η0(参见图3，其中图3为图2的局部视图)。Step 6. Draw the sensitive parameter time-history curve output in step 5.7. If the overall curve is stable and there is no impulse response (see Figure 1), there is no need to call the adaptive unscented Kalman filter algorithm, and the conventional unscented Kalman filter is used. If there is a peak pulse in the ηk time-history curve (see Figure 2), the adaptive unscented Kalman filter algorithm needs to be called (step 7). ) to identify the lesion location and make the maximum sensitive parameter value that occurs before the peak pulse equal to the sensitive parameter threshold η0 (see FIG. 3 , which is a partial view of FIG. 2 ).

步骤7、基于自适应无迹卡尔曼滤波器算法进行识别，过程如下：Step 7. Identify based on the adaptive unscented Kalman filter algorithm. The process is as follows:

步骤7.1、同步骤5.1；Step 7.1, same as step 5.1;

步骤7.2、同步骤5.2；Step 7.2, same as step 5.2;

步骤7.3、同步骤5.3；Step 7.3, same as step 5.3;

步骤7.4、同步骤5.4；Step 7.4, same as step 5.4;

步骤7.5、同步骤5.5；Step 7.5, same as step 5.5;

步骤7.6、同步骤5.6；Step 7.6, same as step 5.6;

步骤7.7、同步骤5.7；Step 7.7, same as step 5.7;

步骤7.8、同步骤5.8；Step 7.8, same as step 5.8;

步骤7.9、同步骤5.9；Step 7.9, same as step 5.9;

步骤7.10、同步骤5.10；Step 7.10, same as step 5.10;

步骤7.11、同步骤5.11；Step 7.11, same as step 5.11;

步骤7.12、判断步骤7.7中计算的ηk与灵敏参数阈值η0的大小，如果ηk＜η0，则时间步变为(k+1)，重复步骤7.1～步骤7.11继续计算；如果ηk≥η0，则继续执行步骤7.13～步骤7.18；Step 7.12. Determine the size of ηk calculated in step 7.7 and the sensitive parameter threshold η0. If ηk <η0 , the time step becomes (k+1), and repeat steps 7.1 to 7.11 to continue the calculation; if ηk ≥η0 , then continue to perform steps 7.13 to 7.18;

步骤7.13、设置一个初始值为0且维度为n的向量，用字母L表示，需要说明的是代表该向量的字符可自选，此处使用字母L表示。Step 7.13. Set a vector with an initial value of 0 and a dimension of n, which is represented by the letter L. It should be noted that the character representing the vector can be selected by yourself, which is represented by the letter L here.

步骤7.14、将状态量协方差

Figure BDA0003598743550000091
的主对角线元素中第一个弹性模量参数对应的协方差值的行号或列号记为m，将状态量协方差
Figure BDA0003598743550000092
的主对角线元素中最后一个弹性模量参数对应的协方差值的行号或列号记为l，则总的弹性模量参数个数为(l-m+1)。Step 7.14, the state quantity covariance
Figure BDA0003598743550000091
The row number or column number of the covariance value corresponding to the first elastic modulus parameter in the main diagonal element of
Figure BDA0003598743550000092
The row number or column number of the covariance value corresponding to the last elastic modulus parameter in the main diagonal element of , is denoted as l, and the total number of elastic modulus parameters is (l-m+1).

取步骤7.11计算的

Figure BDA0003598743550000093
的所有主对角元素组成一个新的对角方阵
Figure BDA0003598743550000094
且保持原对角元素的位置不变。依次扩大
Figure BDA0003598743550000095
中每一个弹性模量参数对应的协方差值，令
Figure BDA0003598743550000096
且每次只对一个协方差值进行扩大，
Figure BDA0003598743550000097
的其余元素值保持不变，其中
Figure BDA0003598743550000098
和
Figure BDA0003598743550000099
分别代表行数列数均为z位置处的协方差值，为一个标量，并且式中λ＝1×10ω，且10ω等于初始状态量协方差P0中最小数量级的倒数；z起始值为m，且z∈[m,l-m+1]。依据z＝m时计算的
Figure BDA00035987435500000910
和步骤7.10计算的
Figure BDA00035987435500000911
执行一步完整的无迹卡尔曼滤波运算，即执行步骤7.1～步骤7.11，并输出ηk，并令L(z)＝ηk。Take the calculated in step 7.11
Figure BDA0003598743550000093
All the main diagonal elements of form a new diagonal square matrix
Figure BDA0003598743550000094
And keep the position of the original diagonal element unchanged. expand in turn
Figure BDA0003598743550000095
The covariance value corresponding to each elastic modulus parameter in , let
Figure BDA0003598743550000096
And only one covariance value is expanded at a time,
Figure BDA0003598743550000097
The remaining element values of , remain unchanged, where
Figure BDA0003598743550000098
and
Figure BDA0003598743550000099
Respectively represent the covariance value at the z position where the number of rows and columns is a scalar, and where λ=1×10ω , and 10ω is equal to the reciprocal of the minimum order of magnitude in the initial state quantity covariance P0 ; z starts The value is m, and z∈[m,l-m+1]. Calculated according to z=m
Figure BDA00035987435500000910
and calculated in step 7.10
Figure BDA00035987435500000911
Perform a complete unscented Kalman filter operation, that is, perform steps 7.1 to 7.11, and output ηk , and let L(z)=ηk .

步骤7.15、令z＝m+1，继续执行步骤7.14，直到z＝l-m+1时结束。Step 7.15, set z=m+1, continue to perform step 7.14, and end when z=l-m+1.

步骤7.16、不计零值，找到L中最小值对应的位置号，记为zmin，即zmin为损伤位置。Step 7.16, ignoring the zero value, find the position number corresponding to the minimum value in L, and record it as zmin , that is, zmin is the damage position.

步骤7.17、仅令

Figure BDA00035987435500000912
中的元素
Figure BDA00035987435500000913
扩大20倍，其中
Figure BDA00035987435500000914
代表行数列数均为z位置处的协方差值，且倍数20来源于数值模拟统计结果。Step 7.17, only make
Figure BDA00035987435500000912
elements in
Figure BDA00035987435500000913
expanded 20 times, where
Figure BDA00035987435500000914
The number of representative rows and columns is the covariance value at the z position, and the multiple of 20 comes from the statistical results of numerical simulation.

步骤7.18、时间步变为(k+1)，重复步骤7.1～步骤7.17直到最大时间步N完成，即直到循环结束，找到所有损伤位置并输出。Step 7.18, the time step becomes (k+1), repeat steps 7.1 to 7.17 until the maximum time step N is completed, that is, until the end of the cycle, find all damage locations and output them.

实施例Example

为了充分说明本发明，本发明以车-桥耦合系统进行实施例的说明。In order to fully explain the present invention, the present invention is described by taking the vehicle-axle coupling system as an embodiment.

为了充分说明本发明，本实施例先对图4所示的车-桥耦合系统进行说明：In order to fully illustrate the present invention, this embodiment first describes the vehicle-axle coupling system shown in FIG. 4 :

所述的车-桥耦合系统，包括四分之一车辆模型的车体质量1、四分之一车辆模型的车体与轮胎之间的悬挂刚度2、四分之一车辆模型的车体与轮胎之间的悬挂阻尼3、四分之一车辆模型的轮胎质量4、四分之一车辆模型的轮胎与桥梁之间的接触刚度5、轮胎与桥梁接触点6、简支边界约束的固定端7、梁单元8、梁单元节点9和简支边界约束的滑动端10；The vehicle-axle coupling system includes thevehicle body mass 1 of the quarter vehicle model, thesuspension stiffness 2 between the vehicle body and the tires of the quarter vehicle model, the vehicle body and the tire of the quarter vehicle model. Suspension damping betweentires 3, tire mass forquarter vehicle model 4, contact stiffness between tire and bridge forquarter vehicle model 5, tire-bridge contact point 6, fixed end for simply supportedboundary constraints 7.Beam element 8, beam element node 9 and simply supported boundary constrained slidingend 10;

所述的轮胎与桥梁接触点6是指在车辆行驶过程中轮胎与桥梁始终密接，不发生分离；Thecontact point 6 between the tire and the bridge refers to that the tire and the bridge are always closely connected during the driving process of the vehicle, and no separation occurs;

所述的梁单元8不局限于图中示意的位置，图中一共有6个梁单元；Thebeam unit 8 is not limited to the position shown in the figure, and there are a total of 6 beam units in the figure;

所述的梁单元节点9不局限于图中示意的位置，其余具有相同形状符号的均是梁单元节点，包括桥梁两端位置也都含有梁单元节点。The beam element node 9 is not limited to the position shown in the figure, and the rest with the same shape and symbol are beam element nodes, including the positions at both ends of the bridge that also contain beam element nodes.

方法实施过程如下：The implementation process of the method is as follows:

1、基于车辆与结构动力相互作用理论，并结合有限元理论，推导车-桥耦合系统的运动控制微分方程，并基于车辆部分的运动控制微分方程推导接触力关系，即轮胎与桥面板间的相互作用力方程，由此当车辆通过桥梁时可计算出桥梁所受的外荷载激励。1. Based on the dynamic interaction theory between the vehicle and the structure, combined with the finite element theory, the motion control differential equation of the vehicle-bridge coupled system is derived, and the contact force relationship between the tire and the bridge deck is derived based on the motion control differential equation of the vehicle part. Interaction force equation, from which the external load excitation of the bridge can be calculated when the vehicle passes the bridge.

2、基于有限元理论，桥梁结构的位移和速度状态可通过各梁单元节点自由度的位移和速度表征，动态荷载作用前，认为各节点自由度的位移和速度都为0，而各梁单元的弹性模量初始值可基于桥梁结构的材料组成求得，由此初始状态量χ0已知，再基于卡尔曼滤波器原理，获得初始状态量协方差P0。2. Based on the finite element theory, the displacement and velocity state of the bridge structure can be characterized by the displacement and velocity of the degrees of freedom of each beam element. Before the dynamic load is applied, the displacement and velocity of each joint degree of freedom are considered to be 0, while the The initial value of the elastic modulus can be obtained based on the material composition of the bridge structure, so the initial state quantity χ0 is known, and then based on the Kalman filter principle, the initial state quantity covariance P0 is obtained.

3、荷载作用过程中，桥梁各梁单元节点的竖向位移可通过传感器或数值模拟手段采集或计算得到，由此观测值已知。3. During the loading process, the vertical displacement of each beam element node of the bridge can be collected or calculated by sensors or numerical simulation methods, and the observed value is known.

4、基于步骤1～7及上述初始信息可进行桥梁结构损伤位置的识别。4. Based onsteps 1 to 7 and the above initial information, the identification of the damage position of the bridge structure can be carried out.

为便于本发明的应用说明，此处对算法重点部分做进一步阐释，分别选择初始状态量χ0和中间向量L作进一步说明，重点选择初始状态量是因为它与初始状态量协方差P0的设置、中间向量L的设置以及最后识别的状态量都密切相关。根据前文阐述，本实施例基于欧拉-伯努利梁单元构建桥梁有限元模型，由于每个梁单元有2个节点，每个节点有2个自由度，而根据图4显示，本实施例共有6个梁单元，故共有7个节点，包含14个自由度。因此，本桥梁的位移和速度状态分别通过14个参数表征，再考虑6个弹性模量参数，故本实施例的状态量维度为(14+14+6＝34)。根据方法介绍，L的维度同样等于34，且其初始值都为0。同样，为方便表达识别效果，假设第30和第31位置处的弹性模量参数发生折减损伤，且状态量中第30和第31位置处的弹性模量分别对应桥梁的两个梁单元，具体编号可事先定义好。通过方法计算，可优先识别出位置30和31两处的状态量参数异常，进而通过梁单元编号可迅速确定损伤位置。In order to facilitate the application description of the present invention, the key parts of the algorithm are further explained here, and the initial state quantity χ0 and the intermediate vector L are respectively selected for furtherexplanation . The settings, the settings of the intermediate vector L, and the final identified state quantities are all closely related. According to the foregoing description, this embodiment builds a finite element model of a bridge based on Euler-Bernoulli beam elements. Since each beam element has 2 nodes, each node has 2 degrees of freedom, and as shown in FIG. 4 , this embodiment There are 6 beam elements in total, so there are 7 nodes in total, including 14 degrees of freedom. Therefore, the displacement and velocity states of the bridge are respectively represented by 14 parameters, and then 6 elastic modulus parameters are considered, so the state quantity dimension of this embodiment is (14+14+6=34). According to the method introduction, the dimension of L is also equal to 34, and its initial value is 0. Similarly, for the convenience of expressing the identification effect, it is assumed that the elastic modulus parameters at the 30th and 31st positions are damaged by reduction, and the elastic moduli at the 30th and 31st positions in the state quantity correspond to the two beam elements of the bridge respectively, The specific number can be defined in advance. Through the method calculation, the abnormal state parameters at positions 30 and 31 can be identified preferentially, and the damage position can be quickly determined by the beam element number.

具体实施方式二：Specific implementation two:

本实施方式为一种存储介质，所述存储介质中存储有至少一条指令，所述至少一条指令由处理器加载并执行以实现一种桥梁损伤位置的无损识别方法。This embodiment is a storage medium, the storage medium stores at least one instruction, and the at least one instruction is loaded and executed by a processor to implement a method for non-destructive identification of bridge damage locations.

应当理解为本实施方式所述的存储介质包括但不限于磁存储介质和光存储介质；所述磁存储介质包括但不限于RAM、ROM，以及其他硬盘、U盘等存储介质。It should be understood that the storage media described in this embodiment include but are not limited to magnetic storage media and optical storage media; the magnetic storage media include but are not limited to RAM, ROM, and other storage media such as hard disks and U disks.

具体实施方式三：Specific implementation three:

本实施方式为一种桥梁损伤位置的无损识别设备，所述设备包括处理器和存储器，所述存储器中存储有至少一条指令，所述至少一条指令由处理器加载并执行以实现一种桥梁损伤位置的无损识别方法。This embodiment is a non-destructive identification device for a bridge damage location, the device includes a processor and a memory, the memory stores at least one instruction, and the at least one instruction is loaded and executed by the processor to realize a bridge damage A non-destructive method for location identification.

应当理解为本实施方式所述的设备包括但不限于包括处理器和存储器的设备，还可以包括其他具有信息采集、信息交互、控制能功能的单元或模块所对应的设备，例如所述设备还可以包括信号采集装置等。所述设备包括但不限于PC机、工作站、移动设备等。It should be understood that the device described in this embodiment includes, but is not limited to, a device including a processor and a memory, and may also include other devices corresponding to units or modules with functions of information collection, information interaction, and control. For example, the device also It may include a signal acquisition device and the like. The devices include, but are not limited to, PCs, workstations, mobile devices, and the like.

以上显示和描述了本发明的基本原理和主要特征和本发明的优点。本行业的技术人员应该了解，本发明不受上述实施例的限制，上述实施例和说明书中描述的只是说明本发明的原理，在不脱离本发明精神和范围的前提下，本发明还会有各种变化和改进，这些变化和改进都落入要求保护的本发明范围内。本发明要求保护范围由所附的权利要求书及其等效物界定。The basic principles and main features of the present invention and the advantages of the present invention have been shown and described above. Those skilled in the art should understand that the present invention is not limited by the above-mentioned embodiments, and the descriptions in the above-mentioned embodiments and the description are only to illustrate the principle of the present invention. Without departing from the spirit and scope of the present invention, the present invention will have Various changes and modifications fall within the scope of the claimed invention. The claimed scope of the present invention is defined by the appended claims and their equivalents.

Claims (10)

1. A nondestructive identification method for a bridge damage position comprises the following steps:
aiming at the bridge structure, determining the initial value of the corresponding state of the bridge structure, and forming an initial state quantity chi0And determining covariance matrix of initial state quantity, called initial state quantity covariance P for short, according to Kalman filtering principle0(ii) a Wherein x0And P0State quantity and state quantity covariance, referred to as 0 th time step, respectively;
performing preliminary identification based on the unscented Kalman filter algorithm, and calculating the observation error epsilon of the kth time step by a measurement updating step based on the unscented Kalman filter algorithm in the process of performing preliminary identification based on the unscented Kalman filter algorithmkAnd measured predicted covariance P at kth time stepyy,kAnd based on epsilonkAnd Pyy,kCalculating and outputting sensitivity parameters corresponding to each step
Figure FDA0003598743540000011
Then drawing an output sensitivity parameter etakTime course curve of (1), if etakWhen the time curve has the peak pulse, identifying the damage position based on the algorithm of the self-adaptive unscented Kalman filterThe device is characterized in that:
in the process of identification based on the adaptive unscented Kalman filter algorithm, eta required to be judged and calculatedkWith a sensitivity parameter threshold η0If η is large or smallk＜η0Continuing to identify the target based on the unscented kalman filter algorithm; if etak≥η0Then, the following steps are continuously executed:
setting an initial value to 0 and dimension equal to initial state quantity chi0The vector of (a), denoted by the letter L; taking state quantity covariance
Figure FDA0003598743540000012
All main diagonal elements form a new diagonal matrix
Figure FDA0003598743540000013
And the position of the original diagonal element is kept unchanged; successively enlarge
Figure FDA0003598743540000014
The covariance value corresponding to each elastic modulus parameter is expanded, only one covariance value is expanded each time, one-step complete unscented Kalman filtering operation is performed, and one-step complete unscented Kalman filtering operation is performed to obtain a sensitive parameter eta aiming at the expansion of a row number or a column number z corresponding to the covariance valuekLet L (z) ═ ηk(ii) a Find the position number z corresponding to the minimum LminThe location is the site of injury; then order
Figure FDA0003598743540000015
Expanding by 20 times, continuing to execute the (k +1) th time step, and judging eta againk+1With a sensitivity parameter threshold η0Until the end of the whole cycle, all the lesion sites are obtained.
2. The method for the nondestructive identification of the damage position of the bridge according to claim 1, characterized in that: and in the process of identification by using the unscented Kalman filter algorithm and the process of identification by using the adaptive unscented Kalman filter algorithm, selecting the vertical displacement response of each beam unit node as an observed value y.
3. The method for non-destructive identification of a damage location of a bridge of claim 2, wherein: the corresponding state of the bridge structure comprises the elastic modulus of each beam unit of the bridge and the displacement and the speed of the bridge structure.
4. The method for the nondestructive identification of the damage position of the bridge as claimed in claim 1,2 or 3, wherein: calculating a sensitive parameter eta in the process of identification based on an adaptive unscented Kalman filter algorithmkComprises the following steps:
step 7.1, using the state quantity chi of the (k-1) th time step based on the UT conversion principle of the unscented Kalman filter algorithmk-1Sum state quantity covariance Pk-1Generating (2n +1) sigma points, and solving the state quantity corresponding to each sigma point through a state equation
Figure FDA0003598743540000021
Where k starts at 1 and k ∈ [1, N ]]N is the total number of time steps, N is the dimension of the state quantity, i is the ith sigma point, and i is an element [1,2N +1 ]]；
Step 7.2, the updating of the state quantity and the state quantity covariance from the (k-1) th time step to the k th time step is completed based on the time updating step of the unscented Kalman filter algorithm, and the updating is respectively recorded as
Figure FDA0003598743540000022
And
Figure FDA0003598743540000023
the formula is as follows:
Figure FDA0003598743540000024
Figure FDA0003598743540000025
in the formula,
Figure FDA0003598743540000026
and
Figure FDA0003598743540000027
respectively the weighted value of the ith sigma point at the kth time step,
Figure FDA0003598743540000028
for the state quantity estimated value, Q, corresponding to the ith sigma point of the kth time stepkNoise at the kth time step;
step 7.3, UT transformation principle based on unscented Kalman filter algorithm, using updated in step 7.2
Figure FDA0003598743540000029
And
Figure FDA00035987435400000210
generating (2n +1) sigma points, and solving an observation estimation value corresponding to each sigma point through an observation equation
Figure FDA00035987435400000211
Step 7.4, calculating and outputting a measurement predicted value of the kth time step in a measurement updating step based on the unscented Kalman filter algorithm
Figure FDA00035987435400000212
And is
Figure FDA00035987435400000213
In the formula,
Figure FDA00035987435400000214
the weight value of the ith sigma point at the kth time step,
Figure FDA00035987435400000215
an observation estimated value corresponding to the ith sigma point at the kth time step;
step 7.5, calculating the observation error epsilon of the kth time step based on the measurement updating step of the unscented Kalman filter algorithmkAnd is made of
Figure FDA00035987435400000216
In the formula, ykIs an observed value at the kth time step,
Figure FDA00035987435400000217
the measured predicted value of the kth time step is taken as a measurement predicted value;
step 7.6, calculating the measurement prediction covariance P of the kth time step based on the measurement updating step of the unscented Kalman filter algorithmyy,kAnd is and
Figure FDA00035987435400000218
in the formula,
Figure FDA00035987435400000219
the weight value of the ith sigma point at the kth time step,
Figure FDA00035987435400000220
is an observed estimated value corresponding to the ith sigma point at the kth time step,
Figure FDA00035987435400000221
is a measured prediction value of the kth time step, RkNoise at the kth time step;
step 7.7, ε calculated based on step 7.5 and step 7.6kAnd Pyy,kConstruction sensitive parameter ηkAnd is and
Figure FDA00035987435400000222
and calculates and outputs eta of each stepkThe value is obtained.
5. The method for the nondestructive identification of the damage position of the bridge according to claim 4, characterized in that: in the process of identification based on the adaptive unscented Kalman filter algorithm, a calculation sensitive parameter eta is obtainedkThen the k time step is calculated
Figure FDA0003598743540000031
And
Figure FDA0003598743540000032
cross covariance of Pxy,kAnd updating data, wherein the specific process comprises the following steps:
step 7.8, calculating the kth time step based on the measurement updating step of the unscented Kalman filter algorithm
Figure FDA0003598743540000033
And
Figure FDA0003598743540000034
cross covariance P ofxy,k，
Figure FDA0003598743540000035
Step 7.9, updating the Kalman gain matrix of the kth time step:
Figure FDA0003598743540000036
and 7.10, updating and outputting the state quantity of the kth time step:
Figure FDA0003598743540000037
and 7.11, updating and outputting the covariance of the state quantity of the kth time step:
Figure FDA0003598743540000038
6. the method for the nondestructive identification of the damage position of the bridge according to claim 5, characterized in that: find the position number z corresponding to the minimum LminThe specific process comprises the following steps:
step 7.14, covariance of state quantity at kth time step
Figure FDA0003598743540000039
The row number or the column number of the covariance value corresponding to the first elastic modulus parameter in the main diagonal elements is marked as m, and the state quantity covariance is marked as
Figure FDA00035987435400000310
The row number or the column number of the covariance value corresponding to the last elastic modulus parameter in the main diagonal elements is marked as l, and the total number of the elastic modulus parameters is (l-m + 1);
taking the current time step
Figure FDA00035987435400000311
All main diagonal elements form a new diagonal matrix
Figure FDA00035987435400000312
And the position of the original diagonal element is kept unchanged; successively enlarge
Figure FDA00035987435400000313
The covariance value corresponding to each elastic modulus parameter in the equation
Figure FDA00035987435400000314
And only one covariance value is expanded at a time,
Figure FDA00035987435400000315
the value of the remaining elements of (1) is kept constant, wherein
Figure FDA00035987435400000316
And
Figure FDA00035987435400000317
respectively representing the covariance values at the z position where the number of rows and the number of columns are respectively a scalar; lambda is the expansion factor;
z is initially m and z is e [ m, l-m +1 ]](ii) a Calculated according to z ═ m hours
Figure FDA00035987435400000318
And state quantity of kth time step
Figure FDA00035987435400000319
Executing one-step complete unscented Kalman filtering operation, namely executing the step 7.1 to the step 7.11, and outputting etakLet L (z) ═ ηk；
Step 7.15, let z be m +1, continue to execute step 7.14, and end until z is l-m + 1;
7.16, counting zero, finding a position number corresponding to the minimum value in the L and recording the position number as zminI.e. zminIs the site of the injury.
7. The method for the nondestructive identification of the damage position of the bridge according to claim 6, characterized in that: the expansion factor λ is 1 × 10ω，10ωEqual to initial state quantity covariance P0The reciprocal of the smallest order of magnitude.
8. The method for non-destructive identification of a damage location of a bridge of claim 7, wherein: the motion control differential equation or the finite element model corresponding to the bridge structure is established based on the Euler-Bernoulli beam unit.
9. A storage medium, characterized by: the storage medium has at least one instruction stored therein, and the at least one instruction is loaded and executed by a processor to implement the method for non-destructive identification of a bridge damage location according to any one of claims 1 to 8.
10. The utility model provides a nondestructive identification equipment of bridge damage position which characterized in that: the apparatus comprises a processor and a memory, the memory having stored therein at least one instruction, the at least one instruction being loaded and executed by the processor to implement a method of non-destructive identification of a bridge damage location according to any of claims 1 to 8.
CN202210398866.7A2022-04-152022-04-15 Non-destructive identification method, storage media and equipment for bridge damage locationActiveCN114692465B (en)

Priority Applications (1)

Application NumberPriority DateFiling DateTitle
CN202210398866.7ACN114692465B (en)2022-04-152022-04-15 Non-destructive identification method, storage media and equipment for bridge damage location

Applications Claiming Priority (1)

Application NumberPriority DateFiling DateTitle
CN202210398866.7ACN114692465B (en)2022-04-152022-04-15 Non-destructive identification method, storage media and equipment for bridge damage location

Publications (2)

Publication NumberPublication Date
CN114692465Atrue CN114692465A (en)2022-07-01
CN114692465B CN114692465B (en)2023-09-08

Family

ID=82142301

Family Applications (1)

Application NumberTitlePriority DateFiling Date
CN202210398866.7AActiveCN114692465B (en)2022-04-152022-04-15 Non-destructive identification method, storage media and equipment for bridge damage location

Country Status (1)

CountryLink
CN (1)CN114692465B (en)

Cited By (1)

* Cited by examiner, † Cited by third party
Publication numberPriority datePublication dateAssigneeTitle
CN116698323A (en)*2023-08-072023-09-05四川华腾公路试验检测有限责任公司Bridge health monitoring method and system based on PCA and extended Kalman filtering

Citations (8)

* Cited by examiner, † Cited by third party
Publication numberPriority datePublication dateAssigneeTitle
US5995620A (en)*1995-02-151999-11-30Telefonaktiebolaget Lm EricssonEcho canceller having Kalman filter for optimal adaptation
US20110084874A1 (en)*2009-09-102011-04-14ThalesHybrid system and device for calculating a position and for monitoring its integrity
CN102184554A (en)*2011-05-282011-09-14西安电子科技大学Infrared target tracking method based on background perception of activation region
CN109982398A (en)*2019-02-252019-07-05广州市香港科大霍英东研究院Indoor and outdoor seamless positioning method, device, equipment and medium
CN110398691A (en)*2019-06-262019-11-01重庆大学 Lithium-ion power battery SoC estimation method based on improved adaptive double unscented Kalman filter
CN112986977A (en)*2021-04-272021-06-18南京雷电信息技术有限公司Method for overcoming radar extended Kalman track filtering divergence
CN114186595A (en)*2021-12-142022-03-15哈尔滨工业大学Time-varying structure parameter identification method, storage medium and device
US20220146264A1 (en)*2020-11-102022-05-12Alpen-Adria-Universität KlagenfurtMethod and system for estimating state variables of a moving object with modular sensor fusion

Patent Citations (8)

* Cited by examiner, † Cited by third party
Publication numberPriority datePublication dateAssigneeTitle
US5995620A (en)*1995-02-151999-11-30Telefonaktiebolaget Lm EricssonEcho canceller having Kalman filter for optimal adaptation
US20110084874A1 (en)*2009-09-102011-04-14ThalesHybrid system and device for calculating a position and for monitoring its integrity
CN102184554A (en)*2011-05-282011-09-14西安电子科技大学Infrared target tracking method based on background perception of activation region
CN109982398A (en)*2019-02-252019-07-05广州市香港科大霍英东研究院Indoor and outdoor seamless positioning method, device, equipment and medium
CN110398691A (en)*2019-06-262019-11-01重庆大学 Lithium-ion power battery SoC estimation method based on improved adaptive double unscented Kalman filter
US20220146264A1 (en)*2020-11-102022-05-12Alpen-Adria-Universität KlagenfurtMethod and system for estimating state variables of a moving object with modular sensor fusion
CN112986977A (en)*2021-04-272021-06-18南京雷电信息技术有限公司Method for overcoming radar extended Kalman track filtering divergence
CN114186595A (en)*2021-12-142022-03-15哈尔滨工业大学Time-varying structure parameter identification method, storage medium and device

Non-Patent Citations (3)

* Cited by examiner, † Cited by third party
Title
周卫东;乔相伟;吉宇人;孟凡彬;: "基于新息和残差的自适应UKF算法", 宇航学报, no. 07*
王硕;宋申民, 于浛: "基于组合导航系统的非线性高斯滤波算法", 控制与决策, no. 09*
齐俊桐;韩建达;: "基于MIT规则的自适应Unscented卡尔曼滤波及其在旋翼飞行机器人容错控制的应用", 机械工程学报, no. 04*

Cited By (2)

* Cited by examiner, † Cited by third party
Publication numberPriority datePublication dateAssigneeTitle
CN116698323A (en)*2023-08-072023-09-05四川华腾公路试验检测有限责任公司Bridge health monitoring method and system based on PCA and extended Kalman filtering
CN116698323B (en)*2023-08-072023-10-13四川华腾公路试验检测有限责任公司Bridge health monitoring method and system based on PCA and extended Kalman filtering

Also Published As

Publication numberPublication date
CN114692465B (en)2023-09-08

Similar Documents

PublicationPublication DateTitle
Zenzen et al.A damage identification technique for beam-like and truss structures based on FRF and Bat Algorithm
CN103116759B (en)A kind of damnification recognition method being applicable to spatial mesh structure
CN102759573B (en)Based on the construction damage positioning of frequency change and the appraisal procedure of degree of injury
CN112001110A (en)Structural damage identification monitoring method based on vibration signal space real-time recursive graph convolutional neural network
Meruane et al.Structural damage assessment using linear approximation with maximum entropy and transmissibility data
Aval et al.A novel optimization algorithm based on modal force information for structural damage identification
CN115455764B (en)Frame structure multi-node damage diagnosis method
CN114048670B (en)Intelligent algorithm for rapidly evaluating damage of regional concrete bridge after earthquake
CN110008520B (en)Structural damage identification method based on displacement response covariance parameters and Bayesian fusion
CN111595541A (en)Multi-dimensional structure damage identification method based on mass vibration transmissibility data convolutional neural network processing
CN116881819B (en) A method for monitoring the working status of cable-stayed cables based on isolated forests
WO2021142901A1 (en)Method of physical mode extraction for engineering structure flexibility identification
CN104750926B (en)Non-uniform beam damnification recognition method based on node curvature and wavelet analysis
Ding et al.Jaya-based long short-term memory neural network for structural damage identification with consideration of measurement uncertainties
CN114638039A (en) A low-rank matrix recovery method for structural health monitoring feature data interpretation
CN114692465A (en)Nondestructive identification method of bridge damage position, storage medium and equipment
ALTABEY et al.Enhancement of structural health monitoring framework on beams based on k-nearest neighbor algorithm
Xiang et al.Vibration-based structural damage identification using P-CNN and time-frequency hybrid index under the conditions of uncertainties and incomplete measurements
Farhadi et al.Damage detection of truss structures using meta-heuristic algorithms and optimized group method of data handling surrogate model
CN105651537B (en)A kind of truss structure damage real-time monitoring system of high susceptibility to damage
Keshmiry et al.AI-based structural health monitoring systems
Wang et al.Structural damage detection using autoregressive-model-incorporating multivariate exponentially weighted moving average control chart
Yan et al.Comparative Analysis of Structural Damage Identification Methods Based on Iterative Reweighted L 1/2 Regularization and Three Optimization Functions
CN112818577B (en)Method for identifying post-fire damage of laminated beam based on deep learning theory
Andreadis et al.Intelligent seismic acceleration signal processing for damage classification in buildings

Legal Events

DateCodeTitleDescription
PB01Publication
PB01Publication
SE01Entry into force of request for substantive examination
SE01Entry into force of request for substantive examination
GR01Patent grant
GR01Patent grant
[8]ページ先頭
©2009-2025 Movatter.jp