梁崗　沈奎雙　馬玉付

摘要：

針對在結構損傷診斷過程中存在的不確定性和對稱結構損傷參數(shù)難以識別的問題，以梁損傷前后兩階頻率變化平方比為基礎損傷指標，采用對基礎損傷指標進行積分處理的方法，構建新的損傷指標?；谪惾~斯結構損傷診斷理論，建立損傷參數(shù)的后驗概率分布。采用馬爾科夫鏈蒙特卡洛（Markov chain Monte Carlo， MCMC）方法解決貝葉斯方程中存在的高維積分問題。仿真和算例分析說明，該方法可實現(xiàn)對損傷參數(shù)的有效估計。

關鍵詞：

貝葉斯推斷； 簡支梁； 裂紋； 損傷識別； 馬爾科夫鏈蒙特卡洛（Markov chain Monte Carlo， MCMC）方法

中圖分類號： U653.921； TH215

文獻標志碼： A

Abstract：

In order to solve the problems of the uncertainty in the process of structural damage diagnosis and the difficulty in identifying damage parameters of symmetrical structure， the two order frequencychange square ratio of beam before and after damage is taken as the basic damage index， and a new damage index is constructed by integrating the basic damage index. The posterior probability distribution of damage parameters is established based on the Bayesianian structural damage diagnosis theory. The Markov chain Monte Carlo （MCMC） method is adopted to solve the problem of high dimensional integral in Bayesianian equations. Simulation and example analysis show that the method can effectively estimate damage parameters.

Key words：

Bayesianian inference； simply supported beam； crackle； damage identification； Markov chain Monte Carlo （MCMC） method

0引言

隨著我國經(jīng)濟的高速發(fā)展，港口碼頭的建設取得了很大的成就，港口起重機在起重運輸領域中的作用日益增強。橋式起重機作為應用廣泛的起重運輸設備，在日常工作中要承受高載荷以及頻繁而劇烈的振動沖擊，其主梁結構經(jīng)常發(fā)生裂紋損傷現(xiàn)象，安全運行存在極大的隱患。因此，為保證起重機的安全運行并對其健康狀況進行監(jiān)測，對主梁結構進行損傷識別研究尤為重要。

結構損傷會使結構的動力特性發(fā)生變化，在結構損傷前后分別測試結構的模態(tài)值，利用損傷指標方法能有效識別出結構的損傷。由于工程中存在著許多不確性定因素，例如溫度的變化、力振幅的變化、動力測試噪聲等，實測值與真實值之間總存在著誤差，從而導致?lián)p傷識別問題成為不確定性問題[1]。因此，必須在確定性損傷識別研究的基礎上，建立一種能比較合理地反映不確定性損傷識別概率的方法。BECK等[2]和KATAFYGIOTIS等[3]于1998年提出了基于貝葉斯模型修正及統(tǒng)計推斷的基本框架，在此框架下VANIK等[4]建立了基于概率統(tǒng)計的健康監(jiān)測方法。易偉建等[5]和李小華等[6]基于貝葉斯基本原理和

馬爾科夫鏈蒙特卡洛（Markov chain Monte Carlo， MCMC）方法，對多層框架結構的局部加強的損傷情況進行了結構損傷診斷研究。房長宇等[7]考慮模型誤差和測試噪聲的影響，在基于動力響應的損傷識別中引入了貝葉斯估計理論，對軸力作用下的帶損傷混凝土梁段進行了損傷識別。

因為結構振動頻率的測量精度遠高于振型的測量精度，而且頻率的測量方法比較簡單，所以本文通過測試結構振動頻率的改變來識別結構的損傷。以單裂紋簡支梁為例，基于傳遞矩陣法[8]推導出裂紋梁的頻率特征方程，由此求出裂紋梁的低階固有頻率。以梁損傷前后兩階頻率變化平方比為基礎損傷指標，對簡支梁結構在對稱位置的裂紋損傷難以正確識別的問題[9]，采用對基礎損傷指標進行積分處理的方法，構建新的損傷指標。基于貝葉斯推斷理論分別建立損傷位置參數(shù)和損傷深度參數(shù)的后驗概率分布，針對貝葉斯方程中分母維數(shù)高、積分困難的問題，采用MCMC方法，得到損傷位置參數(shù)和損傷深度參數(shù)的最優(yōu)估計值，驗證該方法的有效性。裂紋梁損傷參數(shù)識別研究思路框圖見圖1。

