多質(zhì)心有限質(zhì)點(diǎn)法及其在工業(yè)機(jī)器人動(dòng)力學(xué)建模的應(yīng)用

劉目珅， 張飛斌， 王天楊， 褚福磊， 程衛(wèi)東， 劉佑民

(1.清華大學(xué) 機(jī)械工程系，北京 100084； 2. 北京交通大學(xué) 機(jī)械與電子控制工程學(xué)院，北京 100044； 3. 北京航天發(fā)射技術(shù)研究所，北京 100084)

工業(yè)機(jī)器人動(dòng)力學(xué)建模能夠更好地獲取工作過(guò)程中的動(dòng)態(tài)特性，是分析機(jī)器人末端抖動(dòng)，提高定位精度的理論基礎(chǔ)，因此建立一個(gè)能夠精準(zhǔn)滿足實(shí)際工況的動(dòng)力學(xué)模型能為機(jī)械系統(tǒng)動(dòng)態(tài)性能的評(píng)估和優(yōu)化設(shè)計(jì)提供強(qiáng)有力的理論和技術(shù)支撐[1]。在工業(yè)機(jī)器人實(shí)際生產(chǎn)過(guò)程中，工業(yè)機(jī)器人關(guān)節(jié)處電流激勵(lì)頻率與桿件固有頻率重合導(dǎo)致的工業(yè)機(jī)器人末端抖動(dòng)現(xiàn)象時(shí)有發(fā)生，這就需要考慮工業(yè)機(jī)器人桿件柔性的精準(zhǔn)建模，為有效抑制工業(yè)機(jī)器人末端抖動(dòng)提供理論機(jī)理。當(dāng)前工業(yè)機(jī)器人建模方法主要有Lagrange、Newton-Euler以及Kane法，但是以上方法多基于工業(yè)機(jī)器人桿件剛性假設(shè)，若對(duì)機(jī)器人模型考慮桿件柔性進(jìn)行精準(zhǔn)建模，常常會(huì)面臨建模流程推導(dǎo)困難，求解效率低等問(wèn)題[2]。

有限質(zhì)點(diǎn)法[3]是一種基于向量式有限元[4-6]和固體力學(xué)，對(duì)結(jié)構(gòu)進(jìn)行數(shù)值分析的方法，該方法將模型在物理空間上分為有限個(gè)具有質(zhì)量的質(zhì)點(diǎn)，在時(shí)域上，將時(shí)間分解為一系列途徑單元的集合，在每一個(gè)途徑單元內(nèi)結(jié)構(gòu)的分析都是連續(xù)和線性的，對(duì)于質(zhì)點(diǎn)采用基于牛頓運(yùn)動(dòng)學(xué)的控制方程描述運(yùn)動(dòng)狀態(tài)。有限質(zhì)點(diǎn)法通過(guò)逆向運(yùn)動(dòng)的概念分離單元?jiǎng)傮w位移和節(jié)點(diǎn)純變形，進(jìn)而計(jì)算出虛擬的單元節(jié)點(diǎn)內(nèi)力，再通過(guò)正向運(yùn)動(dòng)得到真實(shí)位置的單元節(jié)點(diǎn)內(nèi)力，因此在本質(zhì)上有限質(zhì)點(diǎn)法是將幾何線性和非線性都使用簡(jiǎn)便統(tǒng)一的框架進(jìn)行求解[7]，在求解結(jié)構(gòu)柔性變形上更有優(yōu)勢(shì)。

