基于附加虛擬質量的儲罐損傷識別方法

侯吉林， 景天雨， 王鵬飛， 張青霞， JANKOWSKI ukasz

(1. 大連理工大學 土木工程學院， 遼寧 大連 116023; 2. 大連理工大學 海岸和近海工程國家重點實驗室， 遼寧 大連 116023; 3. 大連民族大學 土木工程學院， 遼寧 大連 116600; 4. 波蘭科學院基礎技術研究所智能技術部， 華沙 波蘭 02-106)

隨著我國石油化工的發展以及國家的石油戰略儲備計劃[1]的實施，服役和興建的儲罐會越來越多，然而這些儲罐一旦發生工程事故，將會對人民生命財產造成很大的損失，因此，對現役儲罐的監測受到越來越多的關注。目前，儲罐的檢測方法主要有射線檢測法、聲發射檢測法[2]、超聲波檢測法等、漏磁檢測法[3]、渦流檢測法[4]等。這些方法的優點是不依賴于結構有限元模型，但往往不能定量的評估損傷程度。近年來,基于振動測試的損傷識別技術得到了廣泛的發展和應用[5],選擇一個合適的損傷指標非常重要[6]。結構模態是結構物理參數的特征函數，結構物理參數的變化必然導致結構模態參數的變化[7]，因此模態被廣泛應用于結構的損傷識別。許成祥等[8]基于試驗模態分析的基本原理, 完成了對立式拱頂鋼制原油儲罐模型的應變模態試驗。楚孝田[9]根據模態分析技術的原理，結合大型原油儲罐振動特性，提出用模態分析技術來對儲罐損傷進行檢測的方法。然而實際的儲罐形狀往往是圓柱型，結構環向對稱，模態密集，且同一高度的不同位置的損傷可能引起的模態改變是相同的，這些都為儲罐結構的損傷識別增加了困難。

楊秋偉等[10]在原結構某些位置附加已知的集中質量，楊志春等[11]在原結構上添加剛度系數已知的剛度元件，形成新的參考結構，然后測試參考結構的頻率信息，聯合原結構和參考結構的頻率信息識別損傷參數。周衛東等[12]通過在結構的某些部位附加已知的集中質量元件，利用廣義柔度靈敏度矩陣建立了損傷識別方程。Holnicki-Szulc[13]提出了虛擬變形法( Virtual Distortion Method，VDM)，利用虛擬變形法來模擬結構發生的損傷，可以有效的識別結構發生的損傷。Suwala等[14]通過 VDM 法在無模型條件下進行了質量識別。張青霞等[15]通過VDM法識別了車-橋耦合系統的移動質量。Hou等[16]通過VDM法在平面框架結構上附加虛擬質量，利用數值模擬驗證了方法的有效性。本文利用附加物理參數增加了測試頻率數據量的思想，來提高儲罐結構損傷識別的精度。并進一步結合虛擬變形法(VDM)的基本思想，通過附加虛擬質量代替附加實際質量，提高方法應用的靈活性?；谔摂M變形法，只需測量一組激勵和一組加速度響應，即可數值地構造出附加虛擬質量后結構的頻率響應，可避免在實際工程中附加真實質量所面臨的一些困難，降低了實際操作的難度[17]。

本文首先介紹如何利用虛擬變形法附加虛擬質量，然后基于附加虛擬質量的子結構模態信息初步判斷損傷位置，之后基于靈敏度推導子結構附加虛擬質量的損傷識別公式，最后通過儲罐模型驗證方法的有效性。

1 附加質量的儲罐損傷識別方法

1.1 附加虛擬質量

根據VDM法的基本思想，通過一組激勵和一組加速度響應可以數值地構造出附加虛擬質量m后結構的頻率響應[18]，實現在原結構上附加虛擬質量。構造的公式為

(1)

