范新亮， 王 彤， 張麗君， 夏遵平

(1. 南京航空航天大學(xué) 機(jī)械結(jié)構(gòu)力學(xué)及控制國(guó)家重點(diǎn)實(shí)驗(yàn)室， 南京 210016；2. 中車唐山機(jī)車車輛有限公司， 唐山 063000)

有限元模型修正首先在航空航天領(lǐng)域提出，發(fā)展至今日已形成了一個(gè)龐大的理論體系，并且在運(yùn)載火箭、衛(wèi)星、航天飛機(jī)、飛機(jī)和直升機(jī)結(jié)構(gòu)的響應(yīng)與載荷預(yù)示、顫振分析、振動(dòng)控制均得到了廣泛的應(yīng)用[1]。模型修正方法按用以修正的試驗(yàn)數(shù)據(jù)不同可分為基于模態(tài)參數(shù)的方法[2-6]和基于頻響函數(shù)的方法[7-9]。Modak等[10-11]通過數(shù)值仿真詳細(xì)對(duì)比了兩種方法的差異及收斂性?；谀B(tài)參數(shù)的方法易受模態(tài)識(shí)別的影響而引入了誤差及不確定性；而基于頻響函數(shù)的方法則避免了提取模態(tài)參數(shù)帶來的誤差，且具有大量數(shù)據(jù)，對(duì)于解決待修正參數(shù)眾多的問題具有優(yōu)勢(shì)。然而，在大噪聲干擾下和初始頻響殘差較大時(shí)收斂困難、規(guī)模龐大的有限元模型修正效率低下等問題均制約著基于頻響函數(shù)的模型修正方法的進(jìn)一步發(fā)展。

為了提高頻域修正方法的收斂性，Pascual等[12-13]基于頻響相關(guān)性準(zhǔn)則，提出一種移頻方法，并利用分析和測(cè)試頻響相匹配的頻率點(diǎn)對(duì)應(yīng)數(shù)據(jù)進(jìn)行修正，提高了參數(shù)的收斂速度及范圍。Gang等[14]在該方法基礎(chǔ)上引入偽自由度的概念解決測(cè)試頻響信息不足而導(dǎo)致收斂困難的問題。為了解決自由度數(shù)目龐大的有限元模型修正問題，常常需要引入模型縮聚方法將有限元自由度縮聚至測(cè)試自由度。而利用模型縮聚難以避免引入較大的近似誤差，Li等[15]通過迭代過程更新縮聚矩陣來減小縮聚引起的誤差。模型縮聚方法僅可利用測(cè)試自由度上的頻響信息的特點(diǎn)不利于收斂，Gang等提出將有限元模型縮聚到實(shí)測(cè)和用戶選擇的偽自由度上，從而增加了用于修正的頻響的自由度數(shù)目。但是利用模型縮聚的坐標(biāo)縮減方法引起的誤差依舊較大，且形成縮聚矩陣本身計(jì)算量也較大。因此又有學(xué)者研究了采用子結(jié)構(gòu)綜合技術(shù)對(duì)系統(tǒng)的自由度縮減后進(jìn)行模型修正的方法。Papadimitriou等[16]在貝葉斯模型修正技術(shù)的基礎(chǔ)上，利用模態(tài)綜合技術(shù)得到縮減的系統(tǒng)，并計(jì)算相應(yīng)的靈敏度，減小了修正的計(jì)算規(guī)模。王陶等[17]提出了基于改進(jìn)自由界面子結(jié)構(gòu)模態(tài)綜合法的結(jié)構(gòu)模型修正方法，通過子結(jié)構(gòu)綜合技術(shù)得到由攝動(dòng)量表示的綜合系統(tǒng)方程和靈敏度方程，僅需對(duì)綜合方程進(jìn)行修正。然而對(duì)于頻響函數(shù)在綜合坐標(biāo)下的靈敏度以及如何對(duì)頻域的系統(tǒng)縮減后進(jìn)行模型修正的方法尚少有研究。

