CN104376231A - 基于改进近似贝叶斯计算的损伤识别方法 - Google Patents

Abstract

本发明涉及一种基于改进近似贝叶斯计算的损伤识别方法。首先，利用近似贝叶斯计算使得参数后验概率分布的求解无需计算参数的似然函数，解决了实际工程应用中似然函数往往无法求解的难题；蒙特卡罗马尔科夫链抽样先建立与求解问题具有相似性的概率模型，再对该模型进行随机抽样，然后利用抽取的样本求出其统计特征估计值，并作为原问题的近似解。但对复杂工程问题来说，蒙特卡罗马尔科夫链抽样方法所需的计算量往往非常大，因此本发明在抽样过程中利用随机响应面快速计算参数样本所对应的结构响应统计特征值，避免了调用有限元模型进行数值求解，因而大幅提高了计算效率，解决了贝叶斯方法在多参数、大样本量下由于计算量过大而无法实现的问题。

Description

基于改进近似贝叶斯计算的损伤识别方法

技术领域

本发明涉及一种基于改进近似贝叶斯计算的损伤识别方法。

背景技术

传统实际工程结构长期处于复杂的运营环境和外荷载作用中，不可避免地会发生不同程度和类型的损伤。随着时间的增长，已有损伤会不断累积，导致结构性能不断退化。若不能及时发现损伤并采取有效的加固措施，在极端条件下就可能导致结构发生灾难性事故。损伤识别作为结构健康监测系统的核心内容^[1]，是近些年来相关领域的研究热点之一。当前已有的损伤识别方法可以归为确定性和不确定性方法两大类^[2-4]。确定性方法中参数和响应都被看作是确定值，所构建的损伤识别指标也是确定性的，仅能有效地应用在实验室试验模型或精确控制的特定结构上，目前还未有能够在实际工程中通用的可靠方法^[5]。因此在现实中，理论上可行的方法在实际应用时的效果往往很不理想，最终导致工程人员不得不重新依赖费时费力的常规检测技术。因而考虑到工程结构的构件本身、所处的工作环境以及所受到的外荷载都不可避免地存在不确定性因素（如构件几何尺寸误差、材料离散性、边界和连接条件的不确定、荷载随机性、环境噪声等），在损伤识别过程中结合概率统计分析方法是目前的发展趋势，其也可以看作是对确定性方法的拓展。

目前已有的概率损伤识别方法主要采用了经典概率统计分析理论^[6]、贝叶斯理论^[7-12]或随机有限元方法^[13-16]。经典概率统计方法基于现有的样本观测值，通过构建适当的估计量和假设检验方法来求取未知参数的统计值。其在选取检验统计量时往往很困难，且无法利用参数的先验知识，也不考虑后续样本所能提供的信息，因此在具体应用上带有一定的局限性。随机有限元法主要通过对实测数据或模型参数进行摄动式的随机模拟来获取参数的概率统计特征，其中模型误差和测量噪声对参数的影响可以通过对实测数据的统计平均来降低。然而当参数摄动范围较大时，摄动法的精度会明显降低。同时对常用的一阶灵敏度摄动而言，其分析结果往往是局部收敛的，并且参数初值的选取对结果有着很大影响，因此在应用上具有较大的局限性。贝叶斯方法结合了参数的先验信息（具体表现为先验概率分布）及当前的实测数据，再基于最优概率模型来确定参数的后验概率分布。其主要优点是可以充分利用历史数据或专家经验（即先验信息），并结合当前实测数据来不断更新参数的概率分布，这一点非常符合在线健康监测和损伤识别的需求。但参数后验概率分布公式中的正则化常数往往无法求解，需要采用马尔可夫链蒙特卡罗方法^[17]来求得后验分布的近似解。当结构模型比较复杂且包含较多未知参数时，贝叶斯修正问题的求解也会变得十分困难且计算量大幅增加。更常见的问题是无法获取似然函数的表达式，导致贝叶斯方法的实用性大打折扣。为此，对传统贝叶斯方法进行改进，以求提出适用于工程实际问题的贝叶斯损伤识别方法，是本发明的主要目的。

总体而言，贝叶斯方法由于能充分利用已有数据和新的测试数据，以此不断修正参数的后验概率分布估计，这点对在线损伤识别是非常有利的。因此，若能对传统贝叶斯方法进行改进，解决似然函数求解困难的问题，并提高样本响应的计算效率，对工程实际应用来说具有重要的理论和实用意义。

