錢震杰， 章定國， 金誠謙

(1. 農(nóng)業(yè)部南京農(nóng)業(yè)機械化研究所，南京 210014； 2. 南京理工大學(xué) 理學(xué)院，南京 210094)

柔性多體系統(tǒng)含摩擦碰撞stick-slip過程動力學(xué)仿真

錢震杰1， 章定國2， 金誠謙1

(1. 農(nóng)業(yè)部南京農(nóng)業(yè)機械化研究所，南京 210014； 2. 南京理工大學(xué) 理學(xué)院，南京 210094)

基于高次剛?cè)狁詈侠碚摵蚅agrange 乘子法，研究了柔性多體含摩擦碰撞stick-slip過程的全局動力學(xué)的精確建模與自動切換仿真問題?；谧兺負渌枷耄鶕?jù)分離、碰撞、黏滯接觸和滑動接觸等狀態(tài)分別構(gòu)造相應(yīng)的約束條件和動力學(xué)方程。運用沖量/動量法求解碰撞初始條件；引入切向滑動摩擦力勢能的概念描述切向滑動接觸力；給出接觸、分離、黏滯、正向/逆向滑動狀態(tài)之間的切換準(zhǔn)則，實現(xiàn)了系統(tǒng)全局動力學(xué)自動切換。通過算例的數(shù)值仿真，分析了滑移/黏滯(微滑動)、正/逆向滑動等復(fù)雜非光滑現(xiàn)象，驗證了該模型和算法的有效性。

含摩擦碰撞；黏滯/滑移；接觸約束；柔性機械臂；動力學(xué)建模

柔性多體系統(tǒng)，如機械臂，操作過程存在著大量的抓取、分離、內(nèi)外碰撞，末端執(zhí)行器和操作控制較復(fù)雜。目前，工程應(yīng)用上大量采用沖量/動量法和連續(xù)碰撞力法來近似求解多體系統(tǒng)的碰撞問題。Pfeiffer等[1]提出了附加約束法求解多體系統(tǒng)變結(jié)構(gòu)問題。由于沖量動量法存在碰撞瞬時完成假設(shè)，無法求解接觸力變化過程，而基于Herze理論的連續(xù)接觸力法存在碰撞體互相嵌入假設(shè)，動力學(xué)響應(yīng)依賴于接觸模型參數(shù)選擇。附加約束法[2-5]可得到碰撞過程接觸力的變化、系統(tǒng)動力學(xué)響應(yīng)，彌補前兩者方法的不足。

一方面，該方法需要解決接觸發(fā)生瞬時兩界面間的速度跳躍的非光滑問題，并保證接觸發(fā)生瞬時和分離之后碰撞力為零，否則可能出現(xiàn)接觸和分離狀態(tài)之間接觸點速度和加速度的高頻“假振蕩”現(xiàn)象，或者數(shù)值不收斂。為此，Taylor等[6]基于Newmark積分格式的碰撞初始接觸點速度和加速度的修正形式引入方程，得到穩(wěn)定數(shù)值解。董富祥等[2]引入連續(xù)介質(zhì)力學(xué)間斷面理論來求解碰撞初始速度跳躍問題。另一方面，附加約束法等約束優(yōu)化算法普遍存在由于數(shù)值積分誤差累加造成的約束方程違約問題。目前主要有兩大思路解決這一問題，即降階轉(zhuǎn)為常微分方程法和引入違約修正因子法[7-8]等。前者雖從根本上解決了違約問題，但降階處理常伴隨著新的數(shù)值求解精度問題；而后者數(shù)值收斂效果依賴于修正因子的選取。

含摩擦碰撞過程存在復(fù)雜的接觸邊界非光滑問題[9]，例如：接觸面切向摩擦力與切向相對速度之間的非光滑性導(dǎo)致的黏滯-滑移現(xiàn)象、正/逆滑動方向轉(zhuǎn)換與法向壓縮/恢復(fù)切換不同步、以及摩擦微滑動等問題[10-11]。針對這些問題，Leine等[12]研究了多剛體系統(tǒng)中的沖擊和摩擦的分支非線性問題。Payr等[13]研究了單桿與非彈性半空間碰撞接觸中的黏-滑、正逆向滑動切換等非光滑問題。Zhao等[14]提出了滿足邊界條件的非線性彈簧阻尼模型和切向接觸模型，考慮了碰撞過程中碰撞體之間的微運動，并分析了碰撞與結(jié)構(gòu)柔性之間的耦合作用。姚文莉等[15]基于滑動模式步進算法研究了含鉸摩擦的多剛體系統(tǒng)的碰撞動力學(xué)問題，考慮了切向模式的復(fù)雜性，解決了碰撞前后能量的不協(xié)調(diào)性。對于在柔性體中的含摩擦碰撞黏滯-滑移動力學(xué)建模和分析還有待進一步探索。