為此，本文提出一種基于縮減基的有限元模型頻響修正方法，將系統(tǒng)的位移表示為模態(tài)綜合縮減基或模態(tài)振型縮減基的線性疊加，從而將自然坐標(biāo)映射至縮減坐標(biāo)下，得到了縮減坐標(biāo)下的頻響函數(shù)殘差的靈敏度矩陣，使得修正過程均在縮減坐標(biāo)上進(jìn)行，相比利用模型縮聚將總體自由度映射至測(cè)試自由度上的方法，在計(jì)算效率及精度上均有所提高。同時(shí)引入移頻方法以增加收斂性，由原先的基于測(cè)試頻響的頻率點(diǎn)處的頻響殘差最小化改進(jìn)為基于所匹配的頻率點(diǎn)處的頻響殘差最小化；并基于模態(tài)展開式推導(dǎo)了精度較高的分段頻響擴(kuò)充方法，減小擴(kuò)充誤差對(duì)修正過程的影響。仿真算例和試驗(yàn)驗(yàn)證均表明該方法對(duì)于規(guī)模巨大的有限元模型能高效地得到穩(wěn)健的修正結(jié)果。

1 理論背景

1.1 模態(tài)綜合方法

圖1 薄板模型子結(jié)構(gòu)劃分

(1)

根據(jù)界面協(xié)調(diào)條件可得到系統(tǒng)自然坐標(biāo)下的總體位移向量u為

(2)

其中Bs為子結(jié)構(gòu)的約束模態(tài)集及主模態(tài)集所組成的Φ經(jīng)坐標(biāo)變換矩陣T、S作用后的綜合基

Bs=TΦS

(3)

(4)

式中，M、K、C、f分別為質(zhì)量矩陣、剛度矩陣、阻尼矩陣、載荷向量。式(2)代入式(4)推導(dǎo)得到系統(tǒng)在模態(tài)綜合坐標(biāo)下的運(yùn)動(dòng)方程為

(5)

1.2 基于模態(tài)綜合縮減基的修正方法

假設(shè)系統(tǒng)的真實(shí)參數(shù)為p*，對(duì)應(yīng)系統(tǒng)矩陣分別為MT、KT、CT，由式(5)在綜合坐標(biāo)下有頻域公式

(6)

類似地，假設(shè)有限元模型的初始參數(shù)為p，有：

(7)

(8)

(9)

取DT=DF+ΔD代入式(8)，其中ΔD為動(dòng)剛度矩陣的修正量可得到

(10)

與文獻(xiàn)[14]類似，可推導(dǎo)得到進(jìn)行頻率點(diǎn)匹配后的迭代修正方程

(11)

(12)

(13)

(14)

假設(shè)有限元模型按材料參數(shù)分為Np個(gè)組，則可以將ΔD(ω)表示為

γeKe)

其中：Me、Ke為單元質(zhì)量、剛度矩陣，αe、γe為單元質(zhì)量、剛度矩陣修正系數(shù)增量；Ce、De為單元阻尼矩陣，βe為單元阻尼矩陣修正系數(shù)。將該式代入式(14)可得到系統(tǒng)參數(shù)識(shí)別方程

(15)

其中θ為αe、βe、γe組成的待識(shí)別的修正參數(shù)增量，系數(shù)矩陣S為

(16)

(17)

式(15)即為利用模態(tài)綜合縮減基轉(zhuǎn)換至模態(tài)綜合坐標(biāo)下的頻域模型修正公式。

1.3 基于模態(tài)振型縮減基的修正方法

由頻響模態(tài)展開式

(18)

式中：Φ為系統(tǒng)的模態(tài)振型;Λ為系統(tǒng)特征值;E為單位矩陣;下標(biāo)l與h分別代表系統(tǒng)的低階模態(tài)與高階模態(tài)。將式(18)右乘ej并展開得到

(19)

由式(19)知任意一列頻響可通過系統(tǒng)的模態(tài)振型疊加而得。即測(cè)試頻響可以用實(shí)際系統(tǒng)的完備低階模態(tài)作為基底進(jìn)行表示。而實(shí)際系統(tǒng)的完備低階模態(tài)無法得到，因此以有限元模型的完備低階模態(tài)進(jìn)行代替。由式(19)，取擬合式

