基于MATLAB/Simulink的四連桿機(jī)構(gòu)動力學(xué)仿真及改進(jìn)

張建軍，黃慶學(xué)，馬麗楠

(太原科技大學(xué)重型機(jī)械教育部工程研究中心，太原 030024)

摘要：針對四連桿機(jī)構(gòu)的動力學(xué)仿真，提出了一種通用化的改進(jìn)方案。只需輸入連桿機(jī)構(gòu)的部分初始參數(shù)，再執(zhí)行一個腳本程序，即可進(jìn)行四連桿機(jī)構(gòu)的動力學(xué)仿真。此外，額外提供了一種求解機(jī)構(gòu)位置解的新方法。輸入一組容易獲取且誤差較小的初始角度參數(shù)，再執(zhí)行一個基于牛頓辛普森算法的m函數(shù)文件，即可得到連桿機(jī)構(gòu)的位移曲線與其相應(yīng)的誤差曲線。這不僅為求解機(jī)構(gòu)位置解開辟了新方法，而且還提供了精度比較高的初始值參數(shù)，進(jìn)一步保證了動力學(xué)仿真的可靠運(yùn)行。

關(guān)鍵詞：四連桿機(jī)構(gòu)；動力學(xué)仿真；通用性改進(jìn)

收稿日期：2015-01-06

基金項目：山西科技創(chuàng)新團(tuán)隊項目(20131006)；山西省自然科學(xué)

作者簡介：張建軍(1987-)，男，碩士研究生，主要研究方向為機(jī)電系統(tǒng)控制及其自動化。

中圖分類號：TG659文獻(xiàn)標(biāo)志碼：B

1四連桿機(jī)構(gòu)動力學(xué)數(shù)學(xué)建模及仿真模型

1.1運(yùn)動學(xué)部分

在對四連桿機(jī)構(gòu)進(jìn)行動力學(xué)分析之前，需先進(jìn)行運(yùn)動學(xué)分析，分析結(jié)果主要包括以下3個矢量方程：位置方程、速度方程及加速度方程[1-3]。在用矢量法建立機(jī)構(gòu)的位置方程時，需將構(gòu)件用矢量來表示，并做出機(jī)構(gòu)的封閉矢量多邊形[4]，如圖1所示。由矢量多邊形：

得其復(fù)數(shù)表達(dá)式：

r2eiθ2+r3eiθ3=r1eiθ1+r4eiθ4

(1)

圖1　四連桿機(jī)構(gòu)運(yùn)動簡圖

利用歐拉公式將式(1)展開，有：

(2)

式(2)是關(guān)于θ3、θ4非線性超越方程組，其解由三角函數(shù)的關(guān)系可求出，但該解法卻不適合計算機(jī)編程實現(xiàn)，現(xiàn)采用數(shù)值解法，并借助MATLAB工具編程實現(xiàn)。對位置方程組(2)，根據(jù)求解非線性超越方程的數(shù)值解法——牛頓辛普森算法[5-6]，編寫相應(yīng)的位置求解函數(shù)m文件，求出四連桿機(jī)構(gòu)的位置初始值。

得到位置初始值后，接下來需求解機(jī)構(gòu)的速度初始值。對式(2)求兩次導(dǎo)數(shù)，依次得到該機(jī)構(gòu)的速度與加速度方程，這里直接給出其相應(yīng)的矩陣表達(dá)式，即式(3)和式(4).

(3)

(4)

把位置初值代入式(3)，求得速度初值。根據(jù)式(4)，編寫m文件并搭建Simulink模型，如圖2所示，crank.m與rrrki.m分別為曲柄和從動桿組的運(yùn)動仿真模塊。用已求得的位置、速度值作為相應(yīng)積分模塊的初始值，設(shè)置好仿真時間，即可進(jìn)行運(yùn)動學(xué)仿真。

圖2　運(yùn)動學(xué)仿真Simulink模型

1.2動力學(xué)部分

在運(yùn)動學(xué)分析的基礎(chǔ)上，對該機(jī)構(gòu)進(jìn)行動力學(xué)分析。先對圖1中四連桿機(jī)構(gòu)中的活動構(gòu)件進(jìn)行受力分析，如圖3所示。列出相應(yīng)的力平衡方程，如式(5)~式(7).

圖3　各活動構(gòu)件的受力分析圖

