李莉 劉鑄永 洪嘉振

（上海交通大學工程力學系，上海 200240）

引言

浮動坐標系方法是柔性多體系統(tǒng)動力學最常用的方法，其中柔性體的運動被分成兩部分：用浮動坐標系描述的大范圍剛體運動和相對于浮動坐標系的小的彈性變形.有限單元法被廣泛的采用來描述柔性體的彈性變形，然而有限元節(jié)點坐標數(shù)目龐大，將會給動力學方程求解帶來巨大的計算負擔.如何降低柔性體的自由度，是當前柔性多體系統(tǒng)動力學研究的一個重要命題.因此，為了提高柔性多體動力學仿真的計算效率，便于控制設(shè)計和實施，就必須要對柔性多體系統(tǒng)動力學模型降階進行研究.降階后的模型既要能真實地反映出系統(tǒng)的動力學特性，階數(shù)也要足夠低.

模型降階一般可以從兩方面予以考慮［1］：一方面是從建模的角度進行降階，即根據(jù)經(jīng)典的假設(shè)模態(tài)法，選擇具有較好展開收斂性的模態(tài)集，以便用較少的模態(tài)對系統(tǒng)進行建模.由于假設(shè)模態(tài)方法沒有考慮作用在物體上的外部載荷的空間特性，其收斂速度很慢.為了提高收斂性，人們提出模態(tài)綜合法等方法如來改進假設(shè)模態(tài)方法.但是，如何選取一組恰當?shù)挠商卣髂B(tài)和假設(shè)模態(tài)構(gòu)成的基矢量集是一件非常困難的事情，因為這取決于工程人員的經(jīng)驗和對此問題的認識程度.模型降階的另一方面是在模型建立之后，對所建立的模型通過提供合適的降階方法，進行合理的自由度減縮，從而進一步對模型降階.從結(jié)構(gòu)動力學或者有限元理論角度出發(fā)，目前結(jié)構(gòu)模型的降階方法主要有：頻率截斷法、慣性完備性準則和參數(shù)匹配法等傳統(tǒng)降階方法.

這些方法僅考慮了模型本身的屬性，而忽略了控制系統(tǒng)設(shè)計要求以及外界干擾等因素.在現(xiàn)代控制理論中，從系統(tǒng)的可控性、可觀性出發(fā)，利用一些降階方法將那些不可觀、不可控的模態(tài)剔除.Skelto［2］提出的模態(tài)價值分析方法和 Moore［3］提出來的平衡準則等，綜合考慮了系統(tǒng)本身屬性、系統(tǒng)的可觀可控性以及外界干擾等因素，是目前比較全面實用的降階方法.國內(nèi)學者曲廣吉等［1，9］等對柔性航天器動力學模型降階問題進行了深入的研究，取得了非常良好的結(jié)果.章敏、蔡國平［4］研究表明對與頻率密集的柔性板，內(nèi)平衡降階方法比模態(tài)價值分析方法更為有效.在狀態(tài)空間系統(tǒng)，近年來提出的Krylov子空間模型降階方法，能夠有效的處理大規(guī)模動力系統(tǒng)的模型降階問題［5］.Eberhard［6，7］等對Krylov子空間模型降階技術(shù)進行改進，對柔性體的動力學模型降價進行了很好的研究.總之，當前對柔性結(jié)構(gòu)的模型降階研究較多，對柔性多體系統(tǒng)動力學的研究相對較少.本文采用Krylov模型降階方法［8］，對柔性梁剛?cè)狁詈蟿恿W模型降階問題進行研究.首先采用有限單元法離散彈性變形，建立系統(tǒng)的剛?cè)狁詈蟿恿W方程.然后將柔性梁作為一個輸入輸出系統(tǒng)，采用Krylov模型降階方法縮減系統(tǒng)自由度，建立系統(tǒng)降階的剛?cè)狁詈蟿恿W方程.再次采用傳統(tǒng)的模態(tài)降階方法，縮減系統(tǒng)自由度，建立系統(tǒng)降階的剛?cè)狁詈蟿恿W方程.最后，分別采用有限元全模型、Krylov降階模型和模態(tài)降階模型進行系統(tǒng)的剛?cè)狁詈蟿恿W仿真.

1 剛?cè)狁詈蟿恿W的模型降階

做大范圍運動的任意一個柔性體的動力學方程可以表示為

其中q表示大范圍剛體運動廣義坐標，p表示有限元節(jié)點坐標，Mr，Me，Mer分別代表剛性質(zhì)量陣、彈性質(zhì)量陣、剛?cè)狁詈腺|(zhì)量陣，Ke代表彈性剛度陣，Qr，Qe分別為剛體運動和彈性運動的廣義力陣.