對(duì)工業(yè)機(jī)器人建模采用實(shí)體單元結(jié)構(gòu)能夠更真實(shí)的反映結(jié)構(gòu)響應(yīng)，有著更好的適應(yīng)性和準(zhǔn)確性。然而基于六面體單元的有限質(zhì)點(diǎn)法[8-9]的逆向運(yùn)動(dòng)步驟較為繁瑣，涉及到多次投影、求模以及向量夾角計(jì)算等。有限質(zhì)點(diǎn)法以單元任一節(jié)點(diǎn)作為逆向平移參考點(diǎn)、逆向轉(zhuǎn)動(dòng)旋轉(zhuǎn)中心，并且逆向轉(zhuǎn)動(dòng)參考平面未能夠充分反映單元節(jié)點(diǎn)信息，不能充分將單元?jiǎng)傮w位移與單元節(jié)點(diǎn)純變形分離，尤其是工業(yè)機(jī)器人這種剛體大轉(zhuǎn)動(dòng)-柔性變形耦合模型，單元?jiǎng)傮w位移對(duì)求解純變形的影響更為嚴(yán)重，使部分節(jié)點(diǎn)內(nèi)力值產(chǎn)生較大偏差，促使單元內(nèi)其余節(jié)點(diǎn)內(nèi)力與之對(duì)應(yīng)，以達(dá)到單元內(nèi)力平衡狀態(tài)，因此在整個(gè)收斂過(guò)程中，單元節(jié)點(diǎn)內(nèi)力值波動(dòng)性增加，收斂性變差，甚至?xí)斐晒?jié)點(diǎn)內(nèi)力值發(fā)散的情況。盡管Hou等[10]已經(jīng)對(duì)基于六面體單元有限質(zhì)點(diǎn)法的逆向運(yùn)動(dòng)步驟進(jìn)行了有效的簡(jiǎn)化，但是未考慮逆向運(yùn)動(dòng)參考點(diǎn)和參考平面的選取對(duì)于節(jié)點(diǎn)內(nèi)力收斂性的影響。

本文提出的多質(zhì)心有限質(zhì)點(diǎn)法與有限質(zhì)點(diǎn)法相比采用了全新的逆向運(yùn)動(dòng)步驟，將平面外、內(nèi)轉(zhuǎn)動(dòng)修改為形式相同、步驟簡(jiǎn)單的兩次逆向轉(zhuǎn)動(dòng)；將單元平面質(zhì)心作為逆向運(yùn)動(dòng)參考點(diǎn)，讓每一個(gè)單元平面質(zhì)心都參與到逆向轉(zhuǎn)動(dòng)所需的參考平面構(gòu)造中，有效降低單元?jiǎng)傮w位移對(duì)節(jié)點(diǎn)純變形求解的影響，提高了對(duì)工業(yè)機(jī)器人這種剛體大轉(zhuǎn)動(dòng)-柔性耦合模型計(jì)算的收斂性。多質(zhì)心有限質(zhì)點(diǎn)法不同于傳統(tǒng)有限單元法以能量法和連續(xù)介質(zhì)力學(xué)作為理論基礎(chǔ)，而是基于向量式結(jié)構(gòu)力學(xué)[11]，以牛頓運(yùn)動(dòng)定律對(duì)質(zhì)點(diǎn)的軌跡進(jìn)行描述，對(duì)結(jié)構(gòu)實(shí)現(xiàn)了真正意義上的離散，在本質(zhì)上更適合處理工業(yè)機(jī)器人動(dòng)力學(xué)問(wèn)題。

基于本文方法對(duì)簡(jiǎn)單的懸臂梁模型編寫(xiě)MATLAB程序，將模型結(jié)果與理論公式以及ABAQUS軟件計(jì)算結(jié)果對(duì)比，驗(yàn)證本文理論和算法的正確性。進(jìn)一步將本文方法結(jié)合正運(yùn)動(dòng)學(xué)引入到工業(yè)機(jī)器人動(dòng)力學(xué)建模中，基于MATLAB將建模程序與GUI工具箱建立聯(lián)系，構(gòu)建工業(yè)機(jī)器人動(dòng)力學(xué)建模軟件，將建模結(jié)果與實(shí)驗(yàn)數(shù)據(jù)相對(duì)比，驗(yàn)證本文方法對(duì)工業(yè)機(jī)器人動(dòng)力學(xué)建模的準(zhǔn)確性和實(shí)用性。

1 多質(zhì)心有限質(zhì)點(diǎn)法基本理論

1.1 單元節(jié)點(diǎn)位移

多質(zhì)心有限質(zhì)點(diǎn)法通過(guò)一系列質(zhì)點(diǎn)和途徑單元的集合，使用基于牛頓運(yùn)動(dòng)學(xué)的控制方程描述質(zhì)點(diǎn)運(yùn)動(dòng)，是一種真正意義上的離散。對(duì)于結(jié)構(gòu)的行為可以使用一組運(yùn)動(dòng)方程描述