因此，本文運用附加約束法建立了柔性多體系統(tǒng)的含摩擦碰撞動力學(xué)模型。將系統(tǒng)分為非碰撞、碰撞初始接觸、切向黏滯碰撞以及切向滑動碰撞四個狀態(tài)，分別建立不同的動力學(xué)方程。碰撞初始條件采用沖量/動量法求解。根據(jù)黏滯/滑移切換條件判斷，采用增減Lagrange乘子的附加約束法求解碰撞過程。并給出了仿真算例，分析了含摩擦碰撞的黏滯/滑移切換現(xiàn)象，驗證了所提出模型和算法的有效性。

1 接觸描述

在重力場下作回轉(zhuǎn)運動的柔性機械臂上一點P與剛性碰撞面發(fā)生點-面斜碰撞，如圖1所示。在柔性機械臂轉(zhuǎn)軸O處建立描述大范圍運動的慣性坐標(biāo)系O-XY，碰撞面坐標(biāo)系為O′-X′Y′，用浮動坐標(biāo)系O-xy描述柔性梁的運動與變形。桿件的總長度為L，大范圍運動轉(zhuǎn)角為θ，F(xiàn)N和FT分別為法向、切向接觸力。

圖1 柔性機械臂碰撞示意圖Fig.1 The flexible beam collision schematic diagram

令WP為慣性基到浮動基的轉(zhuǎn)換矩陣[16]，hP為柔性機械臂運動變形后點P(xP,0)在浮動基中的位置矢量，表示如下：

(1)

式中：ux為軸向高階耦合變形項；uy為橫向彎曲變形項。令碰撞面單位外法向量和切向量分別為nP=[nx,ny]T，τP=[τx,τy]T。

令rQ為Q點在碰撞面中的位置矢量，WQ為慣性基到碰撞面坐標(biāo)系的轉(zhuǎn)換矩陣。經(jīng)過推導(dǎo)歸納，得到碰撞法向間距矢量gN和切向間距矢量gT為

(2)

(3)

式中：Z(xP)為軸向高次耦合形函數(shù)矩陣；N為模態(tài)截斷數(shù)；Φy為橫向模態(tài)函數(shù)行陣；其中B(t)為橫向模態(tài)坐標(biāo)列陣。

2 系統(tǒng)全局動力學(xué)

2.1 系統(tǒng)碰撞初始條件的求解

(4)

(nkxcosθ+nkysinθ)Φy(xP)B

(nkxsinθ-nkycosθ)φyi(xP)

(5)

同理，得到切向列陣dT=[dTθ,dTB1,…dTBN]T其中元素為

(-τkxcosθ-τkysinθ)Φy(xP)B

(-τkxsinθ+τkycosθ)φyi(xP)

(6)

(7)

2.2 系統(tǒng)切向黏滯動力學(xué)方程

若接觸對發(fā)生黏滯，系統(tǒng)添加法向和切向接觸約束，可得到系統(tǒng)碰撞黏滯階段的微分-代數(shù)動力學(xué)方程為

(8)

通過推到歸納，黏滯時法向位移約束Jacobi列陣Φz=[Φθ,ΦB1,…,ΦBN]中的元素為

(nkxcosθ+nkysinθ)Φy(xP)B

(-nkxsinθ+nkycosθ)φyi(xP)

(9)

(10)

同理，通過推到歸納，黏滯時切向位移約束Jacobi列陣Γz=[Γθ,ΓB1,…,ΓBN]中的元素為

(-τkxcosθ-τkysinθ)Φy(xP)B

(-τkxsinθ+τkycosθ)φyi(xP)

(11)

(12)

2.3 系統(tǒng)切向滑動動力學(xué)方程

若接觸對切向滑動，法向約束依然存在，而切向約束解除。引入切向滑動摩擦力勢能的概念VTs，根據(jù)庫倫摩擦定律，VTs表達式如下：