為了提高柔性多體動力學仿真的計算效率，便于控制設(shè)計和實施，就必須要對柔性多體系統(tǒng)動力學模型降階進行研究.如果有限元節(jié)點坐標p（t）∈RN，能夠用一個低維的坐標來近似.他們滿足如下關(guān)系式

其中轉(zhuǎn)換矩陣 V∈RN×s，s?N.把式（2）帶入式（1）降階的動力學方程可以表示為

降階的轉(zhuǎn)換矩陣V∈RN×s可以是模態(tài)矩陣Φ∈RN×s，也可以通過其它方法如Krylov模型降階方法得到.

模型降階過程中暫不考慮大范圍剛體運動的影響，僅對柔性體的彈性變形坐標進行降階.柔性體彈性振動的結(jié)構(gòu)動力學方程為

其中 Be∈RN×m為輸入矩陣，Ce∈Rn×N為輸出矩陣，u∈Rm×1為載荷向量，y∈Rn×1為觀測量，系統(tǒng)的傳遞函數(shù)表示為

不失一般性，下文僅討論Be＝CTe的工況.引入轉(zhuǎn)換矩陣V∈RN×s，降階后系統(tǒng)的結(jié)構(gòu)動力學方程為

降階后系統(tǒng)的傳遞函數(shù)表示為

在相同的輸入條件下，降階前后系統(tǒng)的觀測量相同，即降階前后的傳遞函數(shù)滿足以下關(guān)系式

假設(shè)Ke為非奇異矩陣，傳遞函數(shù)H（s）在σ＝0處泰勒展開為

其中mi∈Rn×m為第i次矩，Krylov模態(tài)降階方法又稱為矩匹配法，即降階后系統(tǒng)傳遞函數(shù)的第 i次矩滿足以下關(guān)系式

轉(zhuǎn)換矩陣的二階Krylov形式為［8］

其中數(shù)學符號 κr（A1，A2，G）定義為

H（s）在σ≠0處泰勒展開與 H（s+σ）在0處的泰勒展開相等

系統(tǒng)轉(zhuǎn)換矩陣形式為

其中

2 剛?cè)狁詈蟿恿W方程

本節(jié)建立了典型的剛?cè)狁詈蟿恿W模型，由一個中心剛體—柔性梁組成.以此為例，說明如何實現(xiàn)柔性多體動力學的模型降階過程.基于以下基本假設(shè)：梁是均勻的，各向同性的，細長的，剪切影響可以忽略；柔性梁的變形為小變形；中心剛體大范圍轉(zhuǎn)動的角速度較低，即零次近似耦合模型即可適用.

2.1 運動學描述

如圖1所示，中心剛體的半徑為r0.柔性梁在固定在剛體的O1點，梁長度為 L.系統(tǒng)在XY平面內(nèi)轉(zhuǎn)動，重力影響忽略.為了描述系統(tǒng)的運動，建立了兩個坐標系：慣性坐標系和浮動坐標系是點O1的位置矢量，是未變形梁上任意一點P0的位置矢量，是變形矢量，梁上任意一點P的位置矢量在坐標系上的位置坐標為

其中 r0＝［r0，0］T，ρ0＝［x，y］T，u＝［u1，u2］T是在坐標系上的位置坐標，A是從到的轉(zhuǎn)換矩陣，給出如下：

其中θ是梁大范圍運動的轉(zhuǎn)角.

對（16）式求導可得

其中的變分可以表示為

對（16）式求t的二次導數(shù)可得

圖1 動力學模型Fig.1 Dynamic Model

2.2 變形描述

采用非笛卡爾變形坐標描述，柔性梁的變形場可用線性梁單元插值函數(shù)描述：

此即為單元變形在浮動坐標系上的插值公式，其中pk為單元節(jié)點坐標，為平面梁的插值函數(shù).使用布爾定位陣BK，（21）式可以表示成全部節(jié)點坐標表示的形式.

其中 p（t）∈RN為梁的整體節(jié)點坐標，為整體函數(shù).

2.3 動力學方程

采用Hamilton原理來建立作大范圍運動平面柔性梁的動力學方程：

其中，T為系統(tǒng)的動能由中心剛體動能TH和柔性梁動能TB組成，U為系統(tǒng)的勢能，W為外力所作的功.

中心剛體和動能有關(guān)的變分項可表示為

柔性梁和動能有關(guān)的變分項可表示為

柔性梁變形能的變分為

假設(shè)作用在柔性梁上的外力為分布參數(shù)形式的外力（例如重力），其在連體坐標系下的坐標陣為在慣性坐標系下的坐標陣為f＝Af′.中心剛體受外力矩為Me.那么外力所做的虛功可表示為

將 （24）、（25）、（26）和（27）式代入（23）式，可以得到作大范圍運動中心剛體-柔性梁的剛?cè)狁詈蟿恿W方程為

