一種利用子結構綜合技術的模型修正方法

王 陶， 何 歡,2， 閆 偉，陳國平,2

(1.機械結構力學及控制國家重點實驗室, 南京 210016；2.南京航空航天大學 振動工程研究所, 南京 210016)

一種利用子結構綜合技術的模型修正方法

王 陶1， 何 歡1,2， 閆 偉1，陳國平1,2

(1.機械結構力學及控制國家重點實驗室, 南京 210016；2.南京航空航天大學 振動工程研究所, 南京 210016)

提出了基于改進自由界面子結構模態綜合法的結構模型修正方法。首先給出了一種改進的自由界面子結構模態綜合法，該方法通過構造一組與系統低階模態加權正交的向量集，有效地解決了含有剛體模態時系統剩余柔度矩陣的求解問題。然后利用攝動法對每個子結構進行攝動，并計算靈敏度，通過子結構綜合技術得到由攝動量表示的綜合系統方程和靈敏度方程，給出了基于靈敏度分析的修正計算方法。該方法僅需對綜合方程進行修正，極大的減縮了待修正問題的計算規模，提高了修正效率。最后，以典型數值算例驗證了該方法的有效性。

自由界面；子結構模態綜合；攝動法；模型修正

有限元模型修正技術是用來提高有限元模型預測準確性的有效方法。縱觀模型修正技術近幾十年的發展，主要可以分為矩陣型修正和參數型修正兩大類，由于采用參數型修正后的修正量物理意義明確，修正后的模型與實際結構關聯程度大，因此是當今主流的模型修正方法。

參數型模型修正問題通常被轉化為優化問題進行處理[1-2]，通過在優化過程中不斷迭代更新靈敏度矩陣和修正參數來完成修正過程。然而，隨著當今工業技術的發展，機械結構日漸趨于大型化、復雜化，進而導致有限元模型的規模也越發龐大。采用傳統方法進行修正，勢必需要在靈敏度分析過程中付出高昂的計算代價。為此，在迭代修正這一環節中，部分學者借鑒在優化設計與可靠性分析中倍受重視的代理模型技術，使用代理模型來代替直接的有限元計算，實踐證明該方法可有效的解決復雜模型的分析效率問題。郭勤濤等[3]利用結構動力學測試的響應結果構建了響應面模型，論證了基于響應面進行結構模型修正的可行性。SHYY等[4]對不同形式的響應面擬合精度和效率進行了研究，證實了雙層徑向基神經網絡響應面具有較好的的擬合效果。YANG等[5]用Kriging響應面和多目標遺傳算法對大展弦比飛機進行了模型修正，計算結果顯示基于響應面的模型修正方法能在保證較高精度的基礎上大幅提高計算效率。

然而基于代理模型的模型修正方法同樣也存在其局限性。代理模型在構建過程中不可避免的會引入近似誤差，而模型修正本身屬于一類不適定性問題，微小的擾動可能造成修正結果與真實結果的偏離，因此代理模型的精度將直接決定修正效果。由于代理模型精度受制因素眾多，建模過程對于研究人員的經驗依賴較高，因此至今仍沒有統一的高精度建模規范可供參考，極大的限制了該模型修正方法的推廣。

基于以上原因，很多學者開始探尋其他高效的修正方法。近些年，部分研究人員將關注目光重新落在了子結構模態綜合法[6-8]這類傳統模型降階方法上，通過將子結構模態綜合法與模型修正問題相結合來提高大型復雜結構模型的修正效率。王毅等[9]將子結構模態綜合法與模型修正方法結合對星箭模型進行了修正研究。張美艷[10]提出了一種新的求解特征靈敏度的結構動力重分析方法，并將求解特征值靈敏度的方法引入固定界面子結構模態綜合法中，推導出了固定界面子結構綜合模型特征參數的一階和高階靈敏度綜合方程。侯吉林等[11]提出了一種以局部子結構為修正對象的約束子結構修正法。通過構造出的約束子結構柔度矩陣，運用迭代優化的方法修正約束子結構模型。WENG等[12]提出了一種子結構模型修正方法，對采用子結構模態綜合得到的模型進行了靈敏度分析。PAPADIMITRIOU等[13]基于子結構模態綜合法推導出了減縮修正方程，用減縮模型進行每個修正迭代步的計算，大大提高了模型修正問題的計算效率。JENSEN等[14]將固定界面法與貝葉斯模型修正方法相結合，通過一個框架結構對所提方法進行了驗證。以上的研究多為通過固定界面模態綜合法與模型修正相結合來提高模型修正的效率。然而在實際操作中，通常采用實驗測得的數據對模型進行修正，而實驗中固定界面的邊界條件實現起來較為困難。于是開始有研究者嘗試將自由界面模態綜合法應用在模型修正問題中。文獻[15-16]針對結構的動力學響應分析及常規的數值迭代算法運算效率較低的問題，提出了在自由界面模態綜合法基礎上的模型修正方法。