原動件AB的力平衡方程：

(5)

連桿BC的力平衡方程：

(6)

從動件CD的力平衡方程：

(7)

式(5)~式(7)中，m2～m4及J2～J4分別為活動構(gòu)件2~4的平動與繞質(zhì)心Si(i=2，3，4)的轉(zhuǎn)動慣量，F(xiàn)xi、Fyi以及Mi(i=2，3，4)分別為活動構(gòu)件質(zhì)心處所受的外力在實虛軸上的分力以及外力矩，ri與rci(i=2，3，4)為各活動構(gòu)件的桿長及質(zhì)心位置值，M1為施加在原動件AB上的平衡力矩，其余符號參考圖3.

為便于計算機(jī)編程，將式(5)~式(7)整理成矩陣形式(8)與式(9)：

式(8)為原動件AB的力平衡方程(5)的矩陣形式，式(9)為從動件BC及連桿CD的力平衡方程(6)、(7)的矩陣形式。結(jié)合該機(jī)構(gòu)的運(yùn)動學(xué)仿真模型，根據(jù)式(8)與式(9)，建立四連桿機(jī)構(gòu)的動力學(xué)仿真模型，如圖4所示，式(8)與式(9)分別用來構(gòu)建曲柄動力學(xué)仿真模塊crankdy.m與從動桿組動力學(xué)仿真模塊rrrrdy.m：

設(shè)置圖4中各模塊的值，即可進(jìn)行仿真。通過仿真發(fā)現(xiàn)，積分模塊的初始值參數(shù)對整個仿真有著決定性的作用。然而，這些初值的獲得，需要之前逐次調(diào)用位置求解函數(shù)以及角速度求解函數(shù)。如果Simulink模型文件設(shè)置有變化(比如常數(shù)模塊改變)，還得重新求取位置與速度的初始值，這將大大增加仿真前期的工作量[7-10]。因此，仿真的通用性不強(qiáng)，有必要對原來的仿真進(jìn)行改進(jìn)，使其函數(shù)調(diào)用靈活且重復(fù)修改量少，本文的下一節(jié)將提供這些思想的MATLAB實現(xiàn)。

(8)

(9)

圖4　動力學(xué)仿真Simulink模型

2四連桿機(jī)構(gòu)運(yùn)動學(xué)仿真通用性的改進(jìn)及其實例驗證

所謂的通用性是指，對任意的四連桿機(jī)構(gòu)，只需在MATLAB命令窗口給出桿長矩陣、角度初值矩陣、原動件運(yùn)動及動力參數(shù)等，通過調(diào)用一段程序代碼，即可實現(xiàn)動力學(xué)仿真。為了提高原仿真的通用性，本節(jié)提供一種改進(jìn)方案，并給出相應(yīng)的實例。

2.1改進(jìn)方案及其MATLAB實現(xiàn)

首先，需先設(shè)置Simulink模型文件中的積分模塊的初始值。在Simulink中有兩種方法建立其仿真的初始條件。第一種，就是直接在相關(guān)屬性欄中填入數(shù)字。遺憾的是，這種方法對于更復(fù)雜的初始條件會變得非常笨拙，且不適合通用化程序設(shè)計的要求。另一種方法是為MATLAB環(huán)境變量設(shè)置初始值，即用變量名代替數(shù)字，該變量的值由MATLAB工作空間提供，因此需在仿真前在MATLAB的命令窗口給同樣變量名的變量賦值。使用第二種方法，很容易在MATLAB環(huán)境中用一系列命令改變初始值，而無需在每次初始條件改變時，返回Simulink模型圖，修改相應(yīng)模塊的設(shè)置。本節(jié)將采取第二種方法。

其次，也是比較重要的，編制一個腳本m文件，以調(diào)用其他m函數(shù)文件。該m文件中用代碼依次調(diào)用了處于同一工作目錄的位置初值求解函數(shù)、速度初值求解函數(shù)以及Simulink模型文件。通過m文件的內(nèi)部調(diào)用，符合模塊化編程思想，調(diào)試方便，減少了重復(fù)編程量。此外，為了給程序調(diào)試帶來便捷，考慮到要單獨調(diào)用某個m函數(shù)文件，故將以上各m文件單獨保存在同一工作目錄，而沒有將其作為子程序整合到腳本文件中。

2.2實例佐證