(1)

結(jié)合顯式中心差分公式，將運(yùn)動(dòng)方程改寫(xiě)為

(2)

式中：d為單元節(jié)點(diǎn)位移向量；n為迭代步數(shù)；h為步長(zhǎng)； 節(jié)點(diǎn)集中載荷P=F+f。

當(dāng)?shù)綌?shù)n=0時(shí)，也就是運(yùn)動(dòng)時(shí)間t=t0=0時(shí)，對(duì)應(yīng)第-1步虛擬的位移向量d-1為

(3)

接下來(lái)的迭代步(n≥1)都是以式(2)作為單元節(jié)點(diǎn)位置計(jì)算的控制方程，每一步計(jì)算所需要的節(jié)點(diǎn)內(nèi)力f可以根據(jù)單元節(jié)點(diǎn)純變形計(jì)算得出。

1.2 單元逆向運(yùn)動(dòng)

接下來(lái)在一個(gè)途徑單元tn～tn+1時(shí)刻內(nèi)對(duì)多質(zhì)心有限質(zhì)點(diǎn)法逆向運(yùn)動(dòng)步驟詳細(xì)描述，對(duì)應(yīng)迭代步數(shù)為第n步～第n+1步。

圖1 六面體單元運(yùn)動(dòng)和變形示意圖Fig.1 Motion and deformation of a hexahedral element

1.2.1 逆向平移

將tn+1時(shí)刻單元以單元質(zhì)心位置矢量un+1為參考點(diǎn)整體逆向平移至tn時(shí)刻單元質(zhì)心位置矢量un上，如圖2所示，經(jīng)過(guò)逆向平移的tn+1時(shí)刻單元節(jié)點(diǎn)記為a′-b′-c′-d′-e′-f′-g′-h′，結(jié)點(diǎn)i位置矢量d′i為

圖2 逆向平移后的單元位置Fig.2 Element position after reverse translation

(4)

1.2.2 第一次逆向轉(zhuǎn)動(dòng)

將經(jīng)過(guò)逆向平移的單元水平方向上的3個(gè)平面質(zhì)心連接，構(gòu)成三角形參考平面1′-4′-6′，其法向量為Vec′1，同理，tn時(shí)刻單元的參考平面1-4-6，其法向量為Vec1。構(gòu)造的單元參考平面及其法向量，如圖3所示。

圖3 單元參考平面及其法向量示意圖Fig.3 Reference plane and its normal vector of an element

Vec1與Vec′1的夾角θA及其法向量nvecA分別為

(5)

(6)

計(jì)算出第一次逆向轉(zhuǎn)動(dòng)的轉(zhuǎn)動(dòng)矩陣R′(-θA)

R′(-θA)=[1-cos(-θA)]V′2+sin(-θA)V′

(7)

以tn時(shí)刻的單元質(zhì)心位置矢量un為旋轉(zhuǎn)中心對(duì)節(jié)點(diǎn)a′-b′-c′-d′-e′-f′-g′-h′進(jìn)行第一次逆向轉(zhuǎn)動(dòng)，如圖4所示。獲得新的單元節(jié)點(diǎn)位置為a″-b″-c″-d″-e″-f″-g″-h″，單元結(jié)點(diǎn)i位置矢量d″i為

圖4 經(jīng)第一次逆向旋轉(zhuǎn)的單元位置Fig.4 Position of the element after the first reverse rotation

d″i=R′(-θA)d′i,i=a,b,c,…,g,h

(8)

1.2.3 第二次逆向轉(zhuǎn)動(dòng)