1裂紋梁頻率特征方程推導

L，截面寬度和高度分別為b和h，其外側邊緣距離梁左端L0處有一深度為a的裂紋。彈性模量為E，截面慣性矩為I，密度為ρ，泊松比為ν。下面用傳遞矩陣法推導單裂紋簡支梁的頻率特征方程。

2貝葉斯推斷和MCMC方法

2.1貝葉斯推斷

貝葉斯推斷是貝葉斯統(tǒng)計的核心，離散型隨機變量的貝葉斯公式如下：

2.2馬爾科夫過程

在式（2）中，

h（x）與樣本觀測值x有關，一般的蒙特卡洛模擬難以有效計算。MCMC方法通過在蒙特卡洛模擬中引進Metropolis準則，將模擬過程看成一個馬爾科夫過程，模擬采樣最終收斂于式（2）定義的概率分布。定義隨機變量x的狀態(tài)轉移概率

為

P（i，j）=P（i→j）=P（Xt+1=sjXt=si）

它只取決于隨機變量取值的當前狀態(tài)，與過程無關。以此轉移概率定義的隨機變量x將最終達到一個靜態(tài)分布π，與x的初始狀態(tài)無關，滿足π=πP，此時存在局部平衡：

P（j→k）π*j=P（k→j）π*k

它是π存在的充分條件。

2.3MetropolisHastings （MH）抽樣算法

MH抽樣算法[11]的基本思想是：假定要從目標概率密度函數(shù)

p（θ）（-∞<θ<+∞）中抽樣，用MH抽樣算法可以得到一條馬爾科夫鏈，產(chǎn)生一條真值序列：

θ（1）→θ（2）→…→θ（t）→…

其中，θ（t）代表馬爾科夫鏈中時刻t的狀態(tài)。樣本開始收斂前的迭代次數(shù)稱為燃燒期，當經(jīng)過k步燃燒期以后，抽樣得到的序列可以近似反映目標分布p（θ），即通過平穩(wěn)的馬爾科夫鏈來近似等效替代未知的后驗分布。具體的步驟如下：

（1）令t=1，產(chǎn)生初始值θ（0），令θ（t）=θ（0）。

（2）利用當前θ（t）值，從建議分布q（θ*|θ（t））中產(chǎn)生一個候選參數(shù)θ*。

（3）計算接受率：

ε=min1，p（θ*）p（θ（t））

（4）從均勻分布U（0，1）中產(chǎn)生一個隨機數(shù)u。

（5）判斷是否接受候選參數(shù)：若ε≥u則接受候選參數(shù)，即θ（t+1）=θ*；若ε

（6）重復步驟（2）～（5），直至產(chǎn)生一個收斂序列，終止循環(huán)：t=T。

MH抽樣算法是一種簡單有效的數(shù)值模擬算法，能夠解決從未知的后驗分布中生成隨機樣本的問題。選擇合理的建議分布對抽樣的收斂效果有重要的影響。

3裂紋梁損傷參數(shù)識別

3.1基于頻率變化平方比積分的損傷指標

文獻[12]表明：當某一單元出現(xiàn)損傷或各損傷單元損傷深度相同時，結構損傷前后兩階頻率變化的平方比

（dωi）2（dωj）2只與損傷位置有關，與損傷深度無關，故可以將其定義為損傷位置指標。不考慮結構損傷引起的質量變化時，（dωi）2ωi2只與損傷深度有關，故可以將其定義為損傷深度指標。

由于簡支梁為對稱結構，如果選擇上述損傷指標，則在對稱位置上損傷指標的值相同，識別結果不唯一，導致?lián)p傷參數(shù)無法識別。采用把原損傷指標對損傷參數(shù)進行積分處理的方法，保證了損傷指標與損傷參數(shù)的對應是唯一的，這樣可以根據(jù)損傷指標有效地識別出損傷位置參數(shù)和損傷深度參數(shù)。損傷位置指標為

DL=

