夏恒，崔鴻知

（長沙理工大學(xué) 土木工程學(xué)院，湖南長沙 410114）

0 引言

目前，部分既有砌體[1]屬于歷史文化保護遺產(chǎn)，其歷史人文價值與日俱增，但材料老化、荷載作用、結(jié)構(gòu)疲勞、超期服役等因素會導(dǎo)致砌體結(jié)構(gòu)出現(xiàn)退化和損壞的現(xiàn)象，由于結(jié)構(gòu)復(fù)雜、損傷位置隱蔽等情況，很難對結(jié)構(gòu)進行準確的損傷定位，假如對其損傷放任不管，隨著損傷的累積，將可能造成重大安全事故和經(jīng)濟損失，故對這些砌體結(jié)構(gòu)進行損傷識別研究具有重要的意義。

在20 世紀50 年代，我國的傳統(tǒng)檢測方法（筒壓法、扁頂法、原位單磚雙剪法、原位軸壓法、點荷法）屬于定期的人工檢測，在進行檢測時還需預(yù)知結(jié)構(gòu)損傷的大致位置，對于一些難以測量的部位，其人工檢測的危險程度也會增加。到了20 世紀80 年代，隨著科學(xué)技術(shù)的發(fā)展，基于動力學(xué)的損傷識別[2]方法逐漸成為主流，其原理是結(jié)構(gòu)發(fā)生損傷時，相關(guān)參數(shù)會發(fā)生改變，通過在結(jié)構(gòu)上布置相應(yīng)的傳感器可對結(jié)構(gòu)的健康狀況進行實時監(jiān)測。

在基于動力學(xué)的損傷定位方法上， Pandey[3]提出絕對變化量指標，根據(jù)損傷前后應(yīng)變模態(tài)差的絕對值進行損傷定位。董聰[4]根據(jù)結(jié)構(gòu)應(yīng)變的一階變分關(guān)系，推導(dǎo)出了應(yīng)變模態(tài)相較于位移振型和自振頻率更能有效、準確地進行損傷定位分析，并且在進行試驗后數(shù)據(jù)表明：損傷位置應(yīng)變模態(tài)的變化十分明顯，且變化值的大小與損傷程度成正比關(guān)系。Vanik[5]提出基于參數(shù)與貝葉斯理論相結(jié)合的在線健康監(jiān)測方法。顧培因[6]提出一種在無法獲取損傷前的模態(tài)數(shù)據(jù)下的一種直接指標法，該方法可以直接利用損傷后的應(yīng)變模態(tài)參數(shù)進行損傷定位。殷志翔[7]利用應(yīng)變模態(tài)差對索網(wǎng)結(jié)構(gòu)實現(xiàn)了損傷識別。鄒云峰[8]提出了基于應(yīng)變模態(tài)響應(yīng)重構(gòu)的損傷識別方法。吳韜[8-10]等人將動力特性用于砌體結(jié)構(gòu)的損傷識別。

研究發(fā)現(xiàn)，目前砌體結(jié)構(gòu)的損傷定位研究大部分屬于確定性分析，即損傷前的參數(shù)都有準確的數(shù)據(jù)。但在現(xiàn)實中，大部分的砌體結(jié)構(gòu)因為各種原因，很難獲取損傷前的參數(shù)。本文以湖南省長沙市湖南大學(xué)財經(jīng)學(xué)院的省級保護性建筑砌體柱為例，考慮測量誤差和模型誤差的影響，以改良的應(yīng)變模態(tài)差值公式嵌入到貝葉斯模型修正的框架中。該方法可以有效地處理損傷定位中的不確定性因素，能準確地識別結(jié)構(gòu)損傷位置，并且還能將此次獲得的后驗信息作為下次更新的先驗信息，從而實現(xiàn)長期的動態(tài)觀測。

1 貝葉斯基本理論及MCMC 方法

1.1 貝葉斯基本理論

貝葉斯統(tǒng)計推斷[11]方法是引入概率學(xué)的知識來做出判斷，它是一種承認生產(chǎn)實踐的繼續(xù)性方法，將相關(guān)先驗信息（如歷史可靠數(shù)據(jù)）參與決策。推斷基礎(chǔ)來自于貝葉斯公式，其公式如下：