針對實驗中固定界面不易實現的問題，本文在文獻[17]的基礎上提出了一種考慮剩余柔度的攝動自由界面子結構模態綜合法。這種方法首先將結構分為若干子結構，采用攝動法對子結構進行攝動，然后根據子結構之間的界面協調條件推導出了含有模態參數攝動量的綜合方程，通過減縮后的綜合方程可直接求解出特征值關于修正變量的靈敏度。由于交界面采用了自由界面的連接方式，所以可以采用實驗室中易于實現的自由邊界條件的實驗數據構造目標函數，將求解出的靈敏度進一步代入到優化算法中對機翼結構有限元模型進行了修正。因為本方法直接針對減縮模型進行修正，所以可以有效提高模型修正的效率。

1 子結構攝動方程

n自由度無阻尼離散振動系統的運動方程可以表示為：

(1)

式中M和K∈Rn×n分別為系統質量矩陣和剛度矩陣，u和f∈Rn分別為廣義位移向量和載荷向量。通過特征值分析可得到系統前l階模態向量Φl∈Rn×l。

根據攝動法的思想，當結構修正變量參數發生小范圍的變化后，結構系統的剛度、質量矩陣也會產生相應的改變，假設結構攝動后的剛度矩陣、質量矩陣表示為：

K=K0+ΔK

(2)

M=M0+ΔM

(3)

式(2)和式(3)中K0和M0是原系統的結構剛度矩陣和質量矩陣，ΔK和ΔM為是原系統結構剛度矩陣和質量矩陣的攝動量。同樣，結構的模態向量也相應改變，即Φl=Φ0l+ΔΦl。

(4)

滿足

(5)

來構造一組與系統低階模態加權正交的向量集。

(6)

并將式(6)與考慮了攝動量的結構模態向量組代入式(4)和(5)得

φ=

化簡可得

(8)

若只保留一階微量，可得

(9)

式中，Φ0l是原始結構系統前l階模態向量，ΔΦl是結構設計參數改變后系統前l階模態向量改變量。

u=Φβ

(10)

將結構模態向量代入到結構振動運動方程，可以將方程化簡為：

(11)

其中：

(12)

(13)

(14)

由于系統低階模態向量和構造的加權正交向量關于質量矩陣正交，則易證明得到：

(15)

再結合式(5)可以將式(11)化簡為：

(16)

令

(17)

(18)

以上各式中的βl表示結構系統所要保留的低階模態坐標，βh則表示結構系統需要被縮減的高階模態坐標。

式(16)中的第二式可以表示為：

(19)

將式(19)經過Laplace變換后可以得到：

(20)

取式(20)中的一階近似，然后再進行Laplace反變換后可以得到：

(21)

將式(21)代入式(10)中，化簡可以得到：

u=(Φ0l+ΔΦl)βl+Gf

(22)

式中G∈Rn×n為剩余柔度矩陣，表達式如下式：

(23)

這里在求解剩余柔度矩陣G時，避免了直接對系統矩陣K求逆的過程，很好地解決了結構在自由-自由狀態下系統矩陣K不可逆而無法計算剩余柔度矩陣的問題。

假設子結構的交界面上有m個自由度，且不考慮外力載荷的作用，將廣義物理坐標按照內部坐標與界面坐標進行分塊處理，寫為如下形式：

(24)

將式(24)中的第二式展開進行化簡，可以將界面坐標表示為：

(25)

