洪嘉振 劉鑄永

(上海交通大學工程力學系，上海 200240)

變拓撲柔性多體系統(tǒng)接觸碰撞動力學研究*

洪嘉振 劉鑄永?

(上海交通大學工程力學系，上海 200240)

在實際工程領域中存在著大量接觸碰撞等非連續(xù)動力學問題，現(xiàn)有的解決柔性多體系統(tǒng)連續(xù)動力學過程的建模理論與方法，已經無法解決或無法很好解決這些問題．本文基于變拓撲思想，提出了附加接觸約束的柔性多體系統(tǒng)碰撞動力學建模理論;通過設計柔性圓柱桿接觸碰撞實驗，驗證了所提出附加約束接觸碰撞模型的有效性;針對柔性多體系統(tǒng)全局動力學仿真面臨時間和空間的多尺度問題，提出多變量的離散方法，從而提高了柔性多體系統(tǒng)非連續(xù)動力學的仿真效率．

柔性多體系統(tǒng)， 接觸碰撞， 變拓撲， 數(shù)值仿真， 實驗研究

引言

近三十年來國內外學者在多體系統(tǒng)連續(xù)過程建模理論、計算方法和實驗研究等方面做了大量的工作，取得了一系列重要科研成果，相繼解決了動力剛化、微分代數(shù)方程數(shù)值求解、最優(yōu)控制等重要問題，建立了程式化的多體系統(tǒng)剛柔耦合動力學建模理論和有效的數(shù)值求解算法，并形成了諸如ADAMS、DADS和RECDYN這樣的大型工程軟件，這些軟件已經在眾多工程問題的實踐中經受了考驗，大大提高了新產品的設計開發(fā)的效率．

然而在現(xiàn)代航空航天、機械制造等領域，高速、輕質的系統(tǒng)面臨的動力學環(huán)境與過程越來越復雜．柔性多體系統(tǒng)動力學中，現(xiàn)有的解決連續(xù)動力學過程的建模理論與方法，已經無法解決或無法很好解決這些問題［1－6］．

例如，目前我國載人航天與未來空間站的航天器交會對接;現(xiàn)代航天器已經成為大型空間的機械系統(tǒng)，各類衛(wèi)星太陽電池陣的展開，到位鎖定時的接觸碰撞;空間大型柔性部件(機械臂、太陽電池陣、大型天線陣等)的非連續(xù)驅動;鉸鏈間隙引起的撞擊等等．在這些非連續(xù)動力學過程中，系統(tǒng)的拓撲會發(fā)生突變，以其作用時間短、撞擊力大給系統(tǒng)造成劇烈擾動，引起系統(tǒng)動力學性態(tài)的突變，對航天器機構動力學性態(tài)的災難性影響已經引起了關注．

柔性多體系統(tǒng)非連續(xù)動力學過程的實質是若干不同拓撲系統(tǒng)連續(xù)過程的不斷切換．也就是具有某個自由度的某種連接關系的一種機構的運動，突然改變?yōu)榫哂胁煌杂啥鹊牧硪环N機構的運動．當前最典型的非連續(xù)問題是接觸碰撞．多體系統(tǒng)接觸碰撞問題在力學本質上是一種非光滑、變邊界的高瞬態(tài)、高非線性動力學過程．對這類問題的研究既要考慮柔性體變形運動與大范圍運動之間的耦合，又要考慮碰撞動力學過程的合理建模，其中對碰撞過程的正確處理是解決多體系統(tǒng)碰撞動力學問題的關鍵．

現(xiàn)有的接觸碰撞模型都還有不足之處．基于剛體假設的沖量動量方法雖然具有計算效率高、直觀等優(yōu)點，且對碰撞過程中出現(xiàn)的能量不協(xié)調等問題進行了很好的解決，仍然難以適應以上工程領域對碰撞力及系統(tǒng)碰撞動力學響應求解提出的需求．基于彈簧阻尼的連續(xù)碰撞力公式的出現(xiàn)在一定程度上克服了碰撞力難以求解的困難，并在工程領域中得到了相當廣泛的應用，長期以來也為兩個問題所困擾:碰撞力模型合理形式的選擇和碰撞力方程中參數(shù)的合理辨識．這兩個問題直接決定了仿真結果的可信性．目前存在的實驗研究成果和方案往往只針對具體的個案，缺乏普適性．非線性有限元法、無網格等方法的優(yōu)點是可以對碰撞動力學過程進行建模，但柔性體的大范圍空間運動及彈性變形不加區(qū)分地采用慣性坐標系下的節(jié)點坐標表示，由此得到的動力學方程的廣義坐標數(shù)目非常龐大，計算效率低，難以應用到復雜的柔性系統(tǒng)動力學分析中．同時由于該方法無法區(qū)分剛體運動和彈性變形，給控制系統(tǒng)的設計帶來難以逾越的困難．因此，迫切需要提出一種柔性多體系統(tǒng)碰撞動力學過程全局仿真方法，以滿足柔性多體系統(tǒng)碰撞動力學對通用性、高效率、準確性和可靠性的要求．

