基于模型縮聚-頻響函數型模型修正的子結構損傷識別方法

方有亮, 婁佳琪, 張 穎*, 李宗嬈, 侯童非

(1.河北大學建筑工程學院， 保定 071002； 2.河北省土木工程監測與評估技術創新中心，保定 071002)

隨著社會的發展進步，實際工程中結構的規模越來越大，結構形式越來越復雜，對大型結構的損傷識別方法的研究是必不可少的[1]。大型復雜結構在實際試驗中面臨著測試信息少、修正參數過多、計算量過大等問題。運用傳統損傷識別的方法對整個結構進行求解，會造成計算難度增加，迭代計算消耗大量時間，而且還可能造成難以收斂等問題。為了解決上述困難，可以將整體結構分成多個子結構進行損傷識別。由于子結構是相互獨立的，因此運用子結構方法可單獨識別各子結構。侯吉林等[2]對約束子結構模型修正方法的提出與研究，使得利用少量信息就可準確識別子結構損傷情況；翁順等[3]利用子結構在主模態的特征信息進行模型修正；Li等[4]提出了一種基于頻域動態響應重構的子結構損傷識別方法；周林仁等[5]提出了針對大跨斜拉橋結構的子結構模型修正方法；雷鷹等[6]對大型結構在部分觀測下的子結構損傷診斷法進行了研究。基于此，采用子結構的Guyan靜力縮聚方法，該方法無須考慮子結構模態坐標的概念，不進行子結構主模態分析，方法簡單，操作性更強。

對于頻響函數型模型修正方法而言，采用直接實測得到的頻響函數而避免了模態法進行模態參數識別的誤差，應用更加廣泛。Hu等[7]提出了不依賴于靈敏度分析的正交正則模態修正方法的思想；徐張明等[8]對一種靈敏度分析的頻響函數型模型修正方法進行研究；劉金梅等[9]利用基于靈敏度分析的模型修正方法對井架結構進行分析識別；Wang等[10]將Kriging模型引入算法，提出了一種新的基于加速度頻響函數的模型修正的方法。對于頻響函數型的模型修正方法，要采用實測頻響函數。而實際實驗中，許多自由度信息是無法得到的，如梁的轉角自由度，這就會出現有限元分析模型中的自由度數目遠大于實測的自由度數目。因此，將理論模型進行自由度縮聚，從而滿足理論自由度與實際的相匹配。張以帥[11]提出將縮聚模型代入到模型修正的方法中進行迭代修正。該方法能夠準確識別參數，同時避免了靈敏度分析帶來的誤差和計算效率。但經過多次計算發現，對于測點數目少的識別，使用上述方法，在迭代計算過程中會出現奇異。

為此，將子結構方法與模型縮聚-頻響函數型模型修正迭代方法結合，運用到損傷識別中。首先，將結構劃分出多個子結構，采用子結構Guyan靜力縮聚方法，縮聚掉不識別子結構的內部自由度，與待識別的子結構組裝形成組裝子結構;然后,再次進行自由度縮聚，僅保留測點自由度信息;最后，運用模型縮聚-頻響函數性模型修正方法識別子結構。

1 子結構的Guyan靜力縮聚

將整體結構劃分成N個子結構，按內部和邊界自由度對子結構的剛度矩陣和質量矩陣進行分區，則第i(1,2,…,N)個子結構的剛度矩陣、質量矩陣如式(1)所示：

(1)

式(1)中：上標I表示第i個子結構的內部自由度；上標B表示第i個子結構對接界面自由度。

在Guyan縮聚中選擇對接界面自由度為主自由度(保留)，內部自由度為副坐標(縮減)，則子結構動力方程為

(2)

忽略副坐標上質量對整體的慣性力，由式(2)可得:

(3)

(4)

式(4)中:Ti為Guyan縮聚的轉換矩陣,其公式為

(5)

(6)

這時，子結構在經過靜力縮聚變換后的動力方程為

(7)

將子結構的縮聚方程組集在一起，利用對接界面位移協調和力的平衡條件就可以形成系統縮聚后的動力方程：

(8)

式(8)中：MG、KG分別為子結構縮聚組集后的系統質量矩陣和剛度矩陣；uB為對接界面位移向量。