发明内容

本发明的目的在于提供一种基于改进近似贝叶斯计算的损伤识别方法。

为实现上述目的，本发明的技术方案是：一种基于改进近似贝叶斯计算的损伤识别方法，首先，根据专家经验或历史数据假设结构随机参数的先验概率分布；其次，基于概率配点法和回归分析建立联系结构随机参数和响应的随机响应面；再次，基于参数先验概率分布随机生成初始值，通过转换函数进行采样，并利用随机响应面快速计算结构响应的统计特征值，根据目标函数及设置的阈值判断是否接受样本；然后，进一步计算样本的接受概率，以确定是否最终接受该样本；随后不断重复前两步，实现循环抽样，直到形成一条稳定的马尔科夫链，并根据包含在该链里的所有参数样本来计算参数的后验概率分布；最后，根据估计的参数后验概率分布来构建损伤指标，判断结构是否发生损伤。

在本发明实施例中，该方法的具体实现步骤如下：

步骤S1：根据历史数据或专家经验假设结构随机参数                                                （）的先验概率分布，以建立参数的初始贝叶斯模型；其间若无法确定先验分布的类型，可以先假设为均匀分布；

步骤S2：将随机参数用标准随机变量表示，并基于概率配点法和回归分析建立随机响应面，表现为联系随机参数和响应R的显式表达式：

                     (1)

式中，为待定系数；n为标准正态随机变量的个数；为p阶多维Hermite多项式：

                 (2)

步骤S3：基于随机生成的一个初始值，并利用该初始值从转换函数随机抽取1个参数样本；然后利用随机响应面快速计算所对应的结构随机响应的统计特征值；最后通过构建目标函数并选取目标函数阈值，判断的目标函数是否小于，若小于则进入步骤S4，否则重新抽取一个的初始值；

所构建的目标函数如下：

                                   (3)

式中，为随机响应面计算的第j阶响应，共有m阶；为实测的第j阶响应；

步骤S4：计算的接受概率，以确定是否最终接受；若满足要求，则接受并修正和；否则，回到步骤S3；

                         (4)

步骤S5：重复步骤S3、S4抽取N个样本，直到最终得到一条稳定的马尔科夫链，然后计算包含在该链中所有参数样本的后验概率分布（），作为各参数的最终分布；

步骤S6：分别估计结构在不同状态下的，然后构建概率损伤指标（damage index），以此判断损伤的位置和程度；

                                            (5)

式中，和分别表示无损和损伤结构参数的后验概率分布值。

相较于现有技术，本发明具有以下有益效果：

（1）避免了计算复杂的似然函数，大大减小了问题求解的复杂性，提高了求解效率并增强了实用性；

（2）抽样过程中，直接通过随机响应面计算随机参数所对应的结构响应统计特征值，无需进行复杂的数值分析计算，大大提高了计算效率。

附图说明

图1为本发明基于改进近似贝叶斯计算的损伤识别流程图。

具体实施方式

下面结合附图，对本发明的技术方案进行具体说明。

本发明结合了近似贝叶斯计算^[18]、蒙特卡罗马尔科夫链抽样^[17]和随机响应面^[19]3种方法，提出了一种改进的贝叶斯损伤识别方法。首先利用近似贝叶斯计算使得参数后验概率分布的求解过程无需计算参数的似然函数，解决了实际工程应用中似然函数往往无法求解的难题。而蒙特卡罗马尔科夫链抽样先建立与所求解问题具有相似性的概率模型，再对该模型进行随机模拟或统计抽样，然后利用抽取的样本求出其统计特征估计值，并作为原问题的近似解。但对复杂工程问题来说，蒙特卡罗马尔科夫链抽样方法所需的计算量往往非常大，因此在抽样过程中利用随机响应面^[12-13]快速计算参数样本所对应的结构响应统计特征值，避免了调用有限元模型进行数值求解，因而大幅提高了计算效率，解决了贝叶斯方法在多参数、大样本量下由于计算量过大而无法实现的问题。

如图1所示，本发明所采用的技术方案主要包括以下几个步骤：