式中：Hpp(ω,m)表示結構在p自由度方向附加虛擬質量m后對應該自由度方向的加速度頻率響應；F(ω)和A(ω)分別表示原結構在p自由度方向上的激勵和加速度的頻域響應；m是p自由度方向上附加的虛擬質量。利用構造的頻率響應識別出頻率，作為損傷判斷的依據。

下面利用單自由度系統(如圖1)，對基于虛擬變形法構造虛擬質量的操作方法進行介紹。

圖1 單自由度系統

如圖1質量、剛度和阻尼為M，K和C。在荷載F(ω)的作用下，其加速度響應A(ω)為

(2)

設M=1、K=2和C=1，ω=1,F(ω)=2,則對應A(ω)=-1+i，對應的頻率響應h(ω)=-0.5+0.5i。下面分別利用式(2)和(1)計算增加質量m=0.5后的頻率響應。

(1) 直接法，若在質量增加m=0.5，則將M=1.5代入式(2),得到A1=-0.8+1.6i，對應的頻率響應h1(ω)=-0.4+0.8i。

(2) 虛擬變形法，將m=0.5，A(ω)=1，F(ω)=2代入式(1),可以得hm(ω)=-0.4+0.8i。

兩種方法得到的頻率響應相同h1(ω)=hm(ω)，驗證了方法的正確性。

1.2 損傷位置的初步定位

一般儲罐為圓柱型，相同高度位置環向壁厚相同。利用這個結構特點，通過環向比較附加質量后的頻率，可以初步判斷損傷的位置和程度。

圖2 子結構劃分示意圖

(3)

(4)

1.3 基于靈敏度的損傷識別

如果利用上述方法初步判斷儲罐存在損傷，則可以結合有限元模型進行損傷程度的識別。

將獲得的虛擬結構的模態稱為試驗模態，儲罐有限元模型上附加相同虛擬質量后對應的模態稱為理論模態。把理論模態與試驗模態之間的誤差范數作為目標函數，通過優化子結構參數使目標函數最小，實現損傷優化識別。構造目標函數T(μ)，見式(5)，為了方便表示，將上節中的按照空間位置排列的頻率矩陣，排列為一個列向量。優化方法選取牛頓法。

(5)

式中：Δωi(μ,m)為損傷因子為μ，子結構i附加虛擬質量m時理論模態ωi(μ,m)與試驗模態ωA,i(m)之間的相對誤差。損傷因子為損傷前后的剛度的比值。計算公式如式(6)

(6)

對T(μ)分別求一階導和二階導可以得到

(7)

(8)

這里，Δω(μ,m)是所有子結構試驗和理論模態的相對誤差的向量形式，P是Δω(μ,m)相對于μ的Jacobi矩陣。我們將P定義為靈敏度矩陣，P，Δω(μ,m)和μ按照式(8)構造

(9)

(10)

(11)

將式(7)和式(8)代入式(11)，整理可以得到

μ=μ0+P+Δω(μ,m)

(12)

式中：P+表示P的廣義逆，可由奇異值分解計算。損傷因子μ可以由式(12)近似求出，然后將求得的μ作為下一步迭代的初始值μ0，利用新的結構參數構造新的靈敏度矩陣P，重復上面的步驟進行迭代，直到收斂。

2 油罐模型算例

下面通過油罐模型的數值仿真，來介紹和驗證附加虛擬質量方法的有效性。

2.1 建立油罐理論有限元模型

2.1.1 模型基本信息

模型采用鋼制立式圓罐形固定頂儲罐(標準序號HG 21502.1-92-208)，直徑7 500 mm，罐壁高度7 500 mm，拱頂高度805 mm，總高度8 305 mm，罐壁厚度6 mm，底板厚6 mm，拱頂板厚5.5 mm。罐壁和底座連接處進行焊接。底板固定在地面上。罐體材料為鋼，彈性模量2.1 Gpa，密度7 850 kg/m3，泊松比0.3。

