机械强度

组号 Group number	第1组 1st group	第2组 2nd group	第3组 3rd group	第4组 4th group	第5组 5th group
杆件编号 Bar number	1、2、3、10、11	1、4、12、13、18	13、14、15、19、21	5、6、16、19、20	6、7、8、9、17
弹性模量 Elastic modulus	0.8E₀	E₀	1.2E₀	E₀	1.6E₀
密度Density	ρ₀	ρ₀	1.4ρ₀	ρ₀	1.8ρ₀
泊松比 Poisson ratio	0.28	0.28	0.28	0.28	0.28
横截面 Cross section /m²	0.002 5	0.002 5	0.002 5	0.002 5	0.002 5

组号 Group number	第1组 1st group	第2组 2nd group	第3组 3rd group	第4组 4th group	第5组 5th group
杆件编号 Bar number	1、2、3、10、11	1、4、12、13、18	13、14、15、19、21	5、6、16、19、20	6、7、8、9、17
弹性模量 Elastic modulus	0.8E₀	E₀	1.2E₀	E₀	1.6E₀
密度Density	ρ₀	ρ₀	1.4ρ₀	ρ₀	1.8ρ₀
泊松比 Poisson ratio	0.28	0.28	0.28	0.28	0.28
横截面 Cross section /m²	0.002 5	0.002 5	0.002 5	0.002 5	0.002 5

模态 Modal	修正前 Before updating	MH算法 MH algorithm	本文所提算法 Algorithm proposed in this paper
1	4.630	4.549	1.784	4.591	0.842	4.632	0.039
2	19.621	21.942	10.581	19.475	0.740	19.654	0.170
3	25.646	27.865	7.965	25.543	0.404	25.699	0.204
4	38.494	40.383	4.678	38.288	0.533	38.574	0.209
5	55.538	58.061	4.345	54.557	1.767	55.546	0.016
6	59.524	63.410	6.128	58.776	1.257	59.603	0.133

模态 Modal	修正前 Before updating	MH算法 MH algorithm	本文所提算法 Algorithm proposed in this paper
1	4.630	4.549	1.784	4.591	0.842	4.632	0.039
2	19.621	21.942	10.581	19.475	0.740	19.654	0.170
3	25.646	27.865	7.965	25.543	0.404	25.699	0.204
4	38.494	40.383	4.678	38.288	0.533	38.574	0.209
5	55.538	58.061	4.345	54.557	1.767	55.546	0.016
6	59.524	63.410	6.128	58.776	1.257	59.603	0.133

参数 Parameter	修正前 Before updating	MH算法MH algorithm	本文所提算法Algorithm proposed in this paper
θ₁	0.8	0.778	0.022	2.700	0.807	0.007	0.913
θ₂	1.0	0.989	0.011	1.070	0.998	0.002	0.190
θ₃	1.2	1.130	0.070	5.842	1.196	0.004	0.300
θ₄	1.4	1.374	0.026	1.850	1.398	0.002	0.179
θ₅	1.6	1.518	0.082	5.100	1.590	0.010	0.631
θ₆	1.8	1.756	0.044	2.422	1.790	0.010	0.533

参数 Parameter	修正前 Before updating	MH算法MH algorithm	本文所提算法Algorithm proposed in this paper
θ₁	0.8	0.778	0.022	2.700	0.807	0.007	0.913
θ₂	1.0	0.989	0.011	1.070	0.998	0.002	0.190
θ₃	1.2	1.130	0.070	5.842	1.196	0.004	0.300
θ₄	1.4	1.374	0.026	1.850	1.398	0.002	0.179
θ₅	1.6	1.518	0.082	5.100	1.590	0.010	0.631
θ₆	1.8	1.756	0.044	2.422	1.790	0.010	0.533

基于马尔科夫链种群竞争的贝叶斯有限元模型修正

PDF下载

叶玲 ¹ , 江宏康 ¹ , 邹雨清 ¹ , 陈华鹏 ¹ , 王力骋 ²

机械强度 | 优化·可靠性 2025,47(2): 85-93

收起

机械强度 | 优化·可靠性 2025, 47(2): 85-93

基于马尔科夫链种群竞争的贝叶斯有限元模型修正

全屏

叶玲¹, 江宏康¹, 邹雨清¹, 陈华鹏¹, 王力骋²

作者信息

^1.华东交通大学轨道交通基础设施性能监测与保障国家重点实验室，南昌　330013

^2.江西交投咨询集团有限公司，南昌　330013

叶玲，女，1989年生，湖北荆门人，博士，副教授，硕士研究生导师；主要研究方向为结构健康监测、轨道动力学等；E-mail：58718070@qq.com。

Bayesian finite element model updating based on Markov chain population competition

Ling YE¹, Hongkang JIANG¹, Yuqing ZOU¹, Huapeng CHEN¹, Licheng WANG²

Affiliations

^1.State Key Laboratory of Performance Monitoring Protecting of Rail Transit Infrastructure, East China Jiaotong University, Nanchang 330013, China

^2.Jiangxi Transport Investment Consulting Group Co., Ltd., Nanchang 330013, China

出版时间: 2025-02-15 doi: 10.16579/j.issn.1001.9669.2025.02.011

文章导航

摘要

收起

