基于有限元Craig縮聚模型的基礎激勵和載荷識別

朱斯巖

(北京航空航天大學 航空科學與工程學院，北京 100191)

朱禮文

(北京宇航系統工程研究所，北京 100076)

基于有限元Craig縮聚模型的基礎激勵和載荷識別

朱斯巖

(北京航空航天大學 航空科學與工程學院，北京 100191)

朱禮文

(北京宇航系統工程研究所，北京 100076)

載荷識別是用有限自由度的加速度測量值推斷外力.基礎激勵是在基礎自由度上加給定的加速度時間歷程激振結構，求整個結構的動響應.利用Craig減縮模型的性質，及界面自由度的已知運動可求出整個結構其余自由度的動力學響應，完成基礎激勵求解.在此基礎上用已知的界面自由運動和求得的非界面自由度的動力學響應，代入結構動力學方程即可實現載荷識別.本方法唯一能給基礎激勵和載荷識別結果帶來的偏差是Craig縮聚模型與原始有限元模型間的差異，對此用頻響分析的方式，直接從頻響曲線的變化來確定縮聚模型的精度.算例表明這一新的基礎激勵和載荷識別途徑是有效的.

結構動力學;基礎激勵;載荷識別

基礎激勵和載荷識別是工程應用中常見的二類結構動力學問題，其共同點是均用有限自由度的加速度時間歷程，實現基礎激勵或載荷識別的目的.基礎激勵類似響應計算，二者的區別在于，響應計算中動力學方程的外力函數項已知，問題轉化為二階微分方程(組)的數值積分，目前已有很多成熟精確的算法［1－3］;基礎激勵是用結構少數基礎自由度的加速度時間歷程，求結構其余自由度的響應［4－5］.文獻［4］介紹了 3 種方法:慣性釋放法，適用于基礎為1個節點6個自由度的情況;大質量法，是一種近似算法;Lagrange乘子法(LMT)，較為繁瑣.載荷識別由于種種實際因素限制，只能用少數自由度的加速度響應測量值推斷作用于結構的激振力F(t)，而對其余自由度響應所做的不同假設，都將導致載荷識別結果的不同，即通常所說的不適定性［6］，也是研究的難點.本文利用Craig減縮模型的性質，由已知的界面自由度加速度時間歷程，可精確簡捷地求得全部非界面自由度的動響應，完成基礎激勵的求解.同時，載荷識別就可以先按上述方式，用測得的界面自由度加速度時間歷程，求非界面自由度動力學響應，再將所有界面和非界面自由度的動響應代入結構動力學方程，即可實現載荷識別，獲得作用于結構界面的外載荷.

以上方法唯一會給基礎激勵和載荷識別結果帶來的偏差，是Craig縮聚模型與原始有限元模型間的差異.有限元模型的Craig縮聚是固定界面模態綜合的基礎［7］，最初它與靜力學子結構方法一樣，只是為了適應當時計算機速度低、容量小的硬件條件.1987年長征系列運載火箭開始衛星發射的國際商業服務時，在星箭耦合載荷計算中運用了這一技術.而此時構造的捆綁火箭結構動力學數學模型，同樣采用了固定界面模態綜合方法［8］.可見Craig縮聚模型已十分廣泛地應用于航空航天結構動力學分析中，因此對Craig縮聚模型與原始有限元模型間的差異進行定量分析是十分必要的.最早王文亮教授在文獻［9］中對此作了有效性分析，而普遍采用的是主模態截止頻率，認為縮聚模型與原始模型在主模態截止頻率以下的模態子空間相當.此外文獻［10］運用模態有效質量(modal effectivemass)的概念.本文提出用頻響分析的方式，直接從頻響曲線的變化來確定縮聚模型的精度.

1 Craig模型縮聚

有限元模型的Craig縮聚，要求把運動參數分為界面自由度和非界面自由度，此時運動方程變為

本文將有外力作用的自由度列為界面自由度，界面模態φc定義為