式（1）中P（θi）是事件θi 發(fā)生的先驗分布，與事件D 無關(guān)，先驗信息主要來源于可靠的歷史數(shù)據(jù)和專家經(jīng)驗，假設(shè)先驗分布為正態(tài)分布N（μ，λ2），可以利用先驗信息確定μ 和λ2的值；P(D|θi)是指在事件θi 發(fā)生的基礎(chǔ)上，事件D 的條件分布；是一個全概率數(shù)，起一個正則化因子的作用；P(θi|D)是在事件D 發(fā)生的基礎(chǔ)上，事件θi 的分布，也被稱為后驗分布。一般情況下后驗概率是難以表達的，對后驗概率的計算可以采用MCMC-馬爾科夫鏈蒙特卡羅法。

1.2 MCMC-馬爾科夫鏈蒙特卡羅法

在馬爾科夫鏈蒙特卡羅法中進行抽樣是至關(guān)重要的，可以降低其樣本的相關(guān)性。不同的抽樣方法，將產(chǎn)生不同的馬爾科夫鏈蒙特卡羅算法。1987 年，Durane[13]將馬爾可夫鏈蒙特卡羅法與分子動力學(xué)相結(jié)合，提出了混合蒙特卡羅算法。混合蒙特卡羅算法結(jié)合了梅特羅波利斯-黑斯挺算法、吉布斯算法以及確定性分子動力學(xué)方法，可以用于高維復(fù)雜分布的抽樣模擬，改進了建議概率分布的產(chǎn)生機制，避免陷入局部解。

混合蒙特卡羅法的算法流程[14]如下：

（1）鏈指標k=0，初始化位置向量x0和動力向量p0。

（2）鏈指標k=k+1，從正態(tài)分布N（0，M）中隨機產(chǎn)生一個新的動力向量pt。

（3）賦予蛙跳算法初值[x(0),p(0)=(xk-1,pt)]并允許L 時間步以獲得一個建議狀態(tài)（x*，p*）=[x(t+Lε)，pt+Lε)]。

（4）計算接受概率α=min[1,exp(-△H)],其中△H=H（ x*，p*）-H（xk-1,pt）。如果接受建議狀態(tài)（xk，pk）=（x*，p*），否則（xk，pk）= （xk-1，pt）。

（5）重復(fù)步驟（2）—（4），直至到達預(yù)設(shè)的鏈長為止。

再根據(jù)后驗樣本可計算后驗分布的各階矩，進行相應(yīng)的統(tǒng)計推斷。

2 有限元模型修正方法

部分既有砌體由于年代久遠，無法準確獲取結(jié)構(gòu)損傷前的參數(shù)，只能根據(jù)現(xiàn)場測量數(shù)據(jù)和相關(guān)歷史資料初步建立無損的有限元模型，但所建立的模型與實際情況存在差異。為了修正有限元模型，本文把結(jié)構(gòu)自身物理參數(shù)，如材料密度、彈性模量、以及泊松比等作為修正參數(shù)，以實際結(jié)構(gòu)的測量值與有限元模型的動力響應(yīng)特征值的相關(guān)信息來構(gòu)建目標函數(shù)[15]，以相關(guān)信息差值作為修正依據(jù)。

式（2）中：y（r）表示為現(xiàn)實結(jié)構(gòu)與有限元模型動力響應(yīng)特征值的差值，可以是應(yīng)力、應(yīng)變、頻率、振型、阻尼、應(yīng)變模態(tài)等；ME表示為現(xiàn)實結(jié)構(gòu)測量的數(shù)據(jù) ；MA（r）表示為有限元模型的特征值，其中r 表示為修正的參數(shù)向量。有限元模型修正流程如圖1 所示。

圖1 有限元模型修正流程圖

3 基于貝葉斯理論的損傷定位方法

3.1 基于應(yīng)變模態(tài)的損傷識別原理

在基于動力學(xué)的結(jié)構(gòu)損傷定位中，用于損傷定位的參數(shù)需要滿足兩個條件：一是對局部損傷敏感；二是位置坐標的單調(diào)函數(shù)[4]。

在文獻[4]中證實了應(yīng)變模態(tài)相較于其他參數(shù)更能有效、準確解決損傷定位問題。

由動力學(xué)方程：

文獻[4]推導(dǎo)出結(jié)構(gòu)應(yīng)變的一階變分函數(shù)為：

式（4）中{Δε}為應(yīng)變的變化；[Δψ]為應(yīng)變模態(tài)的變化；[Уr]結(jié)構(gòu)自振頻率矩陣；[φ]為位移模態(tài)矩陣；[ψ]為應(yīng)變模態(tài)矩陣；[ΔУr]結(jié)構(gòu)自振頻率變化；[Δφ]位移模態(tài)的變化；其中{Δε}和[Δψ]的變化在位置上存在明確的對應(yīng)關(guān)系，是具有進行損傷識別的能力。