其中 θ，p廣義坐標，Mθ，Meθ，Me是分塊質(zhì)量陣，Ke是彈剛度陣，Qθ，Qe是廣義外力陣.篇幅關(guān)系它們的具體形式不再給出.

降階后，式（28）可表示為

通過Krylov降階時W＝V，通過模態(tài)降階時W＝Φ.

3 仿真算例

帶末端質(zhì)量的中心剛體-柔性梁系統(tǒng)的幾何和材料參數(shù)給出如下：質(zhì)量密度 ρ＝2.767×103kg／m3，梁長 L＝1.8m，截面慣性矩 I＝1.302×10-10m4，截面積 A＝2.500×10-4m2，彈性模量 E＝6.895×1010N／m2，中心剛體的轉(zhuǎn)動慣量 JH＝0.25kg·m2，中心剛體的半徑r0＝0.1m；梁末端有一外力F＝［40.2］T作用；初始時刻梁未變形.

圖2 w0＝1rad／s時系統(tǒng)動力學響應(yīng)Fig.2 Response of the system when w0＝1rad／s

圖3 w0＝2rad／s時系統(tǒng)動力學響應(yīng)Fig.3 Response of the system when w0＝2rad／s

圖2，圖3給出了中心剛體初始角速度w0＝1rad／s以及w0＝2rad／s時分別采用有限元全模型、Krylov降階模型和模態(tài)降階模型，進行系統(tǒng)剛?cè)狁詈蟿恿W仿真結(jié)果.數(shù)值仿真結(jié)果表明，與采用模態(tài)降階方法（需要25階模態(tài)）相比，采用Krylov模型降階方法（需要5階模態(tài)）只需要較低的自由度，就可以得到與采用有限元方法吻合較好的結(jié)果.

4 結(jié)論

浮動坐標系方法是柔性多體系統(tǒng)動力學最常用的方法，其中柔性體的運動被分成兩部分：用浮動坐標系描述的大范圍剛體運動和相對于浮動坐標系的小的彈性變形，因此又稱之為剛?cè)狁詈蟿恿W.有限單元法被廣泛的采用來描述柔性體的彈性變形，然而有限元節(jié)點坐標數(shù)目龐大，將會給動力學方程求解帶來巨大的計算負擔.如何降低柔性體的自由度，是當前柔性多體系統(tǒng)動力學研究的一個重要命題.本文以中心剛體-柔性梁系統(tǒng)為例，分別采用Krylov方法和模態(tài)方法進行降價.然后分別采用有限元全模型、Krylov降階模型和模態(tài)降階模型，對中心剛體-柔性梁進行剛-柔耦合動力學仿真.仿真結(jié)果表明，與采用模態(tài)降階方法相比，采用Krylov模型降階方法只需要較低的自由度，就可以得到與采用有限元方法完全一致的結(jié)果.說明Krylov模型降階方法能夠有效的用于柔性多體系統(tǒng)的模型降價研究.

1 繆炳祺，曲廣吉，夏邃勤等.關(guān)于柔性航天器動力學模型降階問題.中國工程科學技術(shù)，2001，3（11）：60～64（Miao BQ，Qu G J，Xia SQ，et al.On the order reduction of dynamicsmodels of flexible spacecrafts.Engineering Science，2001，11（3）：60～64（in Chinese））

2 Skelton R E，Gregory C Z.Measurement feedback and model reduction bymodal costanalysis.In：Proceedings of the Joint Automatic Control Conference.USA：Denver，1979：211～218

3 Moore B C.Principal component analysis in linear system：controllability，observability and model reduction.IEEE Transaction on Automatic Control，1981，26（1）：17～31

4 章敏，蔡國平.密集頻率柔性板的內(nèi)平衡模型降階及其主動控制研究.工程力學，2009，26（11）：161～167（Zhang M，CaiG P.Balanced reduction and controlof flexible plate with dense frequencies.Engineering Mechanics，2009，26（11）：161～167（in Chinese））

5 Bai Z.Krylov subspace techniques for reduced-order modeling of large-scale dynamical systems.Applied Numerical Mathematics，2002，43（1）：9～44

6 Michael L，Eberhard P.A two-step approach formodel reduction in flexible multibody dynamics.Multibody System Dynamics，2007，17：157～176

7 Fehr J，Eberhard P.Simulation process of flexible multibody systems with non-modal model order reduction techniques.Multibody System Dynamics，2011，25（3）：313～334.

8 Salimbahrami B，Lohmann B.Order reduction of large scale second-order systems using Krylov subspacemethods.Linear Algebra and its Applications，2006，415（2）：385～405

9 次永偉，邱大蘆，付樂平等.航天器振動試驗控制技術(shù)進展.動力學與控制學報，2014，03：193～200（Ci YW，Qiu D L，F(xiàn)u L P，et al.Progress in spacecraft vibration testing control technology.Journal of Dynamics and Control，2014，03：193～200（in Chinese） ）