步骤1：首先根据历史数据或专家经验假设结构参数（）的先验概率分布,以建立参数的初始贝叶斯模型。其间若无法确定先验分布的类型，可以先假设为均匀分布；

步骤2：将随机参数用标准随机变量表示，并基于概率配点法和回归分析建立随机响应面，表现为联系随机参数和响应R的显式表达式：

              (1)

式中等为待定系数；n为标准正态随机变量的个数；为p阶多维Hermite多项式：

          (2)

步骤 3：基于随机生成的一个初始值，并利用该初始值从转换函数（transition function, 一般假定为正态分布或均匀分布）随机抽取 1个参数样本；然后利用随机响应面快速计算所对应的结构随机响应的统计特征值；最后通过构建目标函数（即为响应误差函数）并选取合适的目标函数阈值，判断的目标函数是否小于，若小于则进入步骤4，否则重新抽取一个的初始值。

所构建的目标函数如下：

                           (3)

式中为随机响应面计算的第j阶响应（共有m阶）；为实测的第j阶响应。

步骤 4：计算的接受概率，以确定是否最终接受。若满足要求，则接受并修正和。否则，回到步骤3。

                (4)

步骤 5：重复步骤 3、4 以抽取N个样本， 直到最终得到一条稳定的马尔科夫链，然后计算所有样本的后验概率分布（），作为各参数的最终分布。

步骤 6：分别估计结构在不同状态下（无损、损伤）的，然后构建概率损伤指标（damage index），以此判断损伤的位置和程度。

                                   (5)

式中和分别表示无损和损伤结构参数的后验概率分布值。

由于土木工程结构在测量噪声、建模误差等因素的影响下，实测响应数据具有很强的不确定性。因而在进行结构损伤识别时，确定性损伤识别方法的识别结果必然不理想，和结构实际损伤情况会有出入。而当前的一些基于不确定性损伤识别方法大都是计算具有复杂性或计算效率低，实际应用上比较困难。 

本发明提出了一种基于近似贝叶斯方法的结构损伤识别方法，其优点在于（1）避免了计算复杂的似然函数，大大减小了问题求解的复杂性，提高了求解效率并增强了实用性；（2）抽样过程中，直接通过随机响应面计算随机参数所对应的结构响应统计特征值，无需进行复杂的数值分析计算，大大提高了计算效率。

以上是本发明的较佳实施例，凡依本发明技术方案所作的改变，所产生的功能作用未超出本发明技术方案的范围时，均属于本发明的保护范围。

本发明中涉及的参考文献如下：

[1] C.R. Farrar, K. Worden, An introduction to structural health monitoring, Philosophical Transactions of The Royal Society A (2007) 365:303-315.

[2] 宗周红, 任伟新, 阮毅. 土木工程结构损伤诊断研究进展[J]. 《土木工程学报》, 2003, 36(5): 105-110.

[3] 冯新, 李国强, 周晶. 土木工程结构健康诊断中的统计识别方法综述[J]. 《地震工程与工程振动》, 2005, 25(2): 105-113.

[4] Y.J. Yan, L. Cheng, Z.Y. Wu, L.H. Yam, Development in vibration-based structural damage detection technique, Mechanical Systems and Signal Processing (2007) 21(5): 2198-2211.

[5] 段忠东, 闫桂荣, 欧进萍, 土木工程结构振动损伤识别面临的挑战, 《哈尔滨工业大学学报》 (2008) 

40(4) :505-513.

[6] 张清华. 基于概率可靠度的结构损伤识别理论研究及应用 [D]. 成都: 西南交通大学, 2006.

[7] J.L. Beck, L.S. Katafygiotis, Updating models and their uncertainties. Part I: Bayesian statistical 

framework, Journal of Engineering Mechanics, ASCE (1998) 124(4):455-61.

[8] Beck J L, Au S K. Bayesian updating of structural models and reliability using Markov chain Monte Carlo simulation[J]. Journal of Engineering Mechanics, 2002, 128(4): 380-391.

[9] S.H. Cheung, J.L. Beck, Bayesian model updating using hybrid Monte Carlo simulation with application to structural dynamic models with many uncertain parameters, Journal of Engineering Mechanics, ASCE (2009) 135(4):243-255. 

[10] 易伟建, 周云, 李浩. 基于贝叶斯统计推断的框架结构损伤诊断研究[J]. 工程力学, 2009, 26(5): 121-129.