本文基于非連續(xù)動力學過程通過合理切換可分解為若干連續(xù)過程的思路，將復雜機械系統(tǒng)變拓撲過程統(tǒng)一在通用的柔性多體系統(tǒng)模型下，研究通用與高效的柔性多體系統(tǒng)剛柔耦合動力學模型與接觸碰撞動力學模型，解決變拓撲柔性動力學全局仿真系統(tǒng)切換的關鍵問題．采用附加約束的方法能夠真正從物理上對接觸碰撞問題進行準確描述，從而實現(xiàn)柔性多體動力學全局仿真．通過數(shù)值仿真和實驗研究驗證本文所提出的附加約束方法的有效性．針對變拓撲柔性多體系統(tǒng)動力學的非連續(xù)、非光滑特性和存在時間及空間的多尺度問題，提出多變量方法，研究了柔性多體系統(tǒng)非連續(xù)動力學全局仿真計算方法．

1 接觸碰撞動力學模型

1．1 附加約束方法

基于彈簧阻尼的多體系統(tǒng)連續(xù)碰撞力模型目前在工程領域得到廣泛應用．然而該方法將碰撞力作為物體間嵌入量和嵌入速度的函數(shù)來計算，且需要人工設置函數(shù)中的參數(shù)．仿真結果過分依賴于人工參數(shù)的選擇，況且嵌入量本身也不符合物理實際［1－2］．

圖1 不同碰撞階段碰撞Fig．1 Different Stages of impact

如圖1所示，本文將柔性多體碰撞分為碰撞前、碰撞階段、碰撞后三個過程，采用附加接觸約束描述碰撞階段的動力學過程，將上述三個過程的動力學方程形式統(tǒng)一．通過引入碰撞檢測機制，實現(xiàn)拓撲模型的自動切換．將此方法稱為附加接觸約束方法．

柔性多體系統(tǒng)碰撞前和碰撞后動力學方程為［3］:

這里M是柔性多體系統(tǒng)質量矩陣，Φq是多體系統(tǒng)鉸約束的雅克比矩陣，λ是與鉸約束對應的拉格朗日乘子矢量，w，fo和fu分別是柔性體系統(tǒng)受到的慣性力矢量，廣義外力矢量和廣義變形力矢量，γ是加速度約束方程的右端項．

碰撞階段，接觸位置位移相同，故附加接觸位移約束方程為:

將其寫成加速度約束方程的形式，即

結合式(1)，柔性多體系統(tǒng)碰撞階段動力學方程為:

和γC是碰撞過程中附加接觸約束方程的左端項和右端項矩陣．

1．2 與LS－DYNA罰函數(shù)法計算結果比較

如圖2所示，水平面上柔性桿B1與柔性桿B2的中心處發(fā)生碰撞，B1和B2的軸線相互垂直．B1為自由桿，B2兩端受固定端約束．碰撞發(fā)生時B1的速度為0．2m/s，B2處于靜止狀態(tài)．桿的幾何和材料參數(shù)見表1，B1和B2扁平等厚，故該問題可簡化為平面問題處理．對此算例分別用ACCM法和LS－DYNA罰函數(shù)法進行仿真．使用計算機的配置參數(shù)為 AMD Athlon(tm)64 processor 3200+，2．01GHz，1GB 的內存．

表1 桿的幾何和材料參數(shù)Table 1 Geometrical and material data

圖2 柔性桿碰撞示意圖Fig．2 Collision between flexible bars