∫β0（dωi）2（dωj）2dβ，β為損傷位置參數(shù)，DL表征裂紋損傷位置參數(shù)的識別參量；損傷深度指標為DS=∫α0（dωi）2ωi2dα，α為損傷深度參數(shù)，DS

表征裂紋損傷深度參數(shù)的識別參量。

3.2先驗分布

對于裂紋損傷

位置參數(shù)β和損傷深度參數(shù)α的先驗概率分布，其選擇依據(jù)主要是經(jīng)驗和相關文獻資料。α和β在工程實踐中為無信息先驗分布，本文按照貝葉斯假設確定其先驗分布。貝葉斯假設中的無信息先驗分布指當參數(shù)θ的取值區(qū)間Θ和參數(shù)的含義已知、參數(shù)θ的其他信息均未知時，貝葉斯假設把取值區(qū)間Θ上的均勻分布視為θ的先驗分布。按照此假設，損傷位置參數(shù)β和損傷深度參數(shù)α的先驗分布為

3.3似然函數(shù)

由于工程中存在許多不確定性因素，例如溫度的變化、力振幅的變化、動力測試噪聲等，固有頻率的測試值與真實值之間總存在著誤差，測試值與理論計算值的誤差符合線性模型：

4算例分析

4.1損傷指標的建立

設有一單裂紋簡支梁模型，模型的幾何參數(shù)和材料常數(shù)為：跨度L=1 m，截面寬度b=0.05 m，截面高度h=0.005 m，彈性模量E=207 GPa，密度ρ=7 800 kg/m3，泊松比ν=0.3。根據(jù)式（1）計算得出損傷位置參數(shù)、損傷深度參數(shù)與損傷前后兩階頻率變化平方比的關系，見圖3。圖4為去除中間較大數(shù)據(jù)后損傷指標與損傷參數(shù)的關系圖。

從圖3和4可以看出，損傷前后兩階頻率變化平方比在損傷位置參數(shù)0.5處呈現(xiàn)對稱性，這與簡支梁屬于對稱結構有關，對稱結構的對稱位置處發(fā)生相同的損傷將引起相同的頻率變化量。為解決對稱結構引起的損傷參數(shù)難以識別的問題，以損傷位置參數(shù)為自變量，對損傷前后兩階頻率變化平方比進行積分處理，這樣可以實現(xiàn)兩者之間保持單調性，能有效識別對稱結構的損傷情況。

圖5給出了損傷深度參數(shù)α分別為0.1、0.3和

0.5時，損傷指標與損傷位置參數(shù)的關系。從圖5可以看出，不同的損傷深度參數(shù)，其損傷指標與損傷位置參數(shù)的曲線在中間點呈現(xiàn)對稱且均在水平軸上方。圖6是損傷指標的積分值與損傷參數(shù)的關系圖，兩個變量之間保持單調性，因此對損傷指標進行積分處理可以實現(xiàn)對損傷參數(shù)的有效識別。

4.2樣本值獲取

用計算機進行仿真，測試值以生成正態(tài)分布隨

機數(shù)的方式產(chǎn)生，理論值作為均值，標準差設為0.1。設梁上有一裂紋，損傷位置參數(shù)β=0.4，損傷深

度參數(shù)α=0.35，對應的第一階和第二階固有頻率分別為73.558 rad/s和295.121 rad/s，各取10組樣本值，樣本數(shù)據(jù)見表1。

4.3基于MH算法的抽樣過程

根據(jù)MH抽樣算法，選取建議分布為正態(tài)分布，

β～N（β（t），σβ2），α～N（α（t），σα2），其中β和α為候選參數(shù)，β（t）和α（t）為上一次的迭代值，σβ2和σα2

為方差，這里選?。é姚?σα2）=（0.050.1）。在利用MCMC方法進行抽樣前，必須對目標參數(shù)賦予初始值。選取不同的初始值，對后驗均值最終的收斂沒有任何影響，因此初始值的選取是任意的，可以根據(jù)相關的經(jīng)驗確定。隨著迭代次數(shù)的增加，樣本值最終會收斂于目標分布。

當損傷位置參數(shù)β和損傷深度參數(shù)α目標值分別為0.40和0.35時，設置其初始值分別為β=0.10，0.70，0.90和α=0.10，0.60，0.90，迭代次數(shù)N的最大值為2 000，通過仿真得到β和α的識別結果（馬爾科夫鏈），見圖7。