(13)

(14)

2.4 黏滯/滑動狀態(tài)切換準(zhǔn)則

黏滯/滑移狀態(tài)切換準(zhǔn)則：連續(xù)接觸狀態(tài)存在黏滯、正滑移、逆滑移等現(xiàn)象，假設(shè)接觸面材料各向同性、正/逆向摩擦因數(shù)相同μ，切向接觸力為λT，則連續(xù)接觸狀態(tài)存在以下4種黏滯/滑動切換情況：

(15)

3 柔性機械臂與斜面含摩擦碰撞仿真算例

柔性梁與斜面發(fā)生碰撞，如圖2所示。柔性梁密度ρ=2 766.7 kg/m3，彈性模量E=68.952 GPa，長度L=1.0 m，截面慣性矩I=7.85×10-9m4，截面積S=3.14×10-4m2。令初始無變形的梁從θ=0.125π位置在重力作用下自由下落，梁末端與面發(fā)生碰撞，碰撞點處切線為2x+3y=0，摩擦因數(shù)μ=0.3。仿真時間為2.0 s，無碰撞階段的步長為1.0×10-4s，碰撞階段的步長為1.0×10-6s，違約修正參數(shù)取α=20，β=20。

圖2 斜碰撞示意圖Fig.2 The schematic diagram of oblique collision

圖3是柔性梁末端在基座坐標(biāo)系中X橫向和Y豎直位移。圖4是整個仿真過程中系統(tǒng)的能量變化。由圖可見：碰撞時能量損失較大；碰撞引發(fā)了柔性梁的高階變形和高頻振動，系統(tǒng)的彈性變形勢能變化劇烈；碰撞結(jié)束后柔性梁向上反彈，動能又逐漸轉(zhuǎn)化為勢能，梁的機械能又保持恒定；多次碰撞后，重力勢能逐漸降低至平衡位置為負值。

圖5為2 s內(nèi)法向碰撞力時程圖。由仿真數(shù)據(jù)得：法向碰撞力首次高峰值發(fā)生在0.373 791 s，第二次高峰值碰撞發(fā)生在1.25 s附近，之間有多次小峰值碰撞。其原因是：碰撞發(fā)生后梁上靠近碰撞點的區(qū)域發(fā)生變形和振動，遠端受到的碰撞影響較小仍維持原來的運動狀態(tài)；隨著時間的推移，碰撞動力學(xué)效應(yīng)逐漸到達向關(guān)節(jié)端傳播，并返回再傳至碰撞點，造成第二次碰撞力峰值。由圖可見：接近靜止?fàn)顟B(tài)(t=1.36～1.37 s)時，碰撞引發(fā)柔性梁的彈性變形和高頻振動，使得在一次接觸中出現(xiàn)多次的次碰撞。

(a) 軸向位移

(b) 橫向位移圖3 梁末端位移Fig.3 The trajectory of link tip

圖4 柔性梁的機械能變化Fig.4 Energy of the system

圖6是切向碰撞力時程圖。如圖所示：首次切向碰撞力高峰值發(fā)生在0.373 791 s，第二次高峰值發(fā)生在1.375 s附近。同樣，切向也存在碰撞力的兩次高峰值不連續(xù)問題，而與法向碰撞力峰值的出現(xiàn)并不同步。接近系統(tǒng)靜止?fàn)顟B(tài)(t=1.68～1.72 s)時，切向碰撞力的正負號變化頗為頻繁，說明頻繁發(fā)生正向滑動、逆向滑動、黏滯轉(zhuǎn)滑動、滑動轉(zhuǎn)黏滯等多樣性接觸狀態(tài)切換。

圖7為碰撞對的切向滑動位移。如仿真數(shù)據(jù)得到：碰撞發(fā)生后碰撞點的法向和切向相對速度振蕩劇烈。與圖6相對應(yīng)的，切向滑動位移的正負號變化頗為頻繁，存在多次正向/逆向滑動、黏滯轉(zhuǎn)滑動、滑動轉(zhuǎn)黏滯等狀態(tài)的頻繁切換。接近系統(tǒng)靜止?fàn)顟B(tài)(t=1.30～1.38 s) 時，碰撞對的滑動位移在微小(黏滯)范圍內(nèi)反復(fù)振蕩，直至靜止。圖8為碰撞對的相對速度切向絕對滑動速度：多次宏觀碰撞、宏觀碰撞中的次碰撞。如圖所示，一次宏觀的碰撞接觸過程中，出現(xiàn)多次次碰撞接觸。而切向接觸在一次宏觀接觸過程中存在多次黏滯-滑動狀態(tài)的頻繁切換。