在利用LS－DYNA罰函數(shù)方法仿真時，需要人工設定罰因子．對同一個罰因子，計算步長不同，計算結果也不同．隨著計算步長的減小，計算結果將趨于一個穩(wěn)定解．對于本算例，表2列出了利用LS－DYNA罰函數(shù)方法仿真時，對于不同的罰因子獲得穩(wěn)定解所需要的最大時間步長，以及相應的CPU仿真時間．利用本文提出的方法，克服了這種對罰因子的依賴．

表2 各方法獲穩(wěn)定解所需時間步長與CPU計算時間Table 2 Time step required for stable solution and CPU computing time

此外，從圖3－6的罰函數(shù)方法計算結果(虛線)發(fā)現(xiàn)，利用罰函數(shù)方法，對于不同的罰因子(依次分別對應罰因子為 0．01、0．1、1．0 和 5．0 時的情況)，所得到的穩(wěn)定解是不同的．過小的罰因子(圖3)會造成仿真結果嚴重失真．當罰因子增大到一定程度，仿真結果將變化不大(圖5－圖6)，可以認為仿真結果收斂于“精確解”．比較計算結果，本算例能得到“精確解”的罰因子應大于1．0．理論上罰因子，得到的解為精確解．但實際計算中，罰因子取值越大，為了得到穩(wěn)定結果，則計算步長要取得非常小，嚴重影響了計算效率．可見，為了兼顧計算精度和計算效率，罰因子取值不能無限制大，必須要經過多次試算來選擇合適的罰因子．

從圖3－圖6的ACM方法計算結果(實線)發(fā)現(xiàn)，這種計算“精確解”依賴于罰因子的不足，對于本文提出的方法不存在．

圖3 ACCM與DYNA罰因子為0．01的碰撞力曲線Fig．3 Contact force(ACCM and LS － DYNA fs=0．01)

圖4 ACCM與DYNA罰因子為0．1的碰撞力曲線Fig．4 Contact force(ACCM and LS － DYNA fs=0．1)

圖5 ACCM與DYNA罰因子為1．0的碰撞力曲線Fig．5 Contact force(ACCM and LS － DYNA fs=1．0)

圖6 ACCM與DYNA罰因子為5．0的碰撞力曲線Fig．6 Contact force(ACCM and LS － DYNA fs=5．0)

2 接觸碰撞試驗研究

2．1 接觸碰撞試驗裝置

為了驗證本文提出的碰撞動力學過程建模方法的正確性，設計了兩圓柱鋁桿碰撞實驗．如圖7所示為實驗裝置示意圖．該實驗裝置包含以下部件:兩圓柱截面鋁桿，Polytech Gmbh公司制造的OFV－300F型激光測振儀(LDV)、動態(tài)應變儀、1．5V直流電源和金屬實驗臺架．碰撞過程中桿2的速度由激光測振儀進行測量，它由激光掃描頭、控制器、連接箱和工作站組成．

圖7 實驗裝置示意圖Fig．7 Experiment equipments

如圖7所示，一開始將桿B1從限制板上某一確定位置處釋放．然后桿B1自由下擺，當下擺到水平位置處時，B1桿達到最大速度，即碰撞前初始速度，碰撞前B1的速度通過激光測振儀進行測量．經過25次測量，得到碰撞時刻桿 B1的速度為0．21m/s，其中速度平均誤差小于0．5%．在碰撞實驗中碰撞力的時間歷程是一個關注點．在一些實驗中，直接通過壓電傳感器測量碰撞力變化歷程然而在實驗過程中壓電傳感器必須固定在物體碰撞表面上，這樣原有的碰撞界面就會受到損壞，導致壓電傳感器測量得到的碰撞力并非原來的碰撞力．因此，只能通過速度或者應變間接地驗證出仿真結果中碰撞力的正確性．

在本部分，將采用不用的建模方法對兩桿碰撞動力學過程進行數(shù)值仿真，其中碰撞過程中桿的幾何和運動學參數(shù)如表3所示．

表3 柔性桿的材料和運動學參數(shù)Table 3 Geometrical and material data

2．2 試驗結果與仿真結果比較

如圖8－11所示，可以看出，本文提出的附加約束方法仿真結果與實驗結果較為吻合．而彈簧阻尼模型仿真結果與實驗結果之間存在著很大的誤差．如圖10所示，為ACM方法仿真結果與實驗結果在桿B2上不同點處速度曲線比較，其中桿B2上11個測量點沿桿軸向方向均勻分布．由于圖中存在22條曲線，限于空間，沒有給出每條曲線的圖標．從曲線的比較可以看出，幾乎所有仿真結果與實驗結果都吻合地很好，這意味著由ACM方法仿真得到的碰撞過程與真實的碰撞過程相吻合．