由式(2)可見，界面自由度的選取應保證矩陣Kmm可逆，因此必要時須選取無外力作用的自由度為界面自由度，并認為其上作用的外力等于0.對于其余的全部自由度均歸為非界面自由度，解如下廣義特征值問題:

得到主模態φm和主模態特征值λm.

由主模態φm和界面模態φc構造變換矩陣T:

則運動參數的物理空間和模態空間變換關系為

由式(5)可知界面自由度保留了物理空間的性質.

縮聚后式(1)變換為

其中

剛度矩陣的非對角塊，與主模態坐標及界面自由度相關的項為0，即

運動方程右側非界面自由度對應外力為0，坐標變換后對應主模態坐標的廣義力仍為0.同時與主模態對應的質量陣和剛度陣均為對角陣，已解耦，對應阻尼項可以用模態阻尼比方式引入.通常界面自由度數較少，可不計相應的阻尼項.故列入阻尼陣后式(6)變為

式(8)可分解為如下二式

模型縮聚是依靠減少主模態數實現的.若算例規模較小，允許計算出與非界面自由度數相等的主模態數，此時通過變換矩陣式(4)生成的結構模型即Craig變換模型(模型大小與原始模型相同).此時，用Craig變換模型計算的特征值、特征矢量，與原始模型計算得到的特征值、特征矢量完全一致.

界面自由度的加速度值，在基礎激勵問題中是已知的激勵項，可據此求整個結構的動響應;而在載荷識別問題中是已知的測量值，需要用它推斷作用于結構上的外力.

2 基礎激勵

基礎激勵是用加速度時間歷程激振結構基礎，因此無法像響應計算那樣，用外力時間歷程通過積分結構動力學方程獲得結構的動響應.NASTRAN早期提供了“大質量法”的基礎激勵計算模塊，近期的版本包含了直接用加速度時間歷程的計算模塊.

本文基礎激勵問題的求解運用式(9a).由于主模態坐標對應的廣義力為0，因此將已知項移至等號右側并變號，得到響應計算的方程形式如下:

3 載荷識別

載荷識別可依據Craig減縮模型，先按照第2節基礎激勵的方式，用有限的界面自由度的測量值計算主模態坐標參數的動響應，再根據式(9b)完成載荷識別計算.式中，主模態坐標加速度按基礎激勵方式求得為已知的加速度測量值;Xc為界面自由度的位移時間歷程，是未知項，可由加速度測量值數值積分二次求出，積分時初值可簡單取0，同時加速度測量值起點也定為0.初值對計算結果有影響，但因阻尼的存在將只影響計算起始的很小一段時間［11］.將代入式(9b)便能得到通過界面自由度作用于結構的載荷.

上述載荷識別的計算方法同樣適用于模擬振動臺運動的基礎激勵分析.由于振動臺作為基礎只包含1個節點6個自由度，是靜定約束.相應的界面模態是對應結構的剛體運動，因此與剛度矩陣相乘為0，故式(9b)中矩陣，可直接根據計算振動臺對試件的激振力和夾持力.顯然這對振動臺反共振問題分析具有重要價值.

4 Craig縮聚模型精度評估

頻響分析是在結構的給定自由度上，加簡諧形式的單位載荷，計算結構指定自由度的位移(速度、加速度)響應隨加載頻率變化的規律.計算得到的位移(速度、加速度)頻響曲線，是結構固有特性(傳遞特性).目前所有有限元商用軟件都有頻響分析模塊.

結構動力學方程為

其中，響應U和激勵F是簡諧函數，頻率為ω.

通常響應和激勵間存在相位差.將式(12)、

式(13)代入式(11)，得到頻響分析運動方程:

將式(14)分成實部、虛部:

通過式(14)或式(15)求解U1和U2，分別由U1，U2對應的自由度處取值，求該自由度的位移響應、相位差，再乘以ω，ω2，即得到速度、加速度響應.