將經(jīng)過(guò)第一次逆向轉(zhuǎn)動(dòng)的單元豎直方向上其余3個(gè)平面質(zhì)心連接構(gòu)造三角形參考平面2′-3′-5′，其法向量為Vec′2，同理，tn時(shí)刻單元的參考平面2-3-5，其法向量為Vec2。其余逆向轉(zhuǎn)動(dòng)步驟與第一次逆向轉(zhuǎn)動(dòng)相同，Vec2與Vec′2的夾角及其法向量分別為：θB和nvecB，第二次逆向轉(zhuǎn)動(dòng)的轉(zhuǎn)動(dòng)矩陣為：R″(-θB)。經(jīng)過(guò)第二次逆向轉(zhuǎn)動(dòng)的單元位置及以tn時(shí)刻單元質(zhì)心為原點(diǎn)的局部坐標(biāo)系，如圖5所示。單元節(jié)點(diǎn)位置為a?-b?-c?-d?-e?-f?-g?-h?，位置矢量d?i為

圖5 經(jīng)第二次逆向旋轉(zhuǎn)的單元位置及其局部坐標(biāo)系Fig.5 Position of the element after the second reverse rotation and its local coordinate

d?i=R″(-θB)d″i,i=a,b,c,…,g,h

(9)

(10)

需要特別說(shuō)明的是，在整個(gè)逆向運(yùn)動(dòng)過(guò)程中，逆向平移是以單元質(zhì)心為參考點(diǎn)整體逆向平移，逆向轉(zhuǎn)動(dòng)中單元節(jié)點(diǎn)都是繞同一旋轉(zhuǎn)向量旋轉(zhuǎn)相同的角度，保證了六面體單元的各個(gè)節(jié)點(diǎn)之間不會(huì)發(fā)生相對(duì)位移。將本文推導(dǎo)過(guò)程與有限質(zhì)點(diǎn)法進(jìn)行對(duì)比，主要差別如表1所示。由表1可知，本文方法有效降低了計(jì)算的復(fù)雜度，簡(jiǎn)化了計(jì)算流程。

表1 本文方法與有限質(zhì)點(diǎn)法復(fù)雜度對(duì)比Tab.1 The complexity of the proposed method is compared with that of the finite particle method

1.3 單元節(jié)點(diǎn)內(nèi)力

單元節(jié)點(diǎn)內(nèi)力通過(guò)單元節(jié)點(diǎn)純變形進(jìn)行計(jì)算，計(jì)算方法可參考有限單元法[12]。如圖5所示，以tn時(shí)刻單元質(zhì)心作為原點(diǎn)建立局部坐標(biāo)系O′-x1y1z1，坐標(biāo)軸方向與全局坐標(biāo)系保持一致，將全局坐標(biāo)系下的單元節(jié)點(diǎn)位置矢量轉(zhuǎn)化到局部坐標(biāo)系中

(11)

對(duì)時(shí)刻的單元建立形函數(shù)，建立方式和表達(dá)式形式與傳統(tǒng)有限元方法相同，形函數(shù)坐標(biāo)系記為O″-rst，如圖6所示。

圖6 八節(jié)點(diǎn)六面體單元形函數(shù)坐標(biāo)系Fig.6 Eight-node hexahedral element shape function coordinate

八節(jié)點(diǎn)六面體單元的形函數(shù)Ni表達(dá)式為

(12)

式中：r,s,t分別為形函數(shù)表達(dá)式參數(shù)；ri=±1，si=±1，ti=±1為形函數(shù)坐標(biāo)系內(nèi)節(jié)點(diǎn)i的坐標(biāo)，i=a,b,c,…,g,h。

將局部坐標(biāo)系內(nèi)坐標(biāo)轉(zhuǎn)換到形函數(shù)坐標(biāo)系中，坐標(biāo)系的轉(zhuǎn)換通過(guò)雅可比矩陣J來(lái)實(shí)現(xiàn)

(13)

對(duì)雅可比矩陣求逆，得到形函數(shù)關(guān)于局部坐標(biāo)的導(dǎo)數(shù)為

(14)

(15)

(16)

利用單元變形描述的虛功與應(yīng)力應(yīng)變描述的虛功相等這一原理建立等式

(17)

式中：Va為形函數(shù)積分參數(shù)； δ為變分符號(hào)。

(18)

(19)

式中，R′(θA)和R″(θB)分別為第一次和第二次正向運(yùn)動(dòng)轉(zhuǎn)動(dòng)矩陣。

1.4 步 長(zhǎng)