(20)

(21)

其中上標(biāo)*表示測(cè)試自由度對(duì)應(yīng)分量，由式(20)、(21)得到擴(kuò)充后測(cè)試頻響

(22)

與式(19)類似，取測(cè)試頻響及有限元頻響為

(23)

自然坐標(biāo)下基于移頻的頻響直接矩陣方法的公式為

(24)

式(23)代入式(24)簡(jiǎn)化得到

(25)

(26)

式(26)即為利用模態(tài)振型縮減基轉(zhuǎn)換至模態(tài)坐標(biāo)下的頻域模型修正公式，與式(14)形式一致，通過1.2節(jié)類似推導(dǎo)可得其迭代方程。將基于模態(tài)綜合縮減基或模態(tài)振型縮減基的方法統(tǒng)稱為縮減基方法。

1.4 算法實(shí)現(xiàn)相關(guān)問題

1.4.1 縮減基與模型縮聚的對(duì)比

應(yīng)用最廣泛的IRS縮聚(improved reduced system)方法，其縮聚矩陣為[18]

(27)

對(duì)于1.3節(jié)中轉(zhuǎn)換至模態(tài)坐標(biāo)下和通過模型縮聚轉(zhuǎn)換至測(cè)試坐標(biāo)下的修正方法，可以分別理解為以有限元模型的低階模態(tài)矩陣ΦF,l和特定的縮聚方法得到的基Td代替式(2)中的模態(tài)綜合基矩陣Bs，而后進(jìn)行與1.2節(jié)相同的推導(dǎo)得到與式(15)形式完全相同的修正公式，即

(28)

(29)

值得注意的是，三種方法無一例外需要以有限元模型的基矩陣代替相應(yīng)真實(shí)系統(tǒng)的基矩陣，帶來了近似誤差，因此三種基矩陣在近似誤差的差異上的比較也成為判斷其優(yōu)劣的標(biāo)準(zhǔn)。雖在數(shù)學(xué)上難以證明，但在實(shí)際應(yīng)用中發(fā)現(xiàn)多數(shù)情形下基于模態(tài)綜合縮減基或模態(tài)振型縮減基的方法的近似誤差要小于模型縮聚方法，尤其當(dāng)測(cè)點(diǎn)較少時(shí)，模型縮聚方法的誤差會(huì)急劇增加。以某模型為例，圖2表示分別以三種有限元基矩陣代替真實(shí)系統(tǒng)的基矩陣所形成的坐標(biāo)下動(dòng)剛度矩陣求解的頻響與真實(shí)系統(tǒng)的理論頻響的幅值相關(guān)性。可見相同情況下，基于縮減基的方法的近似誤差要遠(yuǎn)小于模型縮聚方法。

圖2 三種方法所計(jì)算頻響與理論頻響幅值相關(guān)性對(duì)比

1.4.2 基矩陣的迭代更新

subject topL≤p≤pU

(30)

根據(jù)上一小節(jié)的分析，以有限元模型的基矩陣代替真實(shí)系統(tǒng)的基矩陣將帶來近似誤差。文獻(xiàn)[15]通過迭代中不斷更新縮聚矩陣Td的方法減小了模型縮聚方法中該誤差對(duì)待修正參數(shù)識(shí)別的影響。因此本文方法中每一個(gè)迭代步以上一步修正后的有限元模型更新的縮減基矩陣來代替真實(shí)系統(tǒng)的基矩陣，隨著有限元模型參數(shù)逼近真實(shí)參數(shù)時(shí)，該近似誤差得以消除。

1.4.3 寬頻帶頻響擴(kuò)充方法

如圖3所示，通過選擇模態(tài)指示函數(shù)的半功率帶寬內(nèi)的頻率點(diǎn)可以得到測(cè)試頻響函數(shù)的分段：

圖3 測(cè)試頻響函數(shù)分段

(31)

(32)

2 算例分析

2.1 仿真算例