用Craig縮聚模型作頻響分析時，式(14)、式(15)中的質量陣、剛度陣和阻尼陣可直接由式(8)得到，但激勵自由度必須列為界面自由度，響應自由度建議也列入界面自由度.否則必須將式(14)或式(15)中的 U1和 U2分別取出，用生成Craig縮聚模型的變換陣T，按式(5)恢復到物理空間，再由相應自由度處取值求該自由度的位移響應、相位差，再乘以ω或ω2，獲得速度響應、加速度響應.

5 算例

算例是長10m，直徑500mm，壁厚3mm的鋁質圓管.數學模型包含11個節點，共66個自由度，由10個等長度的鐵木辛可梁單元組成，不加約束.每個節點集中質量900 kg，集中轉動慣量0.3 kg·m2(繞圓管中心線轉動慣量為 0.6 kg·m2).選取 10 個界面自由度，其序號為:1，2，3，4(節點1);25，26，27(節點5);49，50，51(節點9).頻響計算以自由度序號50作為加載點，自由度2的加速度頻率響應曲線如圖1所示.縮聚模型的主模態階數分別取前8，10，24階3種情況，界面自由度數為10，因此縮聚模型總自由度數分別為18，20，34.

圖1 加速度頻率響應曲線

從圖1曲線可以看出主模態階數變化時對應曲線峰值的位置和高度都會發生變化.實際上，該算例當主模態階數大于等于18后，主模態階數的增加不再引起頻響曲線計算結果的變化.此外，對于根部固支或其他多種約束方式，結構系統的頻響曲線也大體類同.

載荷識別算例也采用了同一結構模型，并在界面自由度2，6(即節點1的Y向和節點5的Y向)加激振力，計算結構動響應.其中界面自由度的加速度時間歷程作為載荷識別的加速度測量值.同時為方便起見，兩個力函數僅在時序上顛倒了一下，其余界面自由度上外力函數皆為0.計算結果顯示識別載荷同激振力十分一致.為突出二者差別，在確定界面自由度運動(作為測量值)的響應計算中，結構模型采用主模態階數等于非界面自由度數的Craig變換模型，而載荷識別時采用主模態階數等于8的Craig縮聚模型.因圖1的頻響曲線顯示這兩個模型的差異主要表現在40Hz以上的頻段，因此激振力也應有40Hz以上的頻率成分.激振力的幅頻特性如圖2所示.圖3給出了自由度2上識別載荷與激振力的比較，二者十分一致.

圖2 原始所加載荷傅里葉幅值譜

圖3 原始所加載荷與載荷識別結果比較

6 結論

本文提出利用有限元Craig縮聚模型進行基礎激勵計算和載荷識別，并根據頻響曲線確定模型精度，算例計算表明方法可行有效.

本文基礎激勵的計算方式在預示振動試驗的分析工作中已實際應用，但載荷識別方法由于要求加力點必須同時是加速度測量點，這就要求事先在加速度測點布置中有所安排考慮，因此在現有工作中進行驗證有一定困難.此外，外力作用處直接安裝加速度測點也有可能存在難度，這時可以將測點設定在外力作用點附近、與著力點間連接剛性較好的結構上，相當于將力移動到測點.若該點處在力作用線上，則對隨后的基礎激勵和載荷識別計算沒有影響.若因外力移動出現附加彎矩，則對應的加速度測量不只是線加速度，還要求增加角加速度測量.由于目前尚無法直接測量角加速度，可以在外力作用點附近尋找剛性較大的結構，用多個線加速度測點測量結果換算該點的角加速度.總之應盡可能避免附加彎矩的出現.

References)

［1］鐘萬勰.結構動力方程的精細時程積分法［J］.大連理工大學學報，1994，34(2):131 －135 Zhong Wanxie.On precise time-integration method for structural dynamics［J］.Journal of Dalian University of Technology，1994，34(2):131－135(in Chinese)