對(duì)結(jié)合顯示中心差分推導(dǎo)出的控制方程來(lái)說(shuō)，步長(zhǎng)的大小尤為重要，一般在問(wèn)題的分析中，對(duì)步長(zhǎng)的要求主要有兩方面：①物理行為方面，在結(jié)構(gòu)分析過(guò)程中，使用足夠的數(shù)據(jù)點(diǎn)來(lái)滿足計(jì)算精度的要求；②數(shù)值計(jì)算方面，差分公式是函數(shù)微分的一種近似計(jì)算，必然有一定誤差，所以要避免迭代計(jì)算中的誤差累計(jì)造成的數(shù)據(jù)發(fā)散。

根據(jù)上述對(duì)步長(zhǎng)的要求得知，步長(zhǎng)必須要小于一個(gè)步長(zhǎng)臨界值hd，才能得到一個(gè)收斂、準(zhǔn)確的計(jì)算結(jié)果。然而模型類型不同其步長(zhǎng)臨界值也各異，并且精確求解臨界步長(zhǎng)也是困難的。丁承先等以簡(jiǎn)單的桿單元為對(duì)象，提出了一種有效的估算臨界步長(zhǎng)值的方法。

1.5 質(zhì)量和加速度

六面體單元質(zhì)量分布遵循一個(gè)原則：將單元質(zhì)量均分到各個(gè)節(jié)點(diǎn)上，再采用集成的方式就可得到單元節(jié)點(diǎn)質(zhì)量，此時(shí)的單元節(jié)點(diǎn)只考慮質(zhì)量，忽略節(jié)點(diǎn)尺寸[13]。

單元節(jié)點(diǎn)的加速度同樣可以通過(guò)差分的形式給出

(20)

2 結(jié)構(gòu)算例

2.1 懸臂梁模型

為驗(yàn)證本文提出的多質(zhì)心有限質(zhì)點(diǎn)法的準(zhǔn)確性，現(xiàn)以簡(jiǎn)單的懸臂梁模型為例，基于本文方法使用MATLAB軟件進(jìn)行編程，對(duì)結(jié)構(gòu)的變形進(jìn)行分析。

使用六面體單元對(duì)懸臂梁模型劃分網(wǎng)格，網(wǎng)格數(shù)量為320個(gè)，單元節(jié)點(diǎn)為525個(gè)，如圖7所示。懸臂梁長(zhǎng)、寬、高分別為0.5 m，0.1 m和0.1 m。材料類型為彈性，彈性模量為E=206 GPa，密度為Rho=7 800 kg/m3，泊松比為rp=0.29。分析時(shí)間TE=0.05 s，步長(zhǎng)為h=1×10-6s，末端集中載荷P=-3 N，采用斜坡-平臺(tái)方式緩慢加載，加載時(shí)間rtime=0.005 s，并使用阻尼系數(shù)收斂計(jì)算結(jié)果，阻尼系數(shù)為zeta=100。當(dāng)前模型變形為小變形，所以依然適用于材料力學(xué)理論公式

圖7 懸臂梁六面體單元模型Fig.7 Hexahedral elements of cantilever model

(21)

式中：L為懸臂梁的長(zhǎng)；I為截面慣性矩， 對(duì)于正方形截面，I=B4/12，B為正方形截面的邊長(zhǎng)。

沿懸臂梁中線各節(jié)點(diǎn)的豎直方向的位移結(jié)果，如圖8所示。橫坐標(biāo)為中線各節(jié)點(diǎn)與固定端的距離，縱坐標(biāo)為中線各節(jié)點(diǎn)豎直方向上產(chǎn)生的位移，可以看出豎直方向上位移與材料力學(xué)理論公式及采用C3D8R六面體單元的ABAQUS軟件的計(jì)算結(jié)果基本保持一致。

圖8 懸臂梁中線各節(jié)點(diǎn)豎直方向位移結(jié)果Fig.8 Vertical displacement results of each node in the middle line of cantilever beam

本文方法和有限質(zhì)點(diǎn)法對(duì)懸臂梁自由端中心節(jié)點(diǎn)豎直方向位移變化過(guò)程的收斂性對(duì)比，如圖9所示。以理論公式結(jié)果為標(biāo)準(zhǔn)值，在相同網(wǎng)格數(shù)量前提下，本文方法、有限質(zhì)點(diǎn)法、理論公式以及ABAQUS軟件的自由端中心節(jié)點(diǎn)豎向方向位移的對(duì)比結(jié)果，如表2所示。