如圖3所示，鉸鏈四連桿機(jī)構(gòu)各構(gòu)件的幾何尺寸為r1=1 200 mm，r2=400 mm，r3=1 000 mm，r4=700 mm；質(zhì)心為rc2=200 mm，rc3=500 mm，rc4=350 mm；質(zhì)量為m2=1.2 kg，m3=3 kg，m4=2.2 kg；繞質(zhì)心的轉(zhuǎn)動慣量J2=0.016 kg/m2，J3=0.25 kg/m2，J4=0.09 kg/m2，構(gòu)件4的工作阻力矩M4=100 N/m，順時針方向，其他構(gòu)件所受外力和外力矩均為零，且不考慮摩擦。構(gòu)件2(即原動件)以等角速度ω2=10 rad/s逆時針旋轉(zhuǎn)。

設(shè)置完成后，在MATLAB命令提示符下輸入：

th=[0，45，100，0]；rs=[0.4，1，0.7，1.2]；

om20=10；

該品種由福建省泉州市農(nóng)業(yè)科學(xué)研究所（福建省泉州市晉江市池店鎮(zhèn)泉安北路，郵編：362212）和云南省農(nóng)業(yè)科學(xué)院經(jīng)濟(jì)作物研究所培育。

r10=[1，0.7，0.5，0.35，3，2.2，0.25，0.09，0，0，0，0，0]；

cr10=[0.4，0.2，1.2，0.016，0，0，0，0，0]；ts=1；Ini；

并回車即可。其動力學(xué)仿真結(jié)果如圖5和圖6所示。

圖5依次為連桿3和從動件4的角位移、角速度和角加速以及鉸接副C的加速度分量隨時間的變化曲線，圖6依次為鉸接副A的水平與豎直方向的分力和曲柄AB上的平衡力矩及其所做的功隨時間的變化曲線。

綜上所述，為提高動力學(xué)仿真程序的通用性，本節(jié)僅通過在外部調(diào)用一個腳本文件Ini.m，實現(xiàn)了對任意四連桿機(jī)構(gòu)的動力學(xué)仿真。整個操作過程簡單且通用性強(qiáng)，如果仿真條件有變化(如仿真時間或者初始位置預(yù)估值變化)，只需在外部命令上進(jìn)行修改，而無需修改相應(yīng)的m函數(shù)文件或者Simulink模型文件。

3四連桿機(jī)構(gòu)的位置解求解新途徑

上一節(jié)介紹了動力學(xué)仿真的通用化改進(jìn)——只要給定連桿長度及初始位置等參數(shù)，即可得到整個機(jī)構(gòu)的動力學(xué)動態(tài)參數(shù)。該機(jī)構(gòu)的動態(tài)位置解，是反饋的加速度值經(jīng)過積分模塊處理后的輸出，該輸出又作為下一步仿真的輸入，如此循環(huán)下去，最終得到全程的位置解曲線。本節(jié)將介紹另外一種求解機(jī)構(gòu)位置信息的途徑，該方案以一組容易獲得且誤差較小的初始位置解為起點，將原動件等間距地旋轉(zhuǎn)一周，通過反復(fù)使用牛頓辛普森方法，進(jìn)而求出若干位置的位置解。

以圖1的四連桿機(jī)構(gòu)為例，其初始條件以容易獲取為選擇標(biāo)準(zhǔn)。編寫一個MATLAB函數(shù)文件fourbar.m，在代碼中建立以下循環(huán)：曲柄角度θ2從0°開始，以2°為步長，遞增至360°；對應(yīng)每個θ2值調(diào)用位置函數(shù)求解位置角。

在MATLAB命令提示符下，輸入以下代碼：>>clear；th=[045100]；rs=[0.41 0.71.2]；n=180；theta=fourbar(th，rs，n)；執(zhí)行函數(shù)文件fourbar.m，將得到一個關(guān)于機(jī)構(gòu)位置解和誤差的矩陣theta，以及從動件角位移及其誤差的曲線圖，如圖7所示。

圖5　四連桿機(jī)構(gòu)運(yùn)動學(xué)仿真結(jié)果

由圖7的誤差曲線可知，從動件角位移誤差在±1.5°范圍。此外，可以改變輸入的n值，修改步長得到更為精確的仿真曲線。通過查詢該仿真曲線或者矩陣theta(可以編寫m文件實現(xiàn))，可以獲取任意位置的位置解，為該機(jī)構(gòu)以其他位置為起點進(jìn)行仿真提供可靠的初值。