從圖7可以看出，損傷位置參數(shù)和損傷深度參數(shù)取不同初始值時其識別結果均能達到預設目標值，預設初始值與目標值的距離將決定燃燒期的長度。損傷位置參數(shù)和損傷深度參數(shù)的馬爾科夫鏈均趨于預設目標值，且當?shù)螖?shù)增加時波動較小，因此兩個參數(shù)的馬爾科夫鏈均是收斂的，達到了損傷識別的效果。

4.4損傷參數(shù)的識別結果

為降低所選初始值的影響，先去除前10%的后驗樣本，然后對后驗樣本進行頻數(shù)統(tǒng)計，所得損傷位

置參數(shù)和損傷深度參數(shù)的后驗分布頻數(shù)直方圖見圖8。損傷位置參數(shù)集中在（0.39，0.41）區(qū)間內，其中0.40對應的頻數(shù)最高；損傷深度參數(shù)集中在（0.34，0.36）區(qū)間內，其中0.35對應的頻數(shù)最高。兩個參數(shù)的后驗分布頻數(shù)直方圖在它們取不同的初始值時都保持“瘦高”形狀，類似正態(tài)分布，這與似然函數(shù)選取為正態(tài)分布的情形基本吻合。

5結論

本文用傳遞矩陣法推導了單裂紋簡支梁頻率特征方程，針對對稱結構損傷參數(shù)難以識別的問題，采用對裂紋梁損傷前后頻率變化平方比的基礎損傷指標進行積分處理的方法，建立新的損傷指標。應用貝葉斯結構損傷診斷理論構建了損傷位置參數(shù)和損傷深度參數(shù)的后驗函數(shù)。采用馬爾科夫鏈蒙特卡洛（MCMC）方法解決了后驗分布中存在高維積分難的問題，實現(xiàn)了對損傷參數(shù)的識別與估計。通過對簡支裂紋梁的損傷估算和分析，驗證了該方法的有效性，從而為梁類工程結構的損傷診斷和健康監(jiān)測提供了參考。

參考文獻：

[1]

MUTO M M. Application of stochastic simulation methods to system identification[D]. California： California Institute of Technology， 2007.

[2]BECK J L， KATAFYGIOTIS L S. Updating models and their uncertainties. I： Bayesianian statistical[J]. Journal of Engineering Mechanics， 1998， 124（4）： 455461.

[3]KATAFYGIOTIS L S， BECK J L. Updating models and their uncertainties. II： model identifiability[J]. Journal of Engineering Mechanics， 1998， 124（4）： 463467.

[4]VANIK M W， BECK J L， AU S K. Bayesianian probabilistic approach to structural health monitoring[J]. Journal of Engineering Mechanics， 2000， 126（7）： 738745.

[5]易偉建， 周云， 李浩. 基于貝葉斯統(tǒng)計推斷的框架結構損傷診斷研究[J]. 工程力學， 2009， 26（5）： 121129.

[6]李小華， 公茂盛， 謝禮立. 基于多分辨率分析的結構物理參數(shù)識別貝葉斯估計方法： 方法推導與驗證[J]. 工程力學， 2011， 28（1）： 1218.

[7]房長宇， 張耀庭， 馬超. 基于貝葉斯理論評估軸力作用下帶損傷段梁[J]. 科學技術與工程， 2014， 14（8）： 6569.

[8]芮筱亭. 多體系統(tǒng)傳遞矩陣法及其應用[M]. 北京： 科學出版社， 2008： 2395.

[9]楊秋偉， 劉濟科. 結構損傷識別的附加質量方法[J]. 工程力學， 2009， 26（5）： 159163.

[10]RIZOS P F， ASPRAGATHOS N， DIMAROGONAS A D. Identification of crack location and magnitude in a cantilever beam from the vibration modes[J]. Journal of Sound & Vibration， 1990， 138（3）： 381388.

[11]康崇祿. 蒙特卡羅方法理論和應用[M]. 北京： 科學出版社， 2015： 86149.

[12]SALAWU O S. Detection of structural damage through changes in frequency： a review[J]. Engineering Structures， 1997， 19（9）： 718723.

（編輯賈裙平）