表2 結(jié)果對(duì)比Tab.2 Results contrast

圖9 本文方法和有限質(zhì)點(diǎn)法對(duì)懸臂梁自由端中心 節(jié)點(diǎn)豎直方向位移變化過(guò)程的收斂性對(duì)比Fig.9 Convergence comparison between the method in this paper and the finite mass method on the vertical displacement change process of the center node of the free end of the cantilever beam

由以上結(jié)果證明本文方法有著良好的準(zhǔn)確性和收斂性。

2.2 工業(yè)機(jī)器人動(dòng)力學(xué)模型

首先對(duì)建模對(duì)象RB08A3機(jī)器人進(jìn)行實(shí)驗(yàn)，發(fā)現(xiàn)存在明顯抖動(dòng)現(xiàn)象，在工業(yè)機(jī)器人連桿4末端安裝加速度傳感器，測(cè)得運(yùn)動(dòng)過(guò)程中執(zhí)行末端加速度數(shù)據(jù)，并記錄J2，J3關(guān)節(jié)電流變化值，如圖10所示。

圖10 工業(yè)機(jī)器人及加速度傳感器布置位置Fig.10 Industrial robot and acceleration sensor placement

為分析抖動(dòng)現(xiàn)象的原因，首先將關(guān)節(jié)電流數(shù)據(jù)分別做頻域分析，結(jié)果如圖11所示。

圖11 關(guān)節(jié)電流頻域圖Fig.11 Frequency domain diagram of joint current

對(duì)實(shí)驗(yàn)測(cè)得的加速度數(shù)據(jù)進(jìn)行頻域分析，可知機(jī)器人抖動(dòng)是因?yàn)樵?3.5 Hz和27.5 Hz左右的電流激勵(lì)頻率和工業(yè)機(jī)器人固有頻率重合。下面采用多質(zhì)心有限質(zhì)點(diǎn)法對(duì)工業(yè)機(jī)器人進(jìn)行動(dòng)力學(xué)建模，建模步驟如圖12所示。

圖12 工業(yè)機(jī)器人動(dòng)力學(xué)建模流程圖Fig.12 Flow chart of industrial robot dynamics modeling

工業(yè)機(jī)器人結(jié)構(gòu)基本參數(shù)，如表3所示。使用六面體單元對(duì)機(jī)器人劃分網(wǎng)格，建立的物理模型如圖13所示，所需單元個(gè)數(shù)為960。

表3 RB08A3機(jī)器人結(jié)構(gòu)參數(shù)Tab.3 RB08A3 robot structure parameters

圖13 工業(yè)機(jī)器人物理模型Fig.13 Physical model of industrial robot

設(shè)定機(jī)器人動(dòng)力學(xué)建?；緟?shù)： 運(yùn)動(dòng)時(shí)間為 2.0 s，迭代步長(zhǎng)為1/3 125 s，單元材料密度為7.8×103，彈性模量為206 GPa，泊松比為0.29，根據(jù)單元的連接類型，確定單元節(jié)點(diǎn)質(zhì)量，運(yùn)動(dòng)狀態(tài)為J2和J3關(guān)節(jié)轉(zhuǎn)動(dòng)，其余關(guān)節(jié)保持固定。

為方便計(jì)算，引入固定的整體坐標(biāo)系和隨機(jī)械臂旋轉(zhuǎn)的坐標(biāo)系[14]，結(jié)合模型的單元與節(jié)點(diǎn)信息，以及關(guān)節(jié)輸出角度變化值，根據(jù)正運(yùn)動(dòng)學(xué)求解出機(jī)器人模型剛體運(yùn)動(dòng)的位置，然后對(duì)模型使用多質(zhì)心有限質(zhì)點(diǎn)法進(jìn)行分析，得到執(zhí)行末端變形量。進(jìn)而將變形量疊加到執(zhí)行末端，重復(fù)上述步驟的計(jì)算，直至分析時(shí)間達(dá)到運(yùn)動(dòng)時(shí)間。使用式(20)可以根據(jù)執(zhí)行末端位置求解出末端加速度變化。