(a) t=0～2 s內(nèi)

(b) t=1.36～1.37 s內(nèi)圖5 法向碰撞力Fig.5 The normal impact forces

(a) t=0～2 s內(nèi)

(b) t=1.68～1.72 s內(nèi)圖6 切向碰撞力Fig.6 The tangential impact forces

(a) 滑動位移

(b) 微滑動(黏滯)位移圖7 碰撞點切向位移Fig.7 The tangential displacement of contact-pair

(a) 多次宏觀碰撞

(b) 宏觀碰撞中的次碰撞圖8 切向絕對滑動速度Fig.8 The tangential slipping absolute velocity of contact-pair

4 結(jié) 論

針對柔性機械臂含摩擦碰撞黏滯-滑移過程，研究了系統(tǒng)變拓撲過程動力學(xué)的精確建模和自動切換問題。采用沖量-動量法求解碰撞初始時刻的速度跳躍，避免了速度不協(xié)調(diào)造成的反復(fù)脫開-接觸的高頻“假振蕩”數(shù)值問題。引入切向滑動摩擦力勢能的概念得到切向滑動約束方程。給出接觸/分離、黏滯/滑切換準(zhǔn)則，建立不同約束條件的動力學(xué)方程，解決了碰撞初始速度跳躍和切向摩擦引起的數(shù)值不穩(wěn)定和能量不一致問題。通過仿真算例，分析了黏滯/滑移、正/逆向滑動等復(fù)雜切換現(xiàn)象，分析了柔性變形與碰撞的耦合效應(yīng)，驗證了本文提出的模型和算法的有效性。

[1] PFEIFFER F, FOERG M, ULBRICH H. Numerical aspects of non-smooth multibody dynamics[J]. Computer Methods in Applied Mechanics & Engineering, 2006, 195(50/51): 6891-6908.

[2] 董富祥, 洪嘉振. 平面柔性多體系統(tǒng)正碰撞動力學(xué)建模理論研究[J]. 計算力學(xué)學(xué)報, 2010, 27(6): 1042-1048.

DONG Fuxiang, HONG Jiazhen. Study on the modeling theory of the normal impact dynamics for the planar flexible multibody system[J]. Chinese Journal of Computational Mechanics, 2010, 27(6):1042-1048.

[3] 段玥晨, 章定國, 洪嘉振. 作大范圍運動柔性梁的一種碰撞動力學(xué)求解方法[J]. 機械工程學(xué)報, 2012, 48(19):95-102.

DUAN Yuechen, ZHANG Dingguo, HONG Jiazhen. Method for solving the impact problem of a flexible beam with large overall motion[J]. Journal of Mechanical Engineering, 2012, 48(19):95-102.

[4] 王檢耀, 洪嘉振, 劉鑄永. 接觸碰撞動力學(xué)的多變量選取方法[J]. 力學(xué)學(xué)報, 2014, 46(2):318-322.

WANG Jianyao, HONG Jiazhen, LIU Zhuyong. Multi-variable selection method in contact/impact dynamics[J]. Acta Mechanica Sinica, 2014, 46(2):318-322.

[5] QI Z H, LUO X M, HUANG Z H. Frictional contact analysis of spatial prismatic joints in multibody systems[J]. Multibody System Dynamics, 2011, 26(4):441-468.

[6] TAYLOR R L, PAPADOPOULOS P. On a finite element method for dynamic contact/impact problems[J]. International Journal for Numerical Methods in Engineering, 1993, 36(12): 2123-2140.

[7] FLORES P, MACHADO M, SEABRA E, et al. A parametric study on the Baumgarte stabilization method for forward dynamics of constrained multibody systems[J]. Journal of Computational and Nonlinear Dynamics, 2011, 6(1): 011019.

[8] 付士慧, 王琪. 多體系統(tǒng)動力學(xué)方程違約修正的數(shù)值計算方法[J]. 計算力學(xué)學(xué)報, 2007, 24(1):44-49.