圖8 B2端部速度響應曲線Fig．8 Velocity response of B2top

圖9 B2上距接觸端上2cm處P點應變Fig．9 Strain of Point P in B2(2cm from top)

圖10 桿B2上11個點的速度Fig．10 Velocity response of 11 points in B2

圖11 碰撞力曲線Fig．11 Contact force

從圖11的仿真結果上看，可以看出在碰撞過程中不同位置處的速度并不同時改變，這說明剛體碰撞假定中各點速度同時突然變化的假定與實際情況之間存在著一定得差異．根據如圖7所示實驗結果，同一速度條件下100cm處和0cm處的時間差為1．9561e－004s，通過實驗數(shù)據可以得到彈性波在鋁桿表面的傳播速度為5．1122e+003m/s，而一維情況下彈性波的理論傳播速度是5．0918e+003，相對誤差是0．4006%．這說明碰撞引起的速度擾動在彈性介質中的確以彈性波的形式進行傳播，同時驗證了碰撞動力學方程的速度跳躍條件的解釋符合真實情況．如圖8所示為彈簧阻尼模型和ACM方法仿真得到碰撞力時間歷程比較．由于碰撞力難以通過實驗的方法直接精確測量，因此只能通過仿真計算得到碰撞力的時間歷程．從圖8－11中ACM仿真結果和實驗結果的對比上看，可以得出由ACM方法仿真得到碰撞過程與實際的物理碰撞過程一致，這樣就可以得出由ACM方法仿真得到的碰撞力與真實碰撞力一致．此外，從圖中可以看出這樣一個復雜的碰撞力變化歷程很難通過簡單的彈簧阻尼模型進行建模．

3 多變量仿真算法

3．1 多變量離散方法

為提高全局仿真的計算效率，提出了多變量離散方法．如圖12所示，兩多體系統(tǒng)1和2發(fā)生接觸/碰撞．其中，系統(tǒng)1包含m個物體，系統(tǒng)2包含n個物體，發(fā)生碰撞的物體為系統(tǒng)1中的第i個物體Bi與系統(tǒng)2中的第j個物體Bj．對發(fā)生碰撞的物體Bi、Bj進行分區(qū)，接觸局部區(qū)域采用非線性有限元方法，其他區(qū)域用浮動參考系方法．對于接觸區(qū)域，以系統(tǒng)1上的接觸域為從物體，系統(tǒng)2上的接觸域為主物體．

圖12 多變量離散方法Fig．12 Multi－ variable discrete method

3．2 多變量方法的仿真結果

圖13所示鉛垂面內雙擺由B1、B2兩個物體組成．B1為剛體，長為0．5m，質量為5kg，對質心的轉動慣量為0．15kgm2．B2為平面彈性體，表4給出了B2的材料參數(shù)，同時給出了雙擺的初始位置和速度．雙擺自水平位置已一定初速度下落，撞擊剛性地面后又反彈．其中，表面與地面靜摩擦系數(shù)為0．7，動摩擦系數(shù)為0．6．B2物體可能發(fā)生碰撞局部的區(qū)域采用有限元節(jié)點坐標，其它區(qū)域采用浮動基位形坐標和模態(tài)坐標，兩部分交界如圖13所示，浮動基原點取在O2點，橫軸附著于B2兩節(jié)點上．

圖13 雙擺與地面碰撞Fig．13 Impact between double － pendulum with ground

表4 雙擺的材料和運動學參數(shù)Table 4 Geometrical and material data of double pendulum

本算例中，對于碰撞過程采用顯式法求解，對于未接觸和接觸分離后的階段采用隱式法求解．隱式法時間步長取h=5E－3s，未接觸階段和接觸分離后階段分別只用到了19和30個時間步．顯式法時間步長取0．2E－6s，碰撞階段用到了3000個時間步．本例約有120個自由度(B1有3個自由度，B2浮動基3個自由度+20階模態(tài)+47有限元節(jié)點)，三個階段的仿真計算均在半分鐘至一分鐘內完成(CPU:AMD Athlon(tm)64 Processor 3200+，2．01GHz)．說明多變量方法具有較高的計算效率．