2 模型縮聚-頻響函數型模型修正的子結構損傷識別方法

模型縮聚-頻響函數型模型修正迭代方法為

(9)

由于在測點較少的情況下，使用模型縮聚-頻響函數型模型修正迭代方法，在迭代計算中較容易出現奇異，所以將子結構與該方法結合，進行損傷識別。結構的節點編號由阿拉伯數字“1,2，…”表示，單元號由“(1)，(2)，…”表示，測點編號為“①，②，…”。如圖1所示，為一個三層平面剛架，該結構共有15個節點，18個單元，節點1、11、12處為固定端。以該三層平面剛架為例，推導模型縮聚-頻響函數型模型修正的子結構損傷識別方法。

圖1 整體結構節點和單元圖Fig.1 The nodes and elements of global structure

假設在實驗中測點為節點5、6、7、8、9、10的水平方向，則只能測得保留這6點水平自由度信息的實驗頻率ωt和實驗頻響函數Ht(ωt)，已知待識別子結構為剛度損傷，進行損傷定位和損傷定量。

(1)首先將整體結構劃分成2個子結構，分別對兩個子結構的節點和單元重新編號，如圖2所示，子結構1中共有11個節點，12個單元，子結構2中有7個節點，6個單元。其次運用子結構Guyan靜力縮聚，對子結構2進行簡化，將內部自由度縮聚到邊界自由度上。子結構2的剛度矩陣、質量矩陣分別用K2、M2表示，Guyan縮聚的轉換矩陣為Ta2，則縮聚后子結構的剛度矩陣、質量矩陣為

(10)

(2)將子結構2的內部自由度縮聚到附加邊界自由度上，即縮聚到圖2(b)中節點5、6、7上。再將縮聚后的子結構2的剛度矩陣和質量矩陣附加到子結構1剛度矩陣和質量矩陣的相應位置，形成組裝子結構，如圖3所示。組裝子結構的質量矩陣、剛度矩陣分別用KK、MM表示。

圖3 組裝子結構Fig.3 Assembley substructure

(11)

設損傷的子結構1剛度矩陣為Kt1，則對于實驗模型質量矩陣不變，剛度矩陣為

(12)

(3)組裝子結構的有限元模型再次簡化，進行Guyan縮聚，僅保留圖3中節點5、6、7、8、9、10水平向自由度信息，得到最終的有限元簡化模型并進行了測點編號，如圖4所示。

圖4 簡化模型及測點編號Fig.4 Simplified model and the number of measuring points

縮聚的轉換矩陣為Ta1，縮聚后的質量矩陣、剛度矩陣為

(13)

實驗模型縮聚的轉換矩陣為Tt1，實驗模型縮聚后的剛度矩陣為

(14)

(4)由式(10)～式(14),可得：

(15)

(16)

式中：m為子結構1單元數(m=1,2,…,12)；Δαj為子結構1各單元剛度矩陣修正量；Kj為子結構1各單元剛度矩陣。

(5)將式(15)、式(16)代入式(9)，忽略阻尼得：

(17)

(6)判斷修正量是否小于一個小值ε。如果不小于，則將各單元的修正量代入式(18)、式(19)，對組裝子結構的剛度矩陣及子結構1中的各單元剛度矩陣進行修正。

(18)

(19)

子結構損傷識別流程圖，如圖5所示。

圖5 子結構損傷識別流程圖Fig.5 The flow chart of substructure damage identification

3 數值模擬

圖6為一平面剛架。剛架總高0.84 m，總寬 0.56 m，每根梁的橫截面積為0.031 m×0.008 m，梁的彈性模量E=2.1×1011Pa，質量密度ρ=7 850 kg/m3。

圖6 三層平面剛架Fig.6 Three-layer plane frame

該剛架結構共有18個單元、15個節點、45個自由度。將該剛架結構拆分成兩個子結構，如圖2所示。子結構1共有12個單元、11個節點。子結構2共有6個單元、7個節點、9個內部自由度、9個界面自由度。子結構2進行子結構Guyan靜力縮聚，得到組裝子結構，如圖3所示。最終進行自由度縮聚得到簡化模型，如圖4所示，保留節點5、6、7、8、9、10水平自由度信息。對比這三個模型的前6階固有頻率，如表1所示。