[11] J.M. Nichols, E.Z. Moore, K.D. Murphy, Bayesian identification of a cracked plate using a population-based Markov Chain Monte Carlo method, Computers & Structures (2011) 89(13-14):1323-1332.

[12] H.Y. Guo, Z.L. Li, Structural damage identification based on Bayesian theory and improved immune genetic algorithm, Expert Systems with Applications (2012) 39(7):6426-6434. 

[13] 安伟光, 朱卫兵. 随机有限元法在不确定性分析中的应用[J]. 哈尔滨工程大学学报, 2002, 23(1): 132-135.

[14] K. D’Souza, B.I. Epureanu, Multiple augmentations of nonlinear systems and generalized minimum rank perturbations for damage detection, Journal of Sound and Vibration (2008) 316(1-5):101-121. 

[15] 陈淮, 杜思义, 魏泽丽, 基于摄动 Riccati传递矩阵的梁桥损伤概率识别方法, 《振动工程学报》 (2010) 23(3):283-289.

[16] H. Xu, L. Cheng, Z.Q. Su, J.L. Guyader, Damage visualization based on local dynamic perturbation: theory and application to characterization of multi-damage in a plane structure, Journal of Sound and Vibration (2013) 332(14): 3438-3462.

[17] W.R. Gilks, S. Richardson, D. Spiegelhalter, Markov Chain Monte Carlo in Practice, CRC Press, 1995.

[18] B.M. Turner, T.V. Zandt, A tutorial on approximate Bayesian computation, Journal of Mathematical Psychology (2012) 56(2):69-85.

[19] S.S. Isukapalli, A. Roy, P.G. Georgopoulos, Stochastic response surface methods (SRSMs) for uncertainty propagation: application to environmental and biological systems, Risk Analysis (1998) 18(3):351-63.。

Claims (2)

1.一种基于改进近似贝叶斯计算的损伤识别方法，其特征在于：首先，根据专家经验或历史数据假设结构随机参数的先验概率分布；其次，基于概率配点法和回归分析建立联系结构随机参数和响应的随机响应面；再次，基于参数先验概率分布随机生成初始值，通过转换函数进行采样，并利用随机响应面快速计算结构响应的统计特征值，根据目标函数及设置的阈值判断是否接受样本；然后，进一步计算样本的接受概率，以确定是否最终接受该样本；随后不断重复前两步，实现循环抽样，直到形成一条稳定的马尔科夫链，并根据包含在该链里的所有参数样本来计算参数的后验概率分布；最后，根据估计的参数后验概率分布来构建损伤指标，判断结构是否发生损伤。

2.根据权利要求1所述的基于改进近似贝叶斯计算的损伤识别方法，其特征在于：该方法具体实现步骤如下：

步骤S1：根据历史数据或专家经验假设结构随机参数                                                （）的先验概率分布，以建立参数的初始贝叶斯模型；其间若无法确定先验分布的类型，可以先假设为均匀分布；

步骤S2：将随机参数用标准随机变量表示，并基于概率配点法和回归分析建立随机响应面，表现为联系随机参数和响应R的显式表达式：

                     (1)

式中，为待定系数；n为标准正态随机变量的个数；为p阶多维Hermite多项式：

                 (2)

步骤S3：基于随机生成的一个初始值，并利用该初始值从转换函数随机抽取1个参数样本；然后利用随机响应面快速计算所对应的结构随机响应的统计特征值；最后通过构建目标函数并选取目标函数阈值，判断的目标函数是否小于，若小于则进入步骤S4，否则重新抽取一个的初始值；

所构建的目标函数如下：

                                   (3)

式中，为随机响应面计算的第j阶响应，共有m阶；为实测的第j阶响应；

步骤S4：计算的接受概率，以确定是否最终接受；若满足要求，则接受并修正和；否则，回到步骤S3；

                         (4)

步骤S5：重复步骤S3、S4抽取N个样本，直到最终得到一条稳定的马尔科夫链，然后计算包含在该链中所有参数样本的后验概率分布（），作为各参数的最终分布；

步骤S6：分别估计结构在不同状态下的，然后构建概率损伤指标（damage index），以此判断损伤的位置和程度；

                                            (5)

式中，和分别表示无损和损伤结构参数的后验概率分布值。