損傷定位的公式以應(yīng)變模態(tài)差的絕對值為依據(jù)，其公式為：

式子（5）中[ψI]和[ψD]表示為無損結(jié)構(gòu)和有損結(jié)構(gòu)的應(yīng)變模態(tài)矩陣；i 表示坐標位置；j 表示模態(tài)階數(shù)。在損傷位置處，應(yīng)變模態(tài)差Δψ會發(fā)生突變，Δψ的大小表示損傷的嚴重程度，Δψ 恒大于0。

本文綜合考慮不確定性因素帶來的影響，實際的應(yīng)變模態(tài)差的絕對值調(diào)整為：

式（6）中,f1、f2分別是測量誤差和模型誤差，假設(shè)f1、f2都服從正態(tài)分布，則其均值分別為μ1、μ2，方差分別為σ21、σ22。

3.2 似然函數(shù)與后驗分布

由改良的定位公式（6），可推出似然函數(shù)表示為：

確定似然函數(shù)和先驗分布后，后驗分布[16]則為：

式（9）中，P(D)為邊緣似然函數(shù)，其值為固定常數(shù)。其中參數(shù)θ 中的每一個θj，其邊緣概率密可以通過積分獲得為：

3.3 損傷識別流程

本文將砌體結(jié)構(gòu)識別流程分為三個主要部分：（1）先基于歷史資料初步建立有限元模型，再激勵砌體結(jié)構(gòu)獲取相關(guān)參數(shù)數(shù)據(jù)，進行動力特性分析，最后以相關(guān)信息差值對有限元模型進行修正；（2）將現(xiàn)場收集的數(shù)據(jù)作為先驗信息，再通過馬爾科夫鏈產(chǎn)生候選樣本數(shù)據(jù)，依據(jù)混合蒙特卡羅法計算接受概率，確定樣本集；（3）通過貝葉斯方法，計算應(yīng)變模態(tài)的后驗分布，判斷是否收斂。具體流程圖如圖2 所示。

圖2 損傷識別流程圖

4 砌體結(jié)構(gòu)損傷定位研究

4.1 有限元模型材料參數(shù)和邊界條件

燒結(jié)磚砌體柱來源于湖南省長沙市湖南大學(xué)財經(jīng)學(xué)院的省級保護性建筑的第一層柱體，結(jié)構(gòu)尺寸為720 mm×470 mm×25 440 mm，如圖3 所示。實體結(jié)構(gòu)多處帶有明顯裂縫，且砌塊間的砂漿風(fēng)化程度嚴重，幾乎散失了強度。損傷前的有限元模型依據(jù)現(xiàn)場測量數(shù)據(jù)、相關(guān)資料和GB 50003—2011《砌體結(jié)構(gòu)設(shè)計規(guī)范》確定燒結(jié)磚為MU15、砂漿為M7.5，砌塊密度設(shè)為1 800 kg/m3、彈性模量設(shè)為1 600（ff 為砌體抗壓強度，單位為MPa），計算為3 312 MPa、體積剪切模量為1 577 MPa、剪切模量1 440 MPa、線性膨脹系數(shù)設(shè)為5×10-6/℃、收縮率設(shè)為-0.1 mm/m、摩擦系數(shù)設(shè)為0.70、熱導(dǎo)率設(shè)為1.4 W/(m·K)、泊松比設(shè)為0.15。

圖3 燒結(jié)磚砌體柱

有限元模型采用ANSYS Workbench 建模，單塊磚的模型為240 mm×115 mm×53 mm，邊界條件設(shè)置為：上部承受荷載為均布荷載，方向垂直向下無偏心，根據(jù)GB 50009—2012《建筑結(jié)構(gòu)荷載規(guī)范》和設(shè)計資料，其均布荷載值設(shè)為20 kN/m2，底部采用 Fix supporied（固定約束），加速度設(shè)為9.8 m/s2，方向豎直向下。由于模型整齊，單元類型采用六面體SOLID186 單元，單元尺寸為50 mm，網(wǎng)格劃分采用程序自動控制。有限元模型如圖4 所示。

圖4 有限元模型及檢測點布置圖

4.2 試驗方案和數(shù)據(jù)結(jié)果