圖14 碰撞力曲線Fig．14 Contact force

4 結論

本文基于非連續(xù)動力學過程通過合理切換可分解為若干連續(xù)過程的思路，將復雜機械系統(tǒng)變拓撲過程統(tǒng)一在通用的柔性多體系統(tǒng)模型下，研究通用與高效的柔性多體系統(tǒng)剛柔耦合動力學模型與接觸碰撞動力學模型，解決變拓撲柔性動力學全局仿真系統(tǒng)切換的關鍵問題．基于變拓撲思想，本文提出了附加接觸約束的柔性多體系統(tǒng)碰撞動力學建模理論;通過設計柔性圓柱桿接觸碰撞實驗，驗證了所提出附加約束接觸碰撞模型的有效性;針對柔性多體系統(tǒng)全局動力學仿真面臨時間和空間的多尺度問題，提出多變量的離散方法，從而提高了柔性多體系統(tǒng)非連續(xù)動力學的仿真效率．

1 Diolaiti N，Melchiorri C，Stramigioli S．Contact impedance estimation for robotic systems．IEEE Transactions on Robotics，2005，21(5):925～935

2 Zhang Y，Sharf I．Experimental validation of nonlinear compliant contact force models．Proceedings of the ASME 2007 International Design Engineering Technical Conferences ＆ Computers and Information in Engineering Conference，September 4–7，2007，Las Vegas，Nevada，USA

3 洪嘉振．計算多體系統(tǒng)動力學．北京:高等教育出版社，1999(Hong Jiazhen．Computational dynamics of multibody systems．Beijing:Higher Education Press，1999(in Chinese))

4 王禮立．應力波基礎．北京:國防工業(yè)出版社，2005(Wang L L．Foundation of stress waves．Beijing:National Defence Industry Press，2005(in Chinese))

5 Thomas Klisch．Contact mechanics in multibody system．Multibody System Dynamics，1998，2(4):335～354

6 黃筑平．連續(xù)介質力學基礎．北京:高等教育出版社，2003(Huang Z P．Fundamentals of continuum mechanics．Beijing:Higher Education Press，2003(in Chinese))

7 董富祥．剛柔耦合多體系統(tǒng)碰撞動力學建模理論與實驗研究［博士學位論文］．上海:上海交通大學，2010(Dong Fuxiang．Study of modeling theory and experiment for impact dynamics of rigid－flexible coupling multibody systems．［PhD Thesis］．Shanghai:Shanghai Jiao Tong U－niversity，2010(in Chinese))

*The project supported by the National Natural Science Foundation of China(11132007，11202126)，F(xiàn)oundation of Doctoral Program of Ministry of Education of China(20090073110009)and Natural Science Foundation of Shanghai(11ZR1417000)

? Corresponding author E－mail:zhuyongliu@sjtu．edu．cn

STUDY ON CONTACT/IMPACT DYNAMICS OF FLEXIBLE MULTIBODY SYSTEM WITH TOPOLOGY VARIABLE*

Hong Jiazhen Liu Zhuyong?
(Department of Engineering Mechanics，Shanghai Jiao Tong University，Shanghai200240，China)

There is an increasing number of discontinuous dynamics problems，such as contact and impact，in the engineering field．The known modelling theories and methods，which are used to solve the continuous dynamics of flexible multibody systems，are not able to solve those problems or to solve those problems perfectly．In this paper，an appending constraint method was proposed to solve the contact/impact dynamics of flexible multibody system based on the idea of topology variable．The impact experiment of flexible bars was designed to verify the validity of the proposed method．At last，a multi－variable method was proposed to solve the multi－scale problems of time and space when the global simulation of the discontinuous dynamics of flexible multibody systems was carried out，and the computational efficiency of the global simulation of the discontinuous dynamics of flexible multibody systems was improved by using this method．

flexible multibody sytem， contact/impact， topology variable， numerical simulation， experimental research

11 April 2012，

15 June 2012．

10．6052/1672－6553－2013－002

2012－04－11 收到第 1 稿，2012－06－15 收到修改稿．

*國家自然科學基金重點項目(11132007，11202126)，博士點基金(20090073110009)、上海市自然科學基金(11ZR1417000)資助

E－mail:zhuyongliu@sjtu．edu．cn