式中，Gjj∈Rm×m是結構剩余柔度矩陣與界面自由度uj對應的部分。

2 攝動方程的綜合

假設這里以兩個子結構a與b的連接為例，根據以上理論推導，由式(25)可以得到兩個子結構的界面坐標表達式：

(26)

(27)

由子結構的位移連續性的協調條件和界面力平衡方程可以求解出界面力：

式中

(29)

綜合式(16)以及式(28)，可得系統的綜合方程為

(30)

式中

(31)

(32)

Mab=Mba=0

(33)

(34)

3 子結構綜合模型靈敏度模型修正

3.1 子結構綜合模型的靈敏度

求解子結構綜合模型特征值靈敏度之前，必須先得到的剛度矩陣和質量矩陣對修正變量參數的導數，表達式如下：

(35)

式中xj表示結構系統中第j個修正變量參數。

由式(34)可以計算子結構綜合模型特征值的靈敏度：

(36)

當求解子結構綜合模型的靈敏度時，每個子結構的修正變量參數是相互獨立的，即結構的修正變量xj對整體系統的影響只與第j個子結構的剛度矩陣、質量矩陣有關，與其他子結構的結構參數無關。所以在求解靈敏度的時候，與修正變量參數xj沒有關系的靈敏度都為零。只需計算與修正變量參數有關系的模態參數靈敏度。

3.2 子結構綜合模型修正

模型修正可以轉化為式(36)所表示的求解修正變量最優解的問題

(37)

目標函數E(x)可以通過結構特征量建立的，體現了待修正模型仿真計算值與實驗測量值之間的差異。

(38)

式中fai和fsi分別為待修正模型與實驗模型的結構特征量，wi是各特征量的權系數，N是特征量數目。

通過靈敏度分析可以求解出目標函數的梯度向量，確定優化迭代的搜索方向，找到在搜索區間的最優步長。

采用優化迭代算法對目標函數進行尋優分析，根據收斂條件判斷是否滿足修正指標要求。

4 數值算例

4.1 某機翼有限元模型及子結構劃分

如圖1所示的某機翼結構有限元模型，圖2所示為其的內部結構圖。

為了驗證本文所提出方法的有效性，本文中將整個機翼結構劃分為兩個子結構，其中子結構a(如圖3所示)包含548個殼單元，501個節點，總計3 006個自由度；子結構b(如圖4所示)包含450個單元，424個節點，總計2 544個自由度。其材料屬性如表1所示。

表1 機翼結構有限元模型材料屬性

圖1 機翼整體結構模型Fig.1 Aircraft wing structure

圖2 機翼模型內部結構Fig.2 Internal wing structure

圖3 機翼子結構a有限元模型Fig.3 Component a of wing structure

圖4 機翼子結構b有限元模型Fig.4 Component b of wing structure

Tab.2 Comparison of the eigenvalues from the presented method and the full FEM Hz

首先采用子結構模態綜合法，對機翼模型進行動特性分析，分別提取子結構前10、20和30階低階模態進行子結構綜合，然后與完全有限元模型(FEM)的計算結果進行對比。

從表2可以看出，提取子結構前20階低階模態進行子結構綜合已經可以達到較高的計算精度。為了使修正結果更加精確，后文采用提取子結構前30階低階模態的綜合模型進行模型修正。

4.2 機翼有限元模型修正

這里采用梯度下降法對模型進行修正，其算法的整體流程如圖5所示。

圖5 模型修正算法流程圖Fig 5.The flow chart of model updating

(1)構造目標函數

采用MSC.Patran/Nastran建模并將其計算的結果作為實驗結果來對減縮后的模型進行模型修正，因此后文將MSC.Patran中建立的模型稱為實驗模型。

這里通過前五階固有頻率來構造優化分析的目標函數，記機翼試驗模型的固有頻率表示為ωs，機翼減縮模型計算得到的固有頻率為ωa，則目標函數為：

(39)

式中ωai是機翼有限元模型的第i階固有頻率，ωsi試驗模型測試的第i階固有頻率。

(2)選擇修正變量

結構修正變量在優化迭代的過程中隨著迭代次數變化，結構的固有特性也會隨著結構修正變量參數發生改變。