表1 各階段固有頻率對比Table 1 Comparison of natural frequencies in different stage

由表1可知，經過縮聚后，前三階頻率基本一致，誤差不到0.2%，越高階頻率則誤差越大，因此，取前三階頻響函數進行計算。頻響函數中加入1%噪聲干擾。

3.1 單一單元損傷

假設子結構1中2單元剛度損傷31%(彈性模量減少模擬剛度損傷)。

如圖7所示，在測點①處激勵，測點③得到的單損傷前后的頻響函數曲線。曲線中的峰值的頻率即為表2中對應的頻率值。由圖7、表2可以看出，單一單元損傷前后頻響函數發生了變化。將每次迭代求出的子結構1的單元剛度矩陣修正量求和，得到子結構1中每個單元的損傷程度為

圖7 單一單元損傷前后頻響函數對比Fig.7 Comparison of FRFs before and after single element damage

表2 單一單元損傷前后頻率對比Table 2 Frequency comparison before and after single element damage

Δαj=[-0.012 2 0.296 7 -0.018 1 -0.028 7

0.064 7 -0.056 2 0.019 1 -0.012 3

0.038 4 -0.023 8 -0.004 0 -0.025 0]T

(20)

單一單元損傷結果如圖8所示，能夠準確識別出損傷的位置及損傷程度，最終識別結果為子結構1中2單元剛度損傷29.67%，識別結果較準確。其他未損傷單元的識別結果最大誤差在5%左右。

圖8 單一單元損傷識別結果Fig.8 Single element damage recognition result

3.2 多單元損傷

3單元損傷19%，5單元損傷25%，11單元損傷10%。

如圖9所示，在測點①處激勵，測點③得到的多損傷前后的頻響函數曲線。曲線中的峰值的頻率即為表3中對應的頻率。可以看出，多單元損傷前后頻響函數發生了變化，前三階頻率分別降低了0.28、1.35、1.73。將每次迭代求出的子結構1的單元剛度矩陣修正量求和，得到子結構1每個單元的損傷程度：

圖9 多單元損傷前后頻響函數對比Fig.9 Comparison of FRFs before and after multi-element damage

表3 多單元損傷前后頻率對比Table 3 Frequency comparison before and after multi-element damage

Δαj=[-0.032 8 0.040 2 0.152 9 -0.026 3

0.254 8 -0.005 6 -0.017 6 0.045 6

-0.043 5 -0.030 3 0.134 6 0.01]T

(21)

多單元損傷結果如圖10所示，能夠準確識別出損傷的位置及損傷程度，最終識別結果為子結構1中3單元損傷15.29%，5單元損傷25.48%，11單元損傷13.46%。能夠較準確識別出單元損傷程度。其他未損傷單元的識別結果最大誤差在4%左右。

圖10 多單元損傷識別結果Fig.10 Multi-element damage recognition result

算例中設定的收斂精度為ε=1×10-5，單一單元損傷識別經過7次迭代后收斂至該精度。如圖11所示，損傷單元從第4次迭代后參數識別進入收斂，直至7次迭代全部單元收斂到設定精度。損傷單元最終識別的剛度損傷為每次迭代結果之和(圖11)。多單元損傷識別經過9次迭代，識別參數收斂至設定結果，如圖12所示。從表4中的數據可以明確每次迭代的參數修正量，最終得出單元剛度的損傷程度。

圖11 單一單元損傷迭代收斂情況Fig.11 Single element damage iteration convergence

圖12 多單元損傷迭代收斂情況Fig.12 Multi-elements damage iterative convergence

表4 多單元損傷損傷單元收斂過程Table 4 The convergence of damaged elements in multi-element damage

4 結論

提出了基于模型縮聚-頻響函數型模型修正的子結構損傷識別方法。并通過數值模擬算例，對該方法進行了分析和驗證，得出以下結論。

(1)子結構的Guyan靜力縮聚方法無須考慮子結構模態坐標的概念，不進行子結構主模態分析，方法簡單，操作方便。

(2)通過對三層剛架數值模擬，在含噪聲干擾下，對本文方法進行了驗證，結果表明該方法簡便，易于計算，結果準確。

(3)由于劃分了子結構，因此可以識別任意子結構，相比于識別整體結構，更加高效，方便。