建立ANSYS有限元模型(如圖3)，罐壁、底板和罐頂均采用shell63單元。網格劃分時，采用二維映射劃分，罐壁劃分320個單元，罐底劃分192個單元，罐頂劃分192個單元，總共704個單元。罐壁每四個單元為一個子結構。

圖3 儲罐結構

2.1.2 模態分析

表1 模型的前5階頻率

第1階陣型第2階陣型

第3階陣型第4階陣型

圖4 前4階陣型

Fig.4 The first 4 mode

通過ANSYS計算得到儲罐模型的前5階頻率(見表1)和前4階陣型(見圖4)。由于儲罐結構屬于對稱的空間結構，會出現重頻率并且各頻率比較接近。為了方便結構參數的優化和運算，利用ANSYS導出剛度矩陣和質量矩陣，然后采用MATLAB進行計算。

2.2 損傷識別

本節模擬小錘激勵計算加速度響應，利用加速度響應識別對應質量下的結構頻率，然后通過靈敏度的方法進行損傷識別。

2.2.1 損傷模型

如圖5，選擇儲罐罐壁三分之一區域為待識別區域，劃分為25個子結構，子結構編號用SV,H表示，其中V和H分別代表子結構在原結構中的豎直位置和水平位置。

圖5 儲罐結構

2.2.2 激勵和響應

利用圖 6模擬小錘激勵，考慮5%的白噪聲，小錘激勵的持續的時間是8.8 ms。采樣頻率為2×104Hz。小錘作用在子結構S3,3的中間位置，計算時間為10 s，考慮5%的白噪聲。圖 6為子結構S3,3中間位置前50 ms的激勵。圖 7為子結構S3,3中間位置前0.5 s的加速度響應。

圖6 小錘激勵

圖7 加速度響應

2.2.3 構造虛擬質量

利用理論未損傷模型，在子結構S3,3的中間位置施加垂直于罐壁方向的激勵，得到對應的加速度頻率響應，其頻率響應的幅值見圖 8(Original)；在子結構S3,3的中間位置中間附加集中質量200 kg后，利用有限元模型直接計算的頻率響應幅值為圖 8(FEM)；利用圖 8(Original)的頻率響應代入式(1)可計算出附加虛擬質量200 kg后的頻率響應的幅值，見圖 8(VDM)。比較兩條曲線，兩者完全吻合，說明利用VDM的方法構造附加虛擬質量后的頻率響應是完全精確的，而其不需要集裝有限元模型，大大提高計算效率。

圖8 利用VDM方法構造頻率響應

另外，從圖8中可以看出，原儲罐在激勵作用下，多階低階模態被激勵出來，且比較密集，不利于模態識別和進一步的損傷識別。在附加虛擬質量之后，儲罐以單個模態振動為主，易于識別，且通過靈敏度分析可知該模態對損傷具有較高敏感度。

2.2.4 附加虛擬質量與靈敏度關系

首先在子結構S3,3的正中間位置附加虛擬質量m。計算前5階頻率關于子結構S3,3的相對靈敏度變化，見圖 9，相對靈敏度的計算公式見式(13)，相對靈敏度為對靈敏度關于頻率歸一化的結果，有利于各階靈敏度的比較和分析。由圖 9可以看出，第1階頻率隨質量的增加相對靈敏度一直在提高；而第2、3、4階頻率的相對靈敏度保持在0.01幾乎沒有變化，第5階頻率隨質量的增加相對靈敏度一直在降低。在其他24個子結構中間附加虛擬質量，構造虛擬結構GV,H(μ,m),i=1,2,…,25,i≠1，綜合所有25個子結構GV,H(μ,m)的靈敏度與附加虛擬質量的關系曲線，V=1的虛擬結構G1,H(μ,m)選定第17階頻率，其余各虛擬結構選定第1階頻率。

(13)

圖9 附加質量后相對靈敏度與附加質量的關系

2.3 損傷位置的確定