针对传统马尔科夫链蒙特卡洛（Markov Chain Monte Carlo，MCMC）模拟方法在高维问题或后验概率密度复杂时采样效率低且难收敛的缺陷，建立了基于马尔科夫（Markov）链种群竞争的贝叶斯有限元模型修正算法。在基于Metropolis-Hastings（MH）随机游走算法实现MCMC模拟的传统方法基础上，引入差分进化算法，利用种群中Markov链之间不同携带信息的相互作用关系，得到优化建议以快速逼近目标函数，解决了高维参数模型修正过程中采样滞留的缺点；引进竞争算法，通过不断的竞争刺激和内置失败者向胜利者学习的机制，采用较少的Markov链获得较高的精度，提高了模型修正效率与精度；最后，通过一个桁架结构的有限元模型修正数值算例验证了所提算法，并与标准MH算法的结果对比，得出该算法可以快速修正高维参数模型，具有较高的精度，且对随机噪声有良好的鲁棒性，为考虑不确定性的大型结构有限元模型修正提供了一种稳定有效的手段。

关键词

模型修正 / 贝叶斯估计 / 马尔科夫链蒙特卡洛 / 种群竞争

Abstract

收起

The traditional Markov chain Monte Carlo (MCMC) simulation method is inefficient and difficult to converge in high dimensional problems and complicated posterior probability density.In order to overcome these shortcomings, a Bayesian finite element model updating algorithm based on Markov chain population competition was proposed. First, the differential evolution algorithm was introduced in the traditional method of Metropolis-Hastings (MH) random walk algorithm.Based on the interaction of different information carried by Markov chains in the population, optimization suggestions were obtained to approach the objective function quickly. It solves the defect of sampling retention in the updating process of high-dimensional parameter model. Then, the competition algorithm was introduced, which has constant competitive incentives and a built-in mechanism for losers to learn from winners. Higher precision was obtained by using fewer Markov chains, which improves the efficiency and precision of model updating. Finally, a numerical example of finite element model updating of a truss structure was used to verify the proposed algorithm.Compared with the results of standard MH algorithm, the proposed algorithm can quickly update the high-dimensional parameter model with high accuracy and good robustness to random noise. It provides a stable and effective method for finite element model updating of large-scale structure considering uncertainty.

Key words

Model updating / Bayesian estimation / Markov chain Monte Carlo / Population competition

引用本文

叶玲, 江宏康, 邹雨清, 陈华鹏, 王力骋. 基于马尔科夫链种群竞争的贝叶斯有限元模型修正. 机械强度, 2025 , 47 (2) : 85 -93 . DOI: 10.16579/j.issn.1001.9669.2025.02.011

Ling YE, Hongkang JIANG, Yuqing ZOU, Huapeng CHEN, Licheng WANG. Bayesian finite element model updating based on Markov chain population competition[J]. Journal of Mechanical Strength, 2025 , 47 (2) : 85 -93 . DOI: 10.16579/j.issn.1001.9669.2025.02.011

正文

收起

0　引言

收起

在结构健康监测领域中，随着建立精确有限元模型评估结构服役状态的方法不断完善，有限元法逐渐成为最广泛、最常用、最基本的数值分析手段。但由于实际工程结构存在材料的参数误差、边界条件的不明确、外部荷载的变化、环境因素以及测量噪声等影响，造成实际工程问题具有来源未知、难以测量的不确定性^[1-2]。因此，已有的一些确定性模型修正方法往往处于理论效果好，但实际应用效果不佳的尴尬位置。为了解决不确定性造成的理论与实际之间的差距，采用概率或统计手段量化与结构相关的不确定性，已逐渐成为建立有限元模型关注的重点。其中，数学统计分析中贝叶斯学派以可充分考虑不确定性而闻名。