FU Shihui, WANG Qi. A numerical method for constraint stabilization of dynamic equations of multi-body systems[J]. Chinese Journal of Computational Mechanics, 2007, 24(1):44-49.

[9] 金棟平, 胡海巖. 碰撞振動及其典型現(xiàn)象[J]. 力學(xué)進展,1999,29(2):155-163.

JIN Dongping, HU Haiyan. Vibro-impacts and their typical behaviors of mechanical systems[J]. Advances in Mechanics, 1999, 29(2):155-163.

[10] 王琪, 莊方方, 郭易圓, 等. 非光滑多體系統(tǒng)動力學(xué)數(shù)值算法的研究進展[J]. 力學(xué)進展, 2013, 43(1): 101-111.

WANG Qi, ZHUANG Fangfang, GUO Yiyuan, et al. Advances in the research on numerical methods for non-smooth dynamics of multibody systems[J]. Advances in Mechanics, 2013, 43(1):101-111.

[11] 丁千, 翟紅梅. 機械系統(tǒng)摩擦動力學(xué)研究進展[J].力學(xué)進展, 2013, 43(1): 112-131.

DING Qian, ZHAI Hongmei. The advance in researches of friction dynamics in mechanics system[J]. Advances in Mechanics, 2013, 43(1):112-131.

[12] LEINE R I, CAMPEN D H V, GLOCKER C H. Nonlinear dynamics and modeling of various wooden toys with impact and friction[J]. Journal of Vibration & Control, 2003, 9(1/2):25-78.

[13] PAYR M, GLOCKER C. Oblique frictional impact of a bar: analysis and comparison of different impact laws[J]. Nonlinear Dynamics, 2005, 41(4):361-383.

[14] ZHAO Z, LIU C S. The analysis and simulation for three-dimensional impact with friction[J]. Multibody System Dynamics, 2007, 18(4):511-530.

[15] 姚文莉, 陳濱, 劉才山,等. 含摩擦的多剛體系統(tǒng)碰撞問題的滑動狀態(tài)步進法[J]. 應(yīng)用數(shù)學(xué)和力學(xué), 2007, 28(12):1448-1454.

YAO Wenli, CHEN Bin, LIU Caishan, et al. Sliding state stepping algorithm for solving impact problems of multi-rigid-body system with joint friction[J]. Applied Mathematics & Mechanics, 2007, 28(12):1448-1454.

[16] 吳勝寶, 章定國. 大范圍運動剛體-柔性梁剛?cè)狁詈蟿恿W(xué)分析[J]. 振動工程學(xué)報, 2011, 24(1): 1-7.

WU Shengbao, ZHANG Dingguo. Rigid-flexible coupling dynamic analysis of hub-flexible beam with large overall motion[J]. Journal of Vibration Engineering, 2011, 24(1):1-7.

Dynamicsimulationforflexiblemultibodysystemscontainingfrictionalimpactandstick-slipprocesses

QIAN Zhenjie1, ZHANG Dingguo2, JIN Chengqian1

(1. Nanjing Research Institute of Agricultural Mechanization of National Ministry of Agriculture, Nanjing 210014, China；2. School of Sciences, Nanjing University of Science and Technology, Nanjing 210094, China)

The dynamic simulation including modeling and automatic switching for flexible multibody systems was investigated using Lagrange multiplier method and the high order rigid-flexible coupled theory. Based on variable topology idea, dynamic equations and the corresponding constraint conditions were established, respectively for detachment, impact, viscous contact and slip contact states. The impulse-momentum method was used to solve the impact initial conditions. The concept of tangential sliding friction force potential energy was introduced to describe the corresponding generalized impact forces with Lagrange equations. The switching criteria among contact, detachment, stick and slip were built to realize the global dynamic automatic switching of the system. Examples were numerically simulated to analyze the complex non-smooth phenomena including forward/backward slip and slip/stick(micro slip) to verify the effectiveness of the proposed model and algorithms.

frictional impact; stick-slip; contact constraints; flexible robot arm; dynamic modeling

國家自然科學(xué)基金青年科學(xué)基金(11602120)；國家自然科學(xué)基金(11272155;11132007)；國家科技支撐計劃(2013BAD08B02)

2016-07-13 修改稿收到日期：2016-09-24

錢震杰 女，博士，助理研究員，1985年生

金誠謙 男，博士，研究員，1973年生

O313；O322

A

10.13465/j.cnki.jvs.2017.23.006