基于MATLAB軟件將工業(yè)機(jī)器人動(dòng)力學(xué)建模程序與MATLAB軟件中GUI工具箱連接，搭建工業(yè)機(jī)器人動(dòng)力學(xué)建模軟件，軟件主界面如圖14所示。

圖14 軟件主界面Fig.14 Main interface of software

進(jìn)而對(duì)仿真測(cè)得的執(zhí)行末端z和y方向上的加速度數(shù)據(jù)頻域分析，并與實(shí)驗(yàn)數(shù)據(jù)對(duì)比，對(duì)比結(jié)果如圖15和圖16所示。

圖15 建模與實(shí)驗(yàn)在z方向上的加速度頻域圖Fig.15 Modeling and experiments in the z direction acceleration frequency domain diagram

圖16 建模與實(shí)驗(yàn)在y方向上的加速度頻域圖Fig.16 Modeling and experiments in the y direction acceleration frequency domain diagram

在z方向上的建模結(jié)果具有13.5 Hz和28.5 Hz，相對(duì)實(shí)驗(yàn)結(jié)果14.0 Hz和27.5 Hz的誤差分別為-3.6%和3.6%。在y方向上的建模結(jié)果具有13.5 Hz，相對(duì)實(shí)驗(yàn)結(jié)果14.0 Hz的誤差為-3.6%。根據(jù)以上結(jié)果可以得知，建模與實(shí)驗(yàn)結(jié)果在低階固有頻率基本重合，高階能量趨勢(shì)基本相同。

3 結(jié) 論

本文提出多質(zhì)心有限質(zhì)點(diǎn)法，將單元質(zhì)心不僅作為逆向平移的參考點(diǎn)，也作為逆向轉(zhuǎn)動(dòng)的旋轉(zhuǎn)中心，使每個(gè)單元平面都參與到參考平面的構(gòu)造中，有效減弱了單元?jiǎng)傮w位移對(duì)節(jié)點(diǎn)純變形求解的影響，進(jìn)而解決了工業(yè)機(jī)器人這種剛體大轉(zhuǎn)動(dòng)-柔性變形耦合模型中因單元節(jié)點(diǎn)純變形計(jì)算偏差，造成的節(jié)點(diǎn)內(nèi)力值收斂性差的問(wèn)題。

與傳統(tǒng)有限質(zhì)點(diǎn)法相比，多質(zhì)心有限質(zhì)點(diǎn)法的兩次逆向轉(zhuǎn)動(dòng)，除選取的參考平面節(jié)點(diǎn)不同以外，都是采用簡(jiǎn)單的平面外轉(zhuǎn)動(dòng)的旋轉(zhuǎn)方式，大幅度減少了逆向運(yùn)動(dòng)過(guò)程中的向量求模、叉乘、投影等步驟次數(shù)，有效簡(jiǎn)化了計(jì)算流程，使得整體框架更為統(tǒng)一。

對(duì)簡(jiǎn)單懸臂梁模型使用多質(zhì)心有限質(zhì)點(diǎn)法自編程序建模分析，并將分析結(jié)果與有限質(zhì)點(diǎn)法、理論公式和ABAQUS軟件分析結(jié)果對(duì)比，結(jié)果表明本文方法有著良好的準(zhǔn)確性和收斂性。

將多質(zhì)心有限質(zhì)點(diǎn)法結(jié)合正運(yùn)動(dòng)學(xué)對(duì)工業(yè)機(jī)器人動(dòng)力學(xué)建模，通過(guò)MATLAB軟件將建模程序與GUI圖形化界面結(jié)合，構(gòu)造基于多質(zhì)心有限質(zhì)點(diǎn)法的工業(yè)機(jī)器人動(dòng)力學(xué)建模軟件，并將建模結(jié)果與實(shí)驗(yàn)結(jié)果相對(duì)比，結(jié)果表明使用本文方法的建模結(jié)果能夠有效模擬出工業(yè)機(jī)器人因電流激勵(lì)頻率與固有頻率重合產(chǎn)生的抖動(dòng)頻率。