在實體結(jié)構(gòu)上取五點分別為a、b、c、d、e，這五點處于結(jié)構(gòu)正面，如圖4 所示。其中a、b、c 三處是帶有明顯裂縫的位置，d、e 兩處無裂縫，通過比較現(xiàn)實結(jié)構(gòu)和有限元模型的應(yīng)變模態(tài)差值判斷結(jié)構(gòu)損傷位置，差值越大說明損傷越嚴重。對實體結(jié)構(gòu)a、b、c、d、e 測量點分別測量3 次，以其平均值作為先驗樣本（表1），再經(jīng)過混合蒙特卡羅法計算，設(shè)置樣本接受的概率，經(jīng)過多次循環(huán)后，使樣本值向目標值收斂。由于前面的樣本值高度依賴于初始值，故舍棄前200 個樣本。本文分析前四階模態(tài)下應(yīng)變模態(tài)差值的后驗分布。

表1 實際結(jié)構(gòu)檢測點應(yīng)變模態(tài)

有限元模型在修正參數(shù)后得到無損傷前應(yīng)變模態(tài)云圖（圖5），檢測點a、b、c、d、e 的應(yīng)變模態(tài)如表2 所示，實際結(jié)構(gòu)和有限元模型檢測點的應(yīng)變模態(tài)如圖6 所示。

表2 有限元模型檢測點應(yīng)變模態(tài)

圖5 前四階應(yīng)變模態(tài)云圖

圖6 實際結(jié)構(gòu)與有限元模型檢測點應(yīng)變模態(tài)

從圖6 中可以明顯發(fā)現(xiàn)：測量點處應(yīng)變模態(tài)值的大小與結(jié)構(gòu)損傷程度成正相關(guān)。測量點（a、b、c）在一階模態(tài)和四階模態(tài)下的應(yīng)變模態(tài)值明顯高于二階模態(tài)和三階模態(tài)下的應(yīng)變模態(tài)值，其應(yīng)變模態(tài)差值在一階模態(tài)和四階模態(tài)下有著更為明顯的峰值，故在損傷定位時，針對此結(jié)構(gòu)采用一階模態(tài)和四階模態(tài)能更準確地識別損傷位置。在檢測點（d、e）處，各階模態(tài)上的應(yīng)變模態(tài)差值相較于測量點（a、b、c），應(yīng)變模態(tài)差值保持在穩(wěn)定低值范圍內(nèi)。將五點的應(yīng)變模態(tài)差值作為先驗信息帶入計算模型中，一階模態(tài)下測量點處的應(yīng)變模態(tài)差值統(tǒng)計結(jié)果如表3 所示。

表3 測量點處一階應(yīng)變模態(tài)差值

考慮測量誤差和模型誤差的影響（式6），其均值和標準差在不同的采樣點取值不同，將修正的應(yīng)變模態(tài)差值的均值帶入馬爾科夫鏈中，生成樣本集，通過貝葉斯算法確定后驗分布圖（圖7），具體數(shù)值如表4 所示。

表4 測量點處一階應(yīng)變模態(tài)差值的修正值

圖7 各測量點應(yīng)變模態(tài)差值的后驗分布圖

從圖7 和表4 中可以發(fā)現(xiàn)經(jīng)過貝葉斯推斷后，其應(yīng)變模態(tài)差值的均值和標準差有相對明顯的修正，在帶有裂縫的a、b、c點應(yīng)變模態(tài)差值與沒帶裂縫d、e 點應(yīng)變模態(tài)差值存在明顯不同，在帶有裂縫處的檢測點其應(yīng)變模態(tài)差值遠遠高于沒有裂縫處的檢測點，如圖8 所示。通過應(yīng)變模態(tài)差值的大小可以判斷出a、b、c 點相較于d、e 點是存在損傷，并且其損傷的程度與應(yīng)變模態(tài)差值成正相關(guān)。

圖8 應(yīng)變模態(tài)差值折線圖

5 結(jié)論

貝葉斯理論屬于一種概率統(tǒng)計理論，該方法可借助先驗信息（包括可靠的歷史數(shù)據(jù)以及專家經(jīng)驗），再基于總體信息和樣本信息得出合理推斷。本文將該方法與有限元方法共同運用于砌體柱的損傷定位研究，結(jié)論如下：

（1）砌體柱的損傷前有限元模型可以通過現(xiàn)場測量數(shù)據(jù)和相關(guān)資料以及經(jīng)驗進行初步的建立，利用實際結(jié)構(gòu)與有限元模型的動力特征值的差值作為修正依據(jù)，以貝葉斯算法作為優(yōu)化方法，做到信息差最小，完成有限元模型的建立。

（2）貝葉斯推斷和有限元相結(jié)合的方法在結(jié)構(gòu)損傷定位研究中具有可行性，貝葉斯理論可以有效地處理損傷定位中不確定性因素的影響，從而得出更為合理的推斷。