圖6 各階模態頻率對修正變量的靈敏度Fig.6 Sensitivity of natural frequencies versus design variables

對于目標函數E(x)的變化，以各個修正變量的取值范圍為約束條件，則其各個修正變量的變化范圍為：

2.0 mm≤x1≤5.5 mm

1.5 mm≤x2≤3.0 mm

2.0 mm≤x3≤5.5 mm

(40)

(3)迭代收斂條件準則

目標函數的前后相對誤差滿足收斂容差：

(41)

(4)機翼模型修正結果分析

圖7給出了無量綱化處理的前五階固有頻率以及目標函數值隨迭代次數的變化曲線。為便于對照，采用式(41)和(42)對頻率以及目標函數值分別進行了歸一化處理。

(42)

(43)

式中ai為要進行歸一化的變量，其原本所處的區間為[a1ia2i]，bi為歸一化后的變量，其所處的區間為[-1,1]。

圖7 固有頻率與目標函數值的歸一化迭代曲線Fig.7 Normalized natural frequencies and E(x)

從圖7中可以看出，隨著迭代次數的增加目標函數逐漸收斂，而各階固有頻率值逐漸趨向于修正值。

表3列舉了所選修正變量在修正前后數值的對比，可以看出翼肋厚度的變化最為明顯。從表4中所列舉的數據可以看出，機翼結構模型修正后的固有頻率與試驗模型的前五階固有頻率誤差很小，很好地說明了本文方法的有效性。

表3 機翼結構模型修正前后修正變量值對比

如圖8所示，通過優化算法進行迭代計算，機翼結構模型的前五階固有頻率與實驗模型間的誤差都明顯減小。另外，運用本文方法對機翼結構模型進行修正的計算時間約為直接對整體機翼結構進行修正時間的45%，說明采用本文的方法可以有效提高模型修正的計算效率。

圖8 機翼結構模型修正前后的固有頻率誤差對比Fig8. Errors of natural frequencies of wing before and after updating

模態階數12345有限元模型/Hz修正前24．45118．24149．74156．34246．04修正后25．18121．81153．41160．93264．38試驗模型/Hz25．19121．56154．15161．23265．05誤差/%修正前2．942．732．863．037．17修正后3．97×10-50．210．480．190．25

5 結 論

本文提出了一種基于改進自由界面子結構模態綜合法的模型修正方法。該方法在進行模態綜合的過程中通過構造一組與保留模態關于系統矩陣加權正交的向量集，使得在求解剩余柔度矩陣時避免了直接對剛度矩陣求逆，很好地解決了含有剛體模態時剩余柔度矩陣的求解問題。從而可以采用較為容易實現的自由邊界條件的實驗數據作為修正數據。

然后采用攝動法的思想推導了考慮攝動量的模態綜合方程，并通過得到的綜合方程構造了用于模型修正的靈敏度矩陣。采用減縮后模型所求得的特征值構造目標函數，并通過優化算法進行迭代求解以達到模型修正的目的。

最后通過數值算例對本文提出的方法進行了驗證。從計算結果可以看出，采用本文提出的方法進行模型修正是有效的。并且由于采用了子結構模態綜合法對模型進行了減縮，使得采用優化算法進行迭代計算時的效率有了較大提升。

[1] CHANG C, CHANG T Y P, ZHANG Q W. Iterative constrained optimization scheme for model updating of long-span bridges[J]. Proceedings of SPIE, 2000, 3995:334-345.

[2] TEUGHELS A,ROECK G D, SUYKENS J A K. Global optimization by coupled minimizers and its application to finite element model updating[J].Computer and Structures, 2003, 81:2337-2351.

[3] 郭勤濤，張令彌，費慶國.用于確定性計算仿真的響應面法及其試驗設計研究[J]. 航空學報，2006,27(1):57-61. GUO Qintao, ZHANG Lingmi, FEI Qingguo. Response surface method and its experimental design for deterministic computer simulation[J]. Acta Aeronautica et Astronautica Sinica, 2006,27(1):57-61.

[4] SHYY W，TUCKER P K. Response surface and neural network techniques for rocket engine injector optimization[J]. Journal of Propulsion and Power，2001,17(2):391-411.