假設儲罐結構發生損傷(如圖10)，損傷因子見表2。分別在所有25個結構的中間施加如圖6所示的小錘激勵，計算對應激勵位置的加速度響應。利用考慮5%噪聲的激勵和響應，代入式(1)，分別計算出各子結構附加虛擬質量300 kg后的頻率響應的幅值，通過傅里葉變換并進行峰值提取，得到各個虛擬結構所選定的頻率。其中V=1時，各虛擬結構選定第17階頻率；V=2,3,4,5時，各虛擬結構選定第1階頻率。然后將各個子結構的頻率匯總到表3，計算出每一個子結構的損傷識別系數(圖11)。

表2 損傷因子

圖10 損傷示意圖

表3 頻率

圖11 損傷識別系數柱狀圖

從表3看出損傷因子μV,H越小，橫向比較對應得到的子結構的頻率越低，對應圖11的損傷識別系數DV,H越大。在儲罐中間區域9個子結構(S2,2～S2,4，S3,2～S3,4，S4,2～S4,4)損傷識別系數相對較大，說明損傷應該在此區域，其中子結構S4,4損傷最大，子結構S3,3次之。

2.4 靈敏度優化方法損傷識別

通過上節的損傷識別系數，找到了儲罐結構損傷的大致位置和損傷的大致程度，但是要得到較為精確的結果，需要結合有限元模型，識別結構的損傷因子。為了方便優化，將上一節找到的損傷區域的9個子結構重新編號，編號為1～9。如圖12。

上一節我們已經計算出損傷區域每個子結構的附加300 kg虛擬質量后所選定的頻率。利用這個結果，通過式(13)計算每個子結構的第一階頻率關于9個子結構的相對靈敏度，并照式(9)組成相對靈敏度矩陣，相對靈敏度矩陣的示意圖為圖 13。由于這個靈敏度矩陣是對角矩陣，條件數比較好，所以基本上第一次迭代(i=1)就已經比較接近實際損傷因子(actual)了，然后只需迭代2次就能收斂，并能很精確地識別出損傷因子，見圖 14。

圖12 損傷示意圖

圖13 相對靈敏度矩陣

圖14 識別的損傷因子

2.5 其他位置損傷對識別的影響

為考慮識別區域以外的損傷對識別結果的影響，設置如下工況：如圖 15，其中子結構1～9的損傷因子如圖14保持不變，為待識別區域；在此區域以外的子結構a和b也存在損傷，其中a的損傷因子為0.6，b的損傷因子為0.6。這里只在子結構1～9上進行操作，包括激勵、計算加速度響應(考慮5%噪聲)和構造虛擬質量，具體方法同2.3和2.4節。利用在這9個子結構上附加虛擬質量后的頻率識別的損傷因子如圖 16。由圖16可知，除了第6個子結構，其余的識別誤差在5%以內，滿足要求。子結構6出現較大誤判是由于子結構a存在損傷且離其太近，子結構a的損傷對斜向相鄰的子結構3和9影響就很小了；而子結構b距離識別區域較遠，其損傷對識別結果的影響可以忽略不計。

圖15 損傷示意圖

圖16 識別的損傷因子

3 結 論

針對儲罐的具有空間對稱、低階模態密集且對局部損傷不敏感的特點，提出基于附加虛擬質量的損傷識別方法，并利用儲罐有限元模型驗證了方法的有效性，得到了以下結論：

(1) 通過在儲罐上附加虛擬質量，改變模態分布特點，使其具有對局部損傷靈敏度較高的頻率信息，并提高了測試模態數據量。

(2) 基于虛擬變形法(VDM)的基本思想，利用原結構的激勵和加速度響應構造附加任意質量后的頻率響應，進一步增加了方法的靈活性和應用性。

(3) 該方法分為兩步，首先通過附加虛擬質量后的頻率可快速初步判斷損傷的位置，然后結合儲罐的有限元模型，利用靈敏度矩陣快速識別儲罐的損傷程度。