［2］ Bathe K J.Finite element procedures in engineering alalysis［M］.Englewood Cliffs，New Jersey:Prentice-Hall Inc，1982:499－532

［3］邱吉寶，向樹紅，張正平.計算結構動力學［M］.合肥:中國科學技術大學出版社，2009:260－274 Qiu Jibao，Xiang Shuhong，Zhang Zhengping.Computational structural dynamics［M］.Heifei:University of Science and Technology of China Press，2009:260 －274(in Chinese)

［4］董益亮，郭鋼，徐宗俊，等.道路激勵作用下的汽車后橋動力響應分析［J］.汽車工程，2002，24(4):339－342 Dong Yiliang，Guo Gang，Xu Zongjun，et al.The dynamic response analysis of the automotive rear axle under road excitation［J］.Automotive Engineering，2002，24(4):339 － 342(in Chinese)

［5］任建峰，仇原鷹，段寶巖，等.一種基礎加速度沖擊的仿真新方法［J］.現代機械，2006，132(2):21 －23 Ren Jianfeng，Qiu Yuanying，Duan Baoyan，et al.A new idea to access the simulation of a shock from the base［J］.Modern Machinery，2006，132(2):21 －23(in Chinese)

［6］王彥飛.反演問題的計算方法及其應用［M］.北京:高等教育出版社，2007:1－5 Wang Yanfei.Computational methods for inverse problems and their application［M］.Beijing:Higher Education Press，2007:1－5(in Chinese)

［7］ Craig R R Jr，Chang C C.Substructure coupling for analysis and testing［R］.NASA CR-2781，1977

［8］于海昌，賈文成，衛國，等.大型捆綁火箭模態試驗/分析的相關性研究［J］.導彈與航天運載技術，1993，202(2):42 －53 Yu Haichang，JiaWencheng，WeiGuo，et al.A correlation investigation ofmodel test/analysis for heavy lift strap-on launch vehicle［J］.Missiles and Space Vehicles，1993，202(2):42 －53(in Chinese)

［9］王文亮，杜作潤.結構振動與動態子結構方法［M］.上海:復旦大學出版社，1985:229－426 Wang Wenliang，Du Zuorun.Structural vibration and dynamic substructure methods［M］.Shanghai:Fudan University Press，1985:229－426(in Chinese)

［10］ López-Díez J，Marco-Gómez V，Luengo P.Modal effective mass applied to error localization of finite element models［R］.AIAA-99-1310，1999

［11］朱斯巖，朱禮文.運載火箭動態載荷識別研究［J］.振動工程學報，2008，21(2):135 －139 Zhu Siyan，Zhu Liwen.Dynamic load identification on launch vehicle［J］.Journal of Vibratioin Engineering，2008，21(2):135－139(in Chinese)

(編 輯:李 晶)

Base excitation and load identification based on Craig reduction model

Zhu Siyan

(School of Aeronautic Science and Engineering，Beijing University of Aeronautics and Astronautics，Beijing 100191，China)

Zhu Liwen

(Beijing Institute of Astronautical Systems Engineering，Beijing 100076，China)

Load identification use measured acceleration of few degree of freedom(DOF)to determine external applied load.Base excitation give movement upon foundation freedom to excite the structure，and then obtain the whole dynamic response.Making use of characteristics of Craig reduction model，and the known motion of interface DOF，the dynamic response of all other DOF could be calculated，thus the base excitation was achieved.The load identification could be performed through the dynamic equation with these motion of interface DOF and the calculated dynamic response of non-interface DOF.The only error of this method for base excitation and load identification is the difference between the Craig reduction model and the original finite-element model.For this reason，frequency response analysis was applied to determine the precision of Craig reduction model through the variation of the frequency response curve.An example was given to demonstrate the validity of the present base excitation and load identification method.

structural dynamics;base excitation;load identification

TB 122

A

1001-5965(2011)06-0665-04

2010-03-10

朱斯巖(1973 －)，女，北京人，講師，ysz4444@buaa.edu.cn.