[5] YANG G W，CHEN D W. Response surface technique for static aeroelastic optimization on a high-aspect-ratio wing[J]. Journal of Aircraft，2009, 46(4):1444-1450.

[6] HURTY W C. Dynamic analysis of structural systems using component modes[J]. AIAA Journal, 1965, 3: 678-785.

[7] GOLDMAN R L. Vibration analysis by dynamic partitioning[J]. AIAA Journal, 1969, 7: 1152-1154.

[8] MACNEAL R H. A hybrid method of component mode synthesis[J]. Computers and Structures, 1971, 1(4): 581-601.

[9] 王毅,朱禮文,王明宇.大型運載火箭動力學關鍵技術及其進展綜述[J].導彈與航天運載技術, 2000,243(1):1-9. WANG Yi, ZHU Liwen, WANG Mingyu. Summary of some key technologies on dynamics of large launch vehicle[J]. Missiles and Space Vehicles, 2000,243(1):1-9.

[10] 張美艷.復雜結構的動力重分析方法研究[D].上海:復旦大學,2007.

[11] 侯吉林,歐進萍.基于局部模態的約束子結構模型修正法[J].力學學報,2009,41(5):748-756. HOU Jilin, OU Jinping. Isolated substurcture model updating based on local mode[J]. Chinese Journal of Theoretical and Applied Mechanics, 2009,41(5):748-756.

[12] WENG Shun, XIA Yong, XU Youlin, et al. Substructure based approach to finite element model updating[J]. Computers and Structures, 2011, 89: 772-782.

[13] PAPADIMITRIOU C, PAPADIOTI D C. Component mode synthesis techniques for finite element model updating[J]. Computers and Structures, 2013, 126: 15-28.

[14] JENSEN H A, MILLAS E, KUSANOVIC D. The use of component mode synthesis techniques for large finite element model updating using dynamic response data[C]//Proceedings of the 9th International Conference on Structural Dynamics.Porto: EURODYN,2014: 2965-2971.

[15] 劉洋,段忠東,周道成.基于模態綜合技術的結構有限元模型修正[J].振動、測試與診斷, 2009,29(3):287-291. LIU Yang, DUAN Zhongdong, ZHOU Daocheng. Finite element model updat ing of complex structure using component mode synthesis technique[J]. Journal of Vibration,Measurement&Diagnosis,2009,29(3): 287-291.

[16] LIU Yang, LI Yan, WANG Dejun, et al. Model updating of complex structures using the combination of component mode synthesis and Kriging predictor[J]. The Scientific World Journal, 2014(2014): 1-13.

[17] HE H, WANG T, CHEN G P, et al. A real decoupled method and free interface component mode synthesis methods for generally damped systems[J]. Journal of Sound and Vibration, 2013, 319: 335-354.

Model updating approach based on improved component mode synthesis

WANG Tao1, HE Huan1,2, YAN Wei1, CHEN Guoping1,2

(1.State Key Lab of Mechanics and Control for Mechanical Structures, Nanjing 210016, China;2. Institute of Vibration Engineering Research, Nanjing University of Aeronautics and Astronautics, Nanjing 210016, China)

A model updating technique based on improved free-interface component mode synthesis was proposed. By constructing weighted-orthogonal vector sets which have weighted-orthogonal relationship with the lower retained modes, it is easier to obtain the residual flexibility attachment matrix. Each of the components was analyzed by perturbation method. Then the modal equation and the sensitivity equation of whole system were derived by using the component mode synthesis method. The model updating method based on sensitivity analysis proposed in the paper greatly reduces the model scale and improves the computational efficiency. The application of the proposed model updating method in a wing model calculation demonstrates that the presented method is effective.

free-interface; component mode synthesis; perturbation method; model updating

國家自然科學基金資助(11472132); 中央高校基本科研業務費資助(NS2014002); 機械結構力學及控制國家重點實驗室(南京航空航天大學)自主研究課題資助(0113Y01); 江蘇高校優勢學科建設工程資助

2015-03-27 修改稿收到日期：2015-12-24

王陶 男，博士，1987年生

何歡 男，副教授，1978年生 E-mail:hehuan@nuaa.edu.cn

O321; TB123

A

10.13465/j.cnki.jvs.2017.02.024