采用Garteur飛機(jī)仿真模型驗(yàn)證本文方法的收斂性及迭代速度。有限元模型劃分為5個(gè)區(qū)域，如圖4所示。每個(gè)區(qū)域的彈性模量、密度及阻尼系數(shù)作為待修正參數(shù)，共計(jì)15個(gè)參數(shù)。分別在無噪聲和20%噪聲水平下采用本文方法進(jìn)行修正。有限元模型劃分為3個(gè)子結(jié)構(gòu)進(jìn)行模態(tài)綜合。自然坐標(biāo)下模型的自由度為99 798，模態(tài)綜合坐標(biāo)下模型自由度僅為582，有效減小了修正過程的矩陣運(yùn)算規(guī)模。

圖4 Garteur飛機(jī)有限元模型子結(jié)構(gòu)劃分

2.1.1 無噪聲下頻響修正

將待修正參數(shù)向量的某組數(shù)值作為真實(shí)參數(shù)值，計(jì)算其有限元模型頻響作為仿真測(cè)試頻響；取上下界為真實(shí)參數(shù)值70%～130%的變化范圍對(duì)該組數(shù)值進(jìn)行攝動(dòng)得到有限元模型修正前的初始參數(shù)，計(jì)算其頻響得到初始分析頻響。初始分析頻響與仿真測(cè)試頻響對(duì)比如圖5所示，其相關(guān)性較小，殘差較大。

圖5 修正前初始分析頻響與測(cè)試頻響對(duì)比

如圖6、圖7所示為經(jīng)過3個(gè)迭代步修正后的分析頻響與測(cè)試頻響，相關(guān)性有了很好的改善。理論上在無噪聲情形下修正后的參數(shù)將與真實(shí)參數(shù)完全一致，但由于擴(kuò)充頻響過程產(chǎn)生的誤差、自由度縮減過程的假設(shè)帶來的誤差是等效的頻響噪聲，對(duì)參數(shù)識(shí)別的過程同樣產(chǎn)生了擾動(dòng)。在模態(tài)綜合坐標(biāo)下，單個(gè)迭代步運(yùn)行時(shí)間僅為220 s左右；而自然坐標(biāo)下則需要2 400 s(硬件信息：CPU-Intel i5 7300HQ，RAM-8G，OS-Windows 10 Pro，SSD-intel 256 G)?？梢娀谀B(tài)綜合縮減基的方法在保證參數(shù)收斂性的前提下，有效提高了模型修正的效率。

圖6 修正后分析頻響與測(cè)試頻響對(duì)比

圖7 修正前后分析頻響與理論測(cè)試頻響幅值相關(guān)性對(duì)比

2.1.2 20%噪聲水平下頻響修正

初始待修正參數(shù)設(shè)置同2.1.1節(jié)。如圖8所示，仿真測(cè)試頻響添加20%的噪聲，以驗(yàn)證基于模態(tài)綜合縮減基的方法在大噪聲情形下的修正精度。

圖8 修正前初始分析頻響與測(cè)試頻響對(duì)比

在大噪聲情形下，經(jīng)過10個(gè)迭代步修正后達(dá)到收斂，除第6階模態(tài)對(duì)應(yīng)頻帶外，其余頻帶的修正后分析頻響與測(cè)試頻響相關(guān)性有了很好的改善，如圖9、圖10所示?？梢娫诖笤肼暩蓴_下，基于模態(tài)綜合縮減基的方法收斂性仍較好。若采用基于模態(tài)振型縮減基的方法，取80階有限元模態(tài)作為坐標(biāo)變換基矩陣，模態(tài)坐標(biāo)下模型自由度為80，單步迭代所需時(shí)間僅為180秒，且修正精度與基于模態(tài)綜合縮減基的方法類似。

圖9 修正后分析頻響與測(cè)試頻響對(duì)比

圖10 修正前后分析頻響與理論測(cè)試頻響幅值相關(guān)性對(duì)比

2.2 試 驗(yàn)