4結(jié)束語

本文依次推導(dǎo)了與四連桿機(jī)構(gòu)動力學(xué)學(xué)相關(guān)的數(shù)學(xué)模型， 并且建立其仿真模型進(jìn)行了相關(guān)分析與改進(jìn)。改進(jìn)后的動力學(xué)仿真模型，操作簡單，調(diào)試方便，通用性強(qiáng)。

圖7　從動件角位移及誤差仿真曲線

此外，本文還提供了四連桿機(jī)構(gòu)位置解的另一種求解途徑，通過執(zhí)行腳本文件fourbar.m發(fā)現(xiàn)，這種方法的結(jié)果誤差也比較小且操作簡單，而且能給仿真前提供可靠的初始值。

總之，與以往的對四連桿機(jī)構(gòu)動力學(xué)仿真的程序相比，本文設(shè)計的程序充分利用了MATLAB軟件的編程特點—m文件之間的內(nèi)部調(diào)用的便捷性，大大提高了程序的通用性與操作簡易的程度，這給方便地修改參數(shù)以適應(yīng)任意四連桿機(jī)構(gòu)的動力學(xué)仿真提供了高效的通道。

參考文獻(xiàn)：

[1]李明磊，賈育秦.基于MATLAB/Simulink的齒輪五桿機(jī)構(gòu)運(yùn)動仿真研究[J].太原科技大學(xué)學(xué)報，2009，30(1)：59-62.

[2]武娟，任家駿，吳鳳林，等.基于MATLAB軟件中SIMULINK的機(jī)構(gòu)運(yùn)動仿真[J].機(jī)械工程與自動化，2004，124(3)：4-6.

[3]徐梓斌.六桿機(jī)構(gòu)運(yùn)動學(xué)仿真的MATLAB實現(xiàn)[J].煤礦機(jī)械，2006，27(4)：617-618.

[4]孫桓，陳作模，葛文杰.機(jī)械原理[M].北京：高等教育出版社，2006.

[5]曲秀全.基于MATLAB/Simulink平面連桿機(jī)構(gòu)的動態(tài)仿真[M].哈爾濱：哈爾濱工業(yè)大學(xué)出版社，2007.

[6]GARDNER J F.機(jī)構(gòu)動態(tài)仿真：使用MATLAB和Simulink[M].周進(jìn)雄，張陵，譯.西安：西安交通大學(xué)出版社，2002.

[7]陸鳳儀，孫占營，徐格寧.基于桿組理論的連桿機(jī)構(gòu)運(yùn)動精度及運(yùn)動仿真研究[J].太原重型機(jī)械學(xué)院學(xué)報，2003，24(1)：57-61.

[8]薛定宇，陳陽泉.基于MATLAB/Simulink的系統(tǒng)仿真技術(shù)與應(yīng)用[M].北京：清華大學(xué)出版社，2002.

[9]王永超.基于MATLAB的機(jī)構(gòu)運(yùn)動學(xué)仿真方法及其比較[J].計算機(jī)仿真，2004，21(8)：81-84.

[10]王芳，張海燕.基于Simulink的連桿機(jī)構(gòu)運(yùn)動學(xué)仿真[J].機(jī)械設(shè)計與研究，2004，20(2)：35-37.

Improvement for Dynamical Simulation of Four-bar Linkage

Mechanism Based on MATLAB/Simulink

ZHANG Jian-jun，HUANG Qing-xue，MA Li-nan

(Heavy Machinery Ministry of Education Engineering Research Center，Taiyuan University of

Science and Technology，Taiyuan 030024，China)

Abstract：To improve dynamical simulation program of four-bar linkage mechanism，a general scheme for better simulation was studied.After providing the required initial parameters，one script will be run to proceed dynamic simulation program.Besides，one more new solution will be introduced to confirm the dynamic positions of mechanism.After inputting some starting angles， a function file based on Newton-Rapson algorithm was executed to draw displacement and its error diagrams.The new method not only solves the location problem of driven links，but also provides more precision initial angles when starting simulation from any other positions to make sure dynamic simulation successful and precise.

Key words：four bar linkage mechanism，dynamical simulation，improvement for generality