近年来，大量国内外学者开始研究基于贝叶斯理论的概率模型修正方法，原因在于贝叶斯推理过程具有可同时考虑结构的先验信息和当前测量值的优势^[3]。因此，其逐步成为模型修正领域的重点研究方向，实现了多种工程结构的模型修正和健康监测。其中BECK等^[4]首次将贝叶斯理论应用于工程结构模型修正，并建立了修正的基本流程体系，通过Laplace渐进法近似得到后验概率分布对两自由度的框架修正；同时，KATAFYGIOTIS等^[5]提出了一种基于Metropolis-Hastings（MH）随机游走算法的马尔科夫链蒙特卡洛（Markov Chain Monte Carlo, MCMC）模拟方法，通过模拟采样近似得到后验概率分布，并以两自由度模型验证了方法的正确性；易伟建等^[6]采用MH法对混凝土框架结构进行修正，编制了MCMC方法的损伤识别程序；WAN等^[7]提出一种延缓拒绝自适应Metropolis(Delayed Rejection Adaptive Metropolis, DRAM)算法，并对一座复杂的斜拉人行桥修正，结果表明DRAM算法拥有较好的样本遍历性；GREEN^[8]针对非线性问题，提出了一种模拟退火思路的新式MCMC模拟方法，有效地提高了模型修正的效率；CHING等^[9]从模型类的选择出发，采取避开复杂概率密度函数的采样方法，提出了一种过渡马尔科夫链蒙特卡洛(Transitional Markov Chain Monte Carlo,TMCMC）算法，通过3个例子说明了TMCMC方法在贝叶斯模型更新、模型类选择和模型平均方面的有效性；ZENG等^[10]利用振动响应或模态参数来估计结构参数和估计本身相关的不确定性，采用差分进化自适应Metropolis (Differential Evolutionary Adaptive Metropolis, DREAM)算法通过并行运行多个Markov链抽取样本提高接受率，并通过十层剪力建筑和三层框架结构的数值示例验证了算法的有效性；WU等^[11]为了解决Expectation‑Maximization（EM）算法因时间复杂导致的计算困难，以MCMC采样用潜在变量的模拟样本逼近EM算法，提出了MCMC-EM算法，并验证了其比传统EM算法的先进性；HUANG等^[12]将贝叶斯概率法和微扰模型相结合的方法，该方法可仅利用少量有噪声干扰的节点振动响应，对结构进行健康监测；刘纲等^[13]针对MCMC方法面对高维待修正参数收敛难效率低问题，引进相关向量机作为回归模型代替有限元计算提高效率；SHERRI等^[14]针对不确定性量化问题提出一种改进的进化MCMC方法，并使用斯诺克更新器改进和扩展更新过程，通过一个F型框架结构验证了方法的有效性；彭珍瑞等^[15]通过新鸟巢更新思想改进标准MH算法，同时使用支持向量机（Support Vector Machines,SVM）建立代理模型，有效地解决了MH算法不收敛、拒绝率高的问题，并在一个三自由度平面桁架上验证了方法的可行性。

由上述研究可以看出，贝叶斯理论的核心内容主要是通过查阅资料或根据已有先前经验得到先验分布，采用当前实测信息来修正先验分布，使其不断接近后验概率分布，从而有效处理参数不确定性。但是后验分布的求解难度较高，原因在于其分母涉及复杂且难以求解的高维积分，通过常规的数值积分手段效果甚微，引进MCMC采样算法，以近似解等效替代精确解，因此采样算法的精确度将直接影响修正结果。目前已有的标准MH采样算法，对单个参数的模型采样效果较好，在多参数的采样过程中容易出现大面积滞留现象，同时样本采样拒绝率过高会导致收敛困难，甚至无法收敛。针对以上缺陷，本文在标准MH算法基础上引入种群竞争算法和差进化算法，可有效地改善高维参数模型修正的不足，并通过一个桁架结构数值算例验证了本文算法的高效性。

1　贝叶斯有限元模型修正理论

收起

采用贝叶斯方法描述有限元模型修正过程中，不确定性参数分布通过后验概率密度函数的数值求解得到，其根据贝叶斯方法具体表述为^[16][17]^115-132

（1）

式中，M为系统的结构模型类别（每个结构类别的定义的修正参数不同）；θ为不确定参数；D为结构体系的实测数据样本（本文采用固有频率和振型）；π（θ，M）为先验概率密度分布，表示不确定参数θ在给定的模型M中且不含D时的初始认知；P（D|θ，M）为似然函数，表示实测数据与有限元解之间的差异；π（θ|D，M）为给定模型M和实测数据D的修正参数的后验概率分布函数；由于本文只研究单个模型，后续公式推导中将M省略，且P（D，M）与修正参数θ不相关，可用常数c表示。因此，π（θ|D）与P（D|θ）π（θ）存在正比例关系为

（2）

（3）

同时，假设结构实测数据与有限元模型的计算值存在以下关系：

（4）

（5）

式中，N为用于修正的结构响应量阶数；θ为修正参数组集的向量（θ∈R^N）；Z^*为实测获得的结构响应量组集的向量（Z^*∈R^N）；Z（θ）为有限元模型计算值组集的向量[Z（θ）∈R^N]；ε为测试波动、环境变化等不确定性引起的差值组集的向量（ε∈R^N）；其中，R^N为对应阶数的实数集范围，R为实数集,上标N为用于修正的结构响应量阶数。故后验概率分布函数表示为

（6）

文中采用测试简易且数值较为精准的结构固有频率和振型作为模型修正的响应量，式（6）可以进一步推导为

（7）

（8）

式中，N_s为实测的有效响应次数；N_n为响应模态对应阶数；ω_ij为第j次i阶模型的实测频率值；φ_ij为第j次i阶模型的实测振型值；ω（θ）、φ（θ）分别为有限元模型得到的频率和振型计算值；J（θ）为实测值与计算值之间的差值函数。

简化结构模型的修正参数后验概率分布函数表达式，计算π（θ|D）最大值可转化为目标函数的最小值为

（9）

贝叶斯修正过程知之非难，行之不易，要想取得良好效果需解决以下两个难题：①结构的复杂性越大，后验概率密度函数解析解的获取越困难；②未知参数的维数以及样本空间的复杂性会显著影响修正过程，且不确定参数的全局解不容易获得。因此，贝叶斯修正过程往往与MCMC方法相结合，通过MCMC方法以采样代替复杂高维积分计算，近似获取π（θ|D），从而得到近似的后验分布函数。

2　基于种群的多链差分竞争算法

收起

传统MCMC方法只能得到单条Markov链，如应用广泛的标准MH算法通过随机游走采样生成一条Markov链，依据Markov链中样本分布近似得到后验分布，而其采样机制的受限于固定的建议分布方差，若方差较大则易导致建议点时常被拒绝沿用上次迭代值，遍历性差出现滞留现象；方差较小则遍历步距小、收敛速度慢，甚至出现迭代结束而样本未更新，导致该方法在面对复杂结构或后验概率密度复杂的情况时，存在采样效率低且精度不佳的缺陷。本节通过建立Markov链种群，以多链角度引入差分算法以克服滞留和难收敛缺陷。在此基础上，加入竞争算法，通过种群中多链差分竞争的方式选择更优修正方案，提出一种基于种群内部竞争的算法以提高采样效率及修正精度，从而高效逼近目标分布。

2.1　Markov链种群

种群竞争需要在多条Markov链中进行，本文采用标准MH算法，对同一目标函数设定不同建议分布方差，采样得到多条Markov链，建立Markov链种群[Y¹，Y²，…，Y^N]。在理论层面Markov链中N越大种群携带测试信息越丰富效果越好，但综合考虑种群竞争计算性价比本文设定N=7，其原因在于刚好能满足完成一次完整的种群竞争，而交叉、变异采用的Markov链不重复使用，建议分布方差在范围[0.1，1.0]依据Markov链数均匀取值。

同时，本方法建立Markov链种群的方法不局限于标准MH算法，其他能生成有效Markov链的采样算法同样适用。建立Markov链种群的具体步骤如下：

步骤1：确定种群Markov链数量N=7及建议分布方差，根据先验概率分布初始化修正参数θ₀和迭代步j=1，且满足P（θ|D_Ns）>0。

步骤2：在每次迭代j时，根据建议分布q（θ^* |θ_j₋₁）生成一个新的样本θ^*。

步骤3：修正有限元模型和计算后验概率分布函数p（θ^*|Z）。

步骤4：计算接受概率。计算式为

（10）

步骤5：从均分布U[0，1]中随机提取数值u，如果u<α（θ^*，θ_j₋₁），则接受θ^*；否则，拒绝θ^*。

步骤6：重复步骤2~5，直到达到最大迭代次数N_s，得到一个收敛序列

。

步骤7：复步骤1~6，直到获得N条Markov链，完成种群建立。

2.2　多链差分算法

多链差分算法类似遗传进化均包含变异、交叉和选择3个操作运算，但与遗传进化不同的是变异这一操作运算通过利用种群中个体与个体之间的差量实现条件分布的更新^[18]。本文中个体为已构建的种群中的Markov链，由于构建每条链的目标函数相同，但每条单链的初始建议分布，链内接受概率不同，以及多链差分过程中链间变异、交叉概率的不同，导致部分交换有助于进一步改善Markov链的质量，且多链相对单链包含更多模态空间信息，通过选择操作运算接受可提升部分，可有效改善传统单链面对高维模态空间能力受限问题。其中多链差分算法变异具体操作方式如下：

（11）

式中，Y为一条Markov链；上标s、a、b为该条Markov链的索引，其中s=1、2、3、…、N，a、b∈｛1，2，…，s+1，N｝；上标p表示经过变异操作得到的向量，即索引为s原生向量

的变异向量记为

，此时第s条Markov链的变异向量

；下标t为当前的进化（迭代）步数；γ为变异比例因子，取值越大突变程度越高，一般在范围[0，0.4]中随机抽取；e为从正态分布N（0，k）^d中提取的向量组，其中d为修正向量的维数，k一般选取1×10^-6。

结合突变后得到的变异向量

和原生向量

，应用交叉生成修正向量

，二项式交叉方式可表达为

（12）

式中，Y为Markov链一组修正向量中的元素；上标s为该条Markov链的索引；上标u表示修正向量；上标w为每条Markov链在当前t迭代步时包含的修正参数的索引，w=1、2、3、…、D，其中D为每条Markov链待修正参数个数；P_CR为交叉概率，一般在范围[0.0，1.0]中随机抽取；w_rand为[1，D]中的一个随机数；下标t为当前的迭代步数。

经过变异和交叉得到的修正向量

和原生向量

结合式（9）代入目标函数F，即可得到接受概率r₁，计算式为

（13）

最后，判别接受概率r₁与1的大小选择是否接受

，如果拒接，则当前迭代步得到的修正向量沿用原生向量

，选择方式为

（14）

因此，以Markov链为个体的差分算法基于进化MCMC方法，利用多链之间的相互作用关系得到新的优化建议，可解决高维参数模型修正过程中的滞留缺陷，提高修正精度和效率。

2.3　种群差分竞争算法

种群竞争算法是一种激励算法，在多链差分算法基础上，采用不断竞争与学习机制提高差分进化效率，充分利用种群中Markov链携带信息，采用较少的Markov链获得较高的精度，是一种基于竞争寻找解决目标函数最优方案的新型改进差分算法。

竞争的具体思路是：Markov链在多链差分算法变异，交叉操作的基础上融入竞争，在选择操作中定义处于拒绝状态的为失败者，反之为胜利者；通过竞争算法提升多链差分的质量，即失败者通过选取胜利者Markov链的部分，试图学习胜利者的优势进而缩小与目标函数的差距得到最优方案；同时这又导致胜利者优势只是暂时的领先，多链差分算法中的变异、交叉、选择3项操作也反过来辅助竞争算法，使胜利者继续突破寻找更佳的解，失败者在快速学习胜利者的基础上仍有超过胜利者的可能，胜利者与失败者的定位在激烈的竞争中随时可能交换。

由于种群中Markov链的前部分处于燃烧阶段，如一开始便加入竞争会影响算法的收敛速度，因此竞争算法是人为假设在m次迭代后出现竞争直到迭代结束的整个过程，因此迭代阶段按该算法思路，可根据是否含有竞争划分为竞争期和正常期（竞争期多链差分算法和竞争算法均存在，正常期只含有多链差分算法）。其中，正常期算法具体思路与操作步骤如下：

步骤1：开始t=1~m时，迭代处于正常期，待修正参数更新方法是在Markov链种群中随机一个链记为

，根据式（11）突变得到

。

步骤2：根据式（12）将变异向量

与原生向量

交叉得到修正向量

。

步骤3：根据式（9）和式（13）计算接受概率r。

步骤4：依据式（14）选择修正向量

。

步骤5：当t=m时，正常期迭代结束。

竞争期的算法是本文算法的核心，具体思路与操作步骤如下：

步骤1：当t=m~T时，迭代从正常期进入竞争期，待修正参数更新方法是在Markov链种群中另随机选取一条Markov链记为

，并根据式（15）、式（16）区分出胜利者

与失败者

：

（15）

（16）

式中，W为胜利者；L为失败者。

步骤2：若得到的是失败者

进入步骤3，若是胜利者

进入步骤4。

步骤3：对

的更新方法：

1）失败者

在胜利者基础上进行变异操作，得到变异向量记为

，变异表达式为

（17）

式中，γ′与上文γ含义相同，取值不同，均为变异比例因子，取值越大突变程度越高，一般在范围[0，0.5]中随机抽取；上标c、d为Markov链的索引。

2）变异向量

结合胜利者

进行交叉操作，交叉表达式为

（18）

3）修正向量

与原生向量

结合式（9）计算接受概率r₂，计算式为

（19）

4）选择经过步骤3优化后的

，失败者修正的选择计算式为

（20）

步骤4：对

的更新方法：

1）胜利者

继续与种群中新选取的Markov链进行变异操作，得到变异向量记为

，胜利者变异表达式为

（21）

2）

结合种群中新选取的Markov链进行交叉操作，交叉表达式为

（22）

3）修正向量

与原生向量

结合式（9）计算接受概率r₃，计算式为

（23）

4）选择经过步骤4优化后的

，胜利者修正的选择表达式为

（24）

步骤5：竞争期每次迭代中，经过变异、交叉、选择后得到的两组修正向量

和

（此时，t=m~T），分别和正常期的修正向量

结合（此时，t=1~m），得到两组完整迭代后Markov链参数优化方案（此时，t=1~T），记为Y^W和Y^L。

步骤6：计算两组优化方案Y^W和Y^L收敛后的参数平均值记为

和

，并根据式（9）比较

与

，选出最优修正方案。

步骤7：重复上述步骤1~6，获得收敛后的最优方案模型，即得到修正参数的后验概率分布。

此外由于Markov链燃烧阶段的样本是未收敛的，若以此为修正模型参数会降低精度，根据经验公式^[17]^115-132可知，只要迭代步数满足T≥300d，得到的样本个数对于修正后参数值的精确几乎不造成影响，因此，本文迭代步数设定为T=1×10⁴>300d，迭代步数充分满足修正要求，其中，正常期与竞争期的m虽是人为假设，但其选值需大于Markov链收敛分界点10%的迭代步数，若取值太小，在Markov链收敛分界点之内，将会影响算法的收敛；若取值太大，超过收敛分界点太多，将会消除竞争作用，延长正常期，影响修正精度和效率。所提方法的具体步骤流程如图1所示。

3　数值算例

收起

3.1　平面桁架

某平面桁架结构如图2所示，以此算例检验所提算法的高效性。该桁架为对称斜三角桁架，由21个杆件铰接组成，共有12个节点、21个单元和21个自由度，如图3所示。按杆件单元位置将其划分为5个小组，其中编号1、2、3、10、11号单元为第1组；编号1、4、12、13、18号单元为第2组；编号13、14、15、19、21号单元为第3组；编号5、6、16、19、20号单元为第4组；编号6、7、8、9、17号单元为第5组。选取第1、3、5组的弹性模量和密度为模型的待修正参数，一共6个待修正参数，分别为θ₁ =E₁ /E₀，θ₂ =ρ₂ /ρ₀，θ₃ =E₃ /E₀，θ₄=ρ₄ /ρ₀，θ₅ =E₁ /E₀，θ₆=ρ₆ /ρ₀，其中桁架基数分别为E₀ =210 GPa，ρ₀ =7 800 kg/m³。桁架组号划分内容及结构参数如表1所示。

本文假定6个试验模型参数服从正态分布N｛[0.8，0.1²]，[1.0，0.1²]，[1.2，0.1²]，[1.4，0.1²]，[1.6，0.1²]，[1.8，0.1²]｝，且参数之间相互独立。通过概率统计手段，随机抽取一组待修正参数代入桁架有限元模型，计算得到前6阶模态，为了充分考虑实际工程测试中的不确定性因素，得到更接近实测情况的测试数据，本文共模拟计算35次，以得到的35次数据均值为结构实测响应。

对样本进行独立检验，由图4可知，6个参数35次模拟得到的样本点绝大部分处于正态分布的假设概率密度山形图的山峰表面，表明随机采样得到样本点满足试验独立假设要求，所得到的样本点具有足够的代表性。

将不确定参数θ₁~θ₆均除以其对应的初始值

（k=1、2、3、4、5、6），得到归一化数据后绘制置信椭圆图（误差椭圆图）如图5所示。图中的点为迭代1×10⁴次后得到的不确定性参数样本，椭圆区域表示95%置信度的样本空间，可以看出不确定性参数基本落位在置信椭圆之内。

图6给出了标准MH算法和本文所提算法的收敛曲线对比图。其中，误差定义为6个待修正参数的绝对误差之和，其计算式为

（25）

如图6所示，在收敛速度上，可明显观察到本文所提算法在第200次迭代步后基本收敛，而MH算法在400次迭代步才逐渐收敛；在收敛精度上，也可看出本文所提算法的误差更小；在收敛路径上，MH算法呈锯齿状，表明采样空间未充分探索，而本文所提算法小范围波动，具备良好的遍历性。

图7给出了MH算法和本文所提算法的修正对比图。由图7分析可知，当完成1×10⁴迭代步时两种算法均已收敛，但MH算法得到的Markov链6个参数均明显出现大块区域的滞留现象，在较高维模型修正中样本的采样效率极低；而本文所提算法融入差分思想有效改善高维参数模型修正滞留缺点，同时收敛速度较快。再根据图7（a）、图7（b）对参数θ₄和θ₅局部对比放大可知（图7右侧图），本文所提算法迭代得到的参数样本滞留区域较小、接受率和采样效率较高、采样具有良好的遍历性，且融入的种群竞争和差分步骤使修正在稳定后仍在收敛值附近小幅度起伏波动，进而提高模型修正精度，相比MH算法可更好地模拟出接近真实情况的修正参数后验分布的样本空间。

提取前800次迭代步的修正过程如图8所示。处于正常期时，以Markov链为个体的差分算法有效地改善了传统MCMC采样效率低且难收敛的问题；200次迭代步加入竞争算法进入竞争期后，差分算法与竞争算法相互提供激励条件，可观察到6个修正参数收敛后的波动幅度更小，进一步提升了模型修正精度。

理论上对比修正后模型得到的响应与实测响应是否一致可判断是否达到修正要求，但在研究中仅以此来判断模型修正效果是片面的，还需结合修正参数值。因此，本文给出了不同算法修正前后模态响应和参数的对比数据表（表2、表3）。从表中可以看出：从频率响应角度分析，首先，初始模型响应与真实响应的前6阶频率均存有5%左右的误差，修正前第2阶误差最大，达到10.58%。标准MH算法修正后前6阶频率误差均值为1%左右，修正后最大误差为1.77%，而本文所提算法修正后频率误差均值在0.21%以内，且修正后最大误差仅为0.21%；从修正前后参数对比分析，MH算法表现远不如本文所提算法，其中θ₃和θ₅修正后相对误差分别为5.84%和5.1%，而本文所提算法的各参数相对误差均在1%以内。对比分析表明本文所提算法识别结构参数精度较高，对高维参数模型修正效果远胜于标准MH算法。

3.2　加噪试验

为了使桁架模型修正更加合乎实际，以及进一步验证所提方法的高效性，本文对试验得到的测试数据加入10%的噪声，用来模拟实际工程问题中测试遇到的不确定性因素（如实测人员操作规范问题、设备装置精度问题等）。为对比不同算法在同一噪声条件下的修正效果，根据式（26）随机加入10%的测试噪声，且假定模态噪声在不同阶和不同参数之间相互独立。

（26）

式中，F₀同式（9）中的F，表示未加噪之前的目标函数；ε为噪声幅值。

图9（b）为Markov链和每个修正参数的概率直方图。首先对采样过程图分析可知，加噪10%对本文所提算法的收敛性影响不大，Markov链均可收敛，且每个参数的概率直方图表明，加噪10%后的参数也大致在真实值±0.04E₀/ρ₀左右，虽略有差距但在可允许接受范围之内，表明了本文所提算法具有良好的鲁棒性。对加噪后得到的修正参数实施正态概率检验，得到结果如图10所示。从图10可知，加噪后计算得到的样本点基本位于假设参数分布线上，表明采样得到的样本点满足正态分布假设，结果可充分接受（对纵坐标刻度进行了调整,以便观察正态率,不影响检验数据是否近似服从正态分布）。

因此，结合上述分析可知，本文所提算法可以在考虑各种不确定性因素下，可精准识别结构模型参数变化，高效地完成有限元模型修正，且算法具备较高的修正精度和良好的鲁棒性。

4　结论

收起

针对不确定性有限元模型修正问题，以贝叶斯统计理论为核心，在标准MH算法基础上引进了竞争算法和差分算法，提出了一种基于Markov链种群竞争差分进化采样的改进算法，通过一个六参数桁架结构模型修正数值算例验证了所提算法的精度和效率，且面对10%的测试噪声干扰下仍保持一定精度，并与传统MH算法计算结果对比，得出以下结论：

1）本文所提算法面对高维参数模型修正采样接受率、遍历性更好，通过提高采样效率的方式极大地提升了模型修正速率，实现了高维参数的不确定性有限元模型修正。

2）本文所提算法通过种群竞争和差分手段，以Markov链之间携带信息差异为建议优化方向，进一步提高了算法对高维参数模型修正精度。

3）本文所提算法在10%的测试噪声干扰下，仍能保证较高精度，得到可靠的识别结果，表明所提算法具备良好的鲁棒性。

本研究目前使用数值算例验证了所提算法的高效性，在保证采样效率的同时具有较高的精度和稳定性。后续将陆续开展结构模态试验研究，为实际结构中考虑不确定性的有限元模型修正提供一种有效稳定的手段。

基金

收起

国家自然科学基金项目(52008168; 52468042)
国家重点研发计划项目(2021YFE015600)

参考文献

收起

文献

收起

参考文献引证文献

排序方式：

[1]

BIONDINI

，FRANGOPOL

D M

. Life-cycle performance of deteriorating structural systems under uncertainty：review[J].Journal of Structural Engineering，2016，142（9）：1-17.

[2]

SIMOEN

，ROECK

G D

，LOMBAERT

. Dealing with uncertainty in model updating for damage assessment：a review[J].Mechanical Systems and Signal Processing，2015（56）：123-149.

[3]

YUEN

K V

. Bayesian methods for structural dynamics and civil engineering[M]. Hoboken：John Wiley and Sons，2010：12-24.

[4]

BECK

J L

，KATAFYGIOTIS

L S

. Updating models and their uncertainties I：Bayesian statistical framework[J].Journal of Engineering Mechanics，1998，124（4）：455-461.

[5]

KATAFYGIOTIS

L S

，BECK

J L

. Updating models and their uncertainties Ⅱ：Model identifiability[J].Journal of Engineering Mechanics，1998，124（4）：463-467.

[6]

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

Weijian

，ZHOU

Yun

，LI

Hao

.Damage assesssment research on frame structure based on Bayesian statistical inference[J].Engineering Mechanics，2009，26（5）：121-129. （In Chinese）

[7]

WAN

H P

，REN

W X

. Stochastic model updating utilizing Bayesian approach and Gaussian process model[J].Mechanical Systems and Signal Processing，2016（70）：245-268.

[8]

GREEN

P L

. Bayesian system identification of a nonlinear dynamical system using a novel variant of simulated annealing[J].Mechanical Systems and Signal Processing，2015（52）：133-146.

[9]

CHING

，WANG

J S

. Discussion of “transitional Markov chain Monte Carlo：observations and improvements” by wolfgang betz，iason papaioannou，and daniel straub[J].Journal of Engineering Mechanics，2017，143（9）：1-3.

[10]

ZENG

J C

，KIM

Y H

. Probabilistic damage detection and identification of coupled structural parameters using Bayesian model updating with added mass[J].Journal of Sound and Vibration，2022（539）：117275.

[11]

，MA

J W

. An effective EM algorithm for mixtures of Gaussian processes via the MCMC sampling and approximation[J].Neurocomputing，2019（331）：366-374.

[12]

HUANG

T X

，SCHRODER

K U

. IWSHM 2019：perturbation-based Bayesian damage identification using responses at vibration nodes[J].Structural Health Monitoring，2021，20（3）：942-959.

[13]

刘纲，罗钧，秦阳，等.基于改进MCMC方法的有限元模型修正研究[J].工程力学，2016，33（6）：138-145.

LIU

Gang

，LUO

Jun

，QIN

Yang

，et al. A finite element model updating method based on improved MCMC method[J].Engineering Mechanics，2016，33（6）：138-145. （In Chinese）

[14]

SHERRI

，BOULKAIBET

，MARWALA

，et al. Bayesian finite element model updating using a population Markov chain Monte Carlo algorithm[M]//Special Topics in Structural Dynamics & Experimental Techniques，Volume 5. Cham：Springer International Publishing，2020：259-269.

[15]

彭珍瑞，郑捷，白钰，等.一种基于改进MCMC算法的模型修正方法[J].振动与冲击，2020，39（4）：236-245.

PENG

Zhenrui

，ZHENG

Jie

，BAI

，et al. A model updating method based on an improved MCMC algorithm[J].Journal of Vibration and Shock，2020，39（4）：236-245. （In Chinese）

[16]

BOULKAIBET

，MTHEMBU

，MARWALA

，et al. Finite element model updating using Hamiltonian Monte Carlo techniques[J].Inverse Problems in Science and Engineering，2017，25（7）：1042-1070.

[17]

BOULKAIBET

，MTHEMBU

，MARWALA

，et al. Finite elment model updating using the shadow hybrid Monte Carlo technique[J].Mechanical Systems and Signal Processing，2015（52）：115-132.

[18]

TER BRAAK

C J F

. A Markov chain Monte Carlo version of the genetic algorithm differential evolution：easy Bayesian computing for real parameter spaces[J].Statistics and Computing，2006，16（3）：239-249.

2025年第47卷第2期

PDF下载

引用本文

BibTeX

文章信息

doi: 10.16579/j.issn.1001.9669.2025.02.011

接收时间：2023-05-24
首发时间：2026-03-18
出版时间：2025-02-15

补充材料

相关文章

文章信息

作者

出版历史

收稿日期：2023-05-24
修回日期：2023-07-17

基金

National Natural Science Foundation of China(52008168; 52468042)

国家自然科学基金项目(52008168; 52468042)

National Key Research and Development Program of China(2021YFE015600)

国家重点研发计划项目(2021YFE015600)

作者信息

^1.华东交通大学轨道交通基础设施性能监测与保障国家重点实验室，南昌　330013

^2.江西交投咨询集团有限公司，南昌　330013

参考文献

分享链接

https://castjournals.cast.org.cn/joweb/jxqd/CN/10.16579/j.issn.1001.9669.2025.02.011

分享至

全文二维码

扫描看全文

引用本文

BibTeX

本文的引用情况

2种不同金属材料的力学参数

科 Family	属数 Number of genus	种数 Number of species	占总种数比例 Percentage of total species (%)	属 Genus	种数 Number of species	占总种数比例 Percentage of total species (%)
鹅膏菌科Amanitaceae	2	11	5.26	鹅膏菌属 Amanita	10	4.78
小菇科 Mycenaceae	2	12	5.74	丝盖伞属 Inocybe	5	2.39
多孔菌科 Polyporaceae	8	14	6.70	蜡蘑属 Laccaria	5	2.39
红菇科 Russulaceae	3	23	11.00	小皮伞属 Marasmius	6	2.87
				小菇属 Mycena	11	5.26
				光柄菇属 Pluteus	5	2.39
				红菇属 Russula	17	8.13
				栓菌属 Trametes	5	2.39

关闭全屏

BibTeX
EndNote
RefWorks
TxT

组号 Group number	第1组 1st group	第2组 2nd group	第3组 3rd group	第4组 4th group	第5组 5th group
杆件编号 Bar number	1、2、3、10、11	1、4、12、13、18	13、14、15、19、21	5、6、16、19、20	6、7、8、9、17
弹性模量 Elastic modulus	0.8E₀	E₀	1.2E₀	E₀	1.6E₀
密度Density	ρ₀	ρ₀	1.4ρ₀	ρ₀	1.8ρ₀
泊松比 Poisson ratio	0.28	0.28	0.28	0.28	0.28
横截面 Cross section /m²	0.002 5	0.002 5	0.002 5	0.002 5	0.002 5

组号
Group number

第1组
1st group

第2组
2nd group

第3组
3rd group

第4组
4th group

第5组
5th group

杆件编号
Bar number

1、2、3、10、11

1、4、12、13、18

13、14、15、19、21

5、6、16、19、20

6、7、8、9、17

弹性模量
Elastic modulus

0.8E₀

E₀

1.2E₀

E₀

1.6E₀

密度Density

ρ₀

1.4ρ₀

ρ₀

1.8ρ₀

泊松比
Poisson ratio

0.28

横截面
Cross section /m²

0.002 5

模态 Modal	修正前 Before updating	MH算法 MH algorithm	本文所提算法 Algorithm proposed in this paper
1	4.630	4.549	1.784	4.591	0.842	4.632	0.039
2	19.621	21.942	10.581	19.475	0.740	19.654	0.170
3	25.646	27.865	7.965	25.543	0.404	25.699	0.204
4	38.494	40.383	4.678	38.288	0.533	38.574	0.209
5	55.538	58.061	4.345	54.557	1.767	55.546	0.016
6	59.524	63.410	6.128	58.776	1.257	59.603	0.133

模态
Modal

修正前
Before updating

MH算法
MH algorithm

本文所提算法
Algorithm proposed in this paper

真实
Actual/Hz

初始
Initial/Hz

误差
Error/%

修正
Updated/Hz

误差
Error/%

修正
Updated/Hz

误差
Error/%

4.630

4.549

1.784

4.591

0.842

4.632

0.039

19.621

21.942

10.581

19.475

0.740

19.654

0.170

25.646

27.865

7.965

25.543

0.404

25.699

0.204

38.494

40.383

4.678

38.288

0.533

38.574

0.209

55.538

58.061

4.345

54.557

1.767

55.546

0.016

59.524

63.410

6.128

58.776

1.257

59.603

0.133

参数 Parameter	修正前 Before updating	MH算法MH algorithm	本文所提算法Algorithm proposed in this paper
θ₁	0.8	0.778	0.022	2.700	0.807	0.007	0.913
θ₂	1.0	0.989	0.011	1.070	0.998	0.002	0.190
θ₃	1.2	1.130	0.070	5.842	1.196	0.004	0.300
θ₄	1.4	1.374	0.026	1.850	1.398	0.002	0.179
θ₅	1.6	1.518	0.082	5.100	1.590	0.010	0.631
θ₆	1.8	1.756	0.044	2.422	1.790	0.010	0.533

参数
Parameter

修正前
Before updating

MH算法MH algorithm

本文所提算法Algorithm proposed in this paper

修正后Updated

误差 Error

修正后Updated

误差 Error

绝对Absolute

相对Relative

绝对Absolute

相对Relative

θ₁

0.8

0.778

0.022

2.700

0.807

0.007

0.913

θ₂

1.0

0.989

0.011

1.070

0.998

0.002

0.190

θ₃

1.2

1.130

0.070

5.842

1.196

0.004

0.300

θ₄

1.4

1.374

0.026

1.850

1.398

0.002

0.179

θ₅

1.6

1.518

0.082

5.100

1.590

0.010

0.631

θ₆

1.8

1.756

0.044

2.422

1.790

0.010

0.533