以某型列車轉(zhuǎn)向架構(gòu)架作為試驗(yàn)對(duì)象，進(jìn)行多參考點(diǎn)錘擊試驗(yàn)，采用N-Modal模態(tài)參數(shù)識(shí)別軟件識(shí)別試驗(yàn)結(jié)果。從測(cè)試角度，測(cè)點(diǎn)布置應(yīng)清晰地反映出各階振動(dòng)形態(tài)，參考點(diǎn)布置需要考慮到避開各階振型的節(jié)點(diǎn)；從模型修正角度，布置盡可能完備的測(cè)點(diǎn)不僅能有效提高擴(kuò)充精度，且能有效反映有限元模型與測(cè)試模型的振動(dòng)形態(tài)差異，所修正的有限元模型亦愈逼近真實(shí)模型。如圖11所示。

圖11 試驗(yàn)結(jié)構(gòu)與測(cè)試設(shè)備

有限元模型劃分為10個(gè)區(qū)域，如圖12所示。每個(gè)區(qū)域的彈性模量、密度及阻尼系數(shù)作為待修正參數(shù)，共計(jì)30個(gè)參數(shù)。由于構(gòu)架為焊接的一體化結(jié)構(gòu)，其各個(gè)區(qū)域之間焊縫厚度的不確定性以及焊接區(qū)域的力學(xué)性質(zhì)的未知導(dǎo)致了有限元模型與實(shí)際結(jié)構(gòu)相差較大。同時(shí)，數(shù)目眾多的待修正參數(shù)均給模型修正帶來了較大的難度。設(shè)置初始待修正參數(shù)后計(jì)算得到初始分析頻響，其中某些測(cè)試自由度的初始分析頻響與測(cè)試頻響對(duì)比如圖13所示，其頻響相關(guān)性較小。

圖12 構(gòu)架有限元模型子結(jié)構(gòu)劃分

圖13 修正前初始分析頻響與測(cè)試頻響對(duì)比

有限元模型劃分為3個(gè)子結(jié)構(gòu)進(jìn)行模態(tài)綜合。自然坐標(biāo)下模型的自由度為181 356，模態(tài)綜合坐標(biāo)下模型自由度僅為1 830，大大減小了修正過程的矩陣運(yùn)算規(guī)模。如圖14、15所示，采用基于模態(tài)綜合縮減基的方法進(jìn)行修正，通過若干迭代步后，修正后分析頻響與測(cè)試頻響相關(guān)性有了很好的改善。

檢驗(yàn)修正前后的模態(tài)頻率變化，如表1所示，修正前后平均頻率誤差由12.72%降低至0.61%。因此，基于模態(tài)綜合縮減基的方法不僅在初始分析頻響與測(cè)試頻響相關(guān)性較差、待修正參數(shù)較多、具有較大測(cè)試噪聲干擾等情況下收斂性較好，且有效提高了模型修正的效率。若采用基于模態(tài)振型縮減基的方法，取80階有限元模態(tài)作為坐標(biāo)變換基矩陣，模態(tài)坐標(biāo)下模型自由度為80，修正精度與基于模態(tài)綜合縮減基的方法類似。

圖14 修正后分析頻響與測(cè)試頻響對(duì)比

圖15 修正前后分析頻響與測(cè)試頻響幅值相關(guān)性對(duì)比

表1 修正后各階頻率對(duì)比

3 結(jié) 論

針對(duì)規(guī)模巨大的有限元模型，本文提出了基于縮減基的頻域有限元模型修正方法，通過數(shù)值算例與試驗(yàn)結(jié)構(gòu)的模型修正表明：

(1) 對(duì)于自由度高達(dá)幾十萬的有限元模型，利用模態(tài)綜合縮減基或模態(tài)振型縮減基可以有效減小矩陣運(yùn)算規(guī)模，提高修正算法的計(jì)算效率。

(2) 在測(cè)試噪聲較大、初始分析頻響與測(cè)試頻響殘差較大、待修正參數(shù)較多等復(fù)雜情形下，本文方法仍具有較好的收斂性。

(3) 轉(zhuǎn)向架模型修正中，本文方法修正后的分析頻響與測(cè)試頻響的相關(guān)性有較大的提高，這對(duì)實(shí)際結(jié)構(gòu)在較復(fù)雜情形下的大型有限元模型修正具有一定的意義。