六桿級機構的動態(tài)仿真

1、中文摘要：運動學仿真及7R六桿機制的自由程度與Matlab動態(tài)仿真 學生姓名：施小康 班級： 指導老師：朱保利摘要：用復數(shù)推導了曲柄、6R級桿組的運動學和動力學數(shù)學模型，并將其轉化為適用于MATLAB仿真的矩陣數(shù)學模型，以該矩陣數(shù)學模型編制了相應的M函數(shù)仿真模塊，對給定的7R六桿級機構為例說明如何使用這二個仿真模塊建立MAllAB仿真模型，并對其仿真結果的正確性加以分析其主要目的是以組成機構的桿組為仿真模塊，搭建各種平面連桿機構MATLAB仿真模型，可以對各種低副機構進行動力學仿真和分析關鍵詞：運動學 動力學 機構 仿真 桿組指導老師簽名：英文摘要：Kinematic simulation &

2、 dynamic simulation of 7R six bar degree of freedom mechanism with MATLABStudent name:Shi Xiaokang Class: Supervisor:Zhu Baoli Abstract：The complex formula of displacement,velocity and acceleration of crank and 6R degree basic bar group is derived by complex.And so it is the dynamic complex formula

3、of the same mechanism.The complex formula is changed to matrix fomlula for MATIAB simulationAccording to the matrix formula，Mfunction of crank and 6Rdegree basic bar group is made up for simulation model blockMATLAB simulation modeI of 7R six bardegree mechanism is givenMain idea of the paper is to

4、build MATLAB simulation model of aU kinds lower pair mechanism by cell of basic bar groupsKey words：kinematics；mechanism；simulation；bar groupSignature of Supervisor:目 錄1 緒論(2)1.1引言(2)1.2平面連桿機構及桿組概述(2)1.3 進行桿組系統(tǒng)仿真的意義(3)1.4 仿真軟件的發(fā)展狀況與應用(3)1.5 MATLAB概述(3)2 7R六桿級機構運動學仿真(5)2.1 曲柄原動件運動學分析(5)2.2 6R級桿組運動學分析

5、(6)2.3 7R六桿級機構MATLAB仿真積分模塊初值的確定(11)2.4 7R六桿級機構運動學仿真模型及結果(16)3 7R六桿級機構動力學仿真(23)3.1 曲柄原動件動力學數(shù)學模型的建立(23)3.2 6R級桿組動力學數(shù)學模型的建立(25)3.3 需要引用的函數(shù)(30)3.4 7R六桿級機構運動學仿真模型及結果(32)4 結論(39)參考文獻(40)致謝(41)1 緒論1.1引言大學的四年生活，通過老師的講解和我自己的學習，我收獲了很多，我也深深的喜歡上了機械這個行業(yè)，對機械加工和制造方面尤為感興趣，我覺得通過自己的努力和思考來改變工藝規(guī)程來提高生產(chǎn)效率，提高經(jīng)濟效益很有成就感。我所研

6、究的課題就是給了這樣的機會我可以通過我的努力來優(yōu)化工藝規(guī)程，提高經(jīng)濟效益。此次畢業(yè)設計，是在我們學完了機械制造工藝學、工藝裝備設計等課程，進行了生產(chǎn)實習之后，進行的一個重要的實踐性環(huán)節(jié)。這要求我們把所學的工藝理論和實踐知識，在實際的工藝、夾具設計中綜合地加以運用，這有助與提高了我們分析和解決生產(chǎn)實際問題的能力，為以后從事相關的技術工作奠定的基礎。1.2平面連桿機構及桿組概述平面連桿機構是將各構件用轉動副或移動副聯(lián)接而成的平面機構。最簡單的平面連桿機構是由四個構件組成的，簡稱平面四桿機構。它的應用非常廣泛，而且是組成多桿機構的基礎。全部用回轉副組成的平面四桿機構稱為鉸鏈四桿機構。僅能在某一角度擺

7、動的連架桿，稱為搖桿。對于鉸鏈四桿機構來說，機架和連桿總是存在的，因此可按照連架桿是曲柄還是搖桿，將鉸鏈四桿機構分為三種基本型式：曲柄搖桿機構、雙曲柄機構和雙搖桿機構。在實際機械中，平面連桿機構的型式是多種多樣的，但其中絕大多數(shù)是在鉸鏈四桿機構的基礎上發(fā)展和演化而成。如曲柄滑塊機構、導桿機構等。任何機構都是由原動件、機架和從動件構成的系統(tǒng)。由于機架的自由度為零，一般每個原動件的自由度為1，且根據(jù)運動鏈成為機構的條件可知，機構的自由度與原動件為應相等，所以，從動件系統(tǒng)的自由度數(shù)必為零。機構的從動件系統(tǒng)還可以進一步分解成若干個不可再分的自由度為零的構件組合，這種組合稱為桿組。設n表示活動構件數(shù)，P

8、L表示低副個數(shù)，根據(jù)n的取值不同，村級可分為級桿組和級桿組。其中級桿組分為5種：RRR級桿組、RRP級桿組、RPR級桿組、PRP級桿組以及RPP級桿組。任何機構都可以看作是由若干個基本村級依次聯(lián)接于原動件和機架而構成的，這就是所謂機構的組成原理。通常，把由最高級別為級桿組的基本桿組構成的機構稱為級機構；把最高級為級桿組的基本桿組構成的機構稱為級機構。1.3 進行桿組系統(tǒng)仿真的意義系統(tǒng)仿真是建立在控制理論、相似理論、信息處理技術和計算機初等理論基礎之上的，以計算機和其他專用物理效應設備為工具，利用系統(tǒng)模型對真實或假設的系統(tǒng)進行試驗，并借助于專家的經(jīng)驗知識、統(tǒng)計數(shù)據(jù)和信息資料對實驗結果進行分析研究

9、，進而做出決策的一門綜合的實驗性學科。仿真技術是分析、研究各種系統(tǒng)，尤其是復雜系統(tǒng)的重要工具。隨著機械行業(yè)的迅速發(fā)展，對研究、設計的機械設備越來越復雜，用于制造各種零件的材料價格越來越昂貴，不可能每一步都采取試制再修改的方法進行設計，采用仿真的方法可以在一定程度上克服這種不足的不足，降低研究成本，提高效率。而連桿機構作為常見的傳動機構，對其進行運動學和動力學仿真，建立起基本桿組模塊的仿真模型，無疑對日后的設計大有裨益。一般機構的運動分析，使用Quik BASIC語言或Fortran語言編寫程序進行計算，其缺點“透明性”差，修改麻煩等而用MATLAB對機構進行運動仿真，利用MATLAB的simu

10、link仿真模型的數(shù)據(jù)可視化的特點，就可以很容易觀察到運動參數(shù)是如何變化的，極其簡便同時，用MATLAB建立和修改仿真模型具有方便、快捷、很容易擴展等優(yōu)點MATLAB仿真求解器提供很多解不同微分方程的方法，可以根據(jù)不同的微分方程類型選擇相應的求解方法機構的動力學分析，由已知工作阻力，求出運動副的約束反力和驅動力（或力矩），為選擇和設計軸承和零部件強度的計算及選擇原動機提供理論依據(jù)。1.4 仿真軟件的發(fā)展狀況與應用早期的計算機仿真技術大致經(jīng)歷了幾個階段：20世紀40年代模擬計算機仿真；50年代初數(shù)字仿真；60年代早期仿真語言的出現(xiàn)等。80年代出現(xiàn)的面向對象仿真技術為系統(tǒng)仿真方法注入了活力。我國早

11、在50年代就開始研究仿真技術了，當時主要用于國防領域，以模擬計算機的仿真為主。70年代初開始應用數(shù)字計算機進行仿真4。隨著數(shù)字計算機的普及，近20年以來，國際、國內(nèi)出現(xiàn)了許多專門用于計算機數(shù)字仿真的仿真語言與工具，如CSMP，ACSL， SIMNOM， MATLAB/Simulink， Matrix/System Build， CSMP-C等。1.5 MATLAB概述MATLAB是國際上仿真領域最權威、最實用的計算機工具。它是MathWork公司于1982年推出的一套高性能的數(shù)值計算和可視化數(shù)學軟件，被譽為“巨人肩上的工具”。MATLAB是一種應用于計算技術的高性能語言。它將計算，可視化和編程

12、結合在一個易于使用的環(huán)境中，此而將問題解決方案表示成我們所熟悉的數(shù)學符號，其典型的使用包括:.數(shù)學計算.運算法則的推導.模型仿真和還原.數(shù)據(jù)分析，采集及可視化.科技和工程制圖.開發(fā)軟件，包括圖形用戶界面的建立MATLAB是一個交互式系統(tǒng)，它的基本數(shù)據(jù)元素是矩陣，且不需要指定大小。通過它可以解決很多技術計算問題，尤其是帶有矩陣和矢量公式推導的問題，有時還能寫入非交互式語言如C和Fortran等。MATLAB的名字象征著矩陣庫。它最初被開發(fā)出來是為了方便訪問由LINPACK和EISPAK開發(fā)的矩陣軟件，其代表著藝術級的矩陣計算軟件。 MATLAB在擁有很多用戶的同時經(jīng)歷了許多年的發(fā)展時期。在大學環(huán)

13、境中，它作為介紹性的教育工具，以及在進階課程中應用于數(shù)學，工程和科學。在工業(yè)上它是用于高生產(chǎn)力研究，開發(fā)，分析的工具之一。Simulink概述 Simulink是用于仿真建模及分析動態(tài)系統(tǒng)的一組程序包，它支持線形和非線性系統(tǒng)，能在連續(xù)時間，離散時間或兩者的復合情況下建模。系統(tǒng)也能采用復合速率，也就是用不同的部分用不同的速率來采樣和更新。Simulink提供一個圖形化用戶界面用于建模，用鼠標拖拉塊狀圖表即可完成建模。在此界面下能像用鉛筆在紙上一樣畫模型。相對于以前的仿真需要用語言和程序來表明不同的方程式而言有了極大的進步。Simulink擁有全面的庫，如接收器，信號源，線形及非線形組塊和連接器。

14、同時也能自己定義和建立自己的塊。模塊有等級之分，因此可以由頂層往下的步驟也可以選擇從底層往上建模?？梢栽诟邔由辖y(tǒng)觀系統(tǒng)，然后雙擊模塊來觀看下一層的模型細節(jié)。這種途徑可以深入了解模型的組織和模塊之間的相互作用。2 7R六桿級機構運動學仿真2.1 曲柄原動件運動學分析2.1.1 曲柄原動件運動學數(shù)學模型的建立圖1 曲柄的復數(shù)坐標系如圖1所示，在復數(shù)坐標系中，曲柄AB復向量的模rj為常數(shù)、幅角j為變量，通過轉動副A與機架連接，轉動副A的復向量的模ri為常量、幅角i為常量，曲柄AB端點B的位移、速度和加速度的推導如下： （2.1）將方程2.1兩邊對時間t求兩次導數(shù)得：（2.2）由式2.2寫成矩陣形式有

15、 （2.3）2.1.2 曲柄MATLAB運動學仿真模塊M函數(shù)根據(jù)式(2.3)編寫曲柄原動件MATLAB的M函數(shù)如下：function y=crank(x)%x(1)=rj 曲柄桿長%x(2)=thetaj曲柄與水平方向夾角%x(3)=dthetaj曲柄角速度%x(4)=ddthetaj曲柄角加速度%y(1)=ReddB轉動副B加速度實軸分量%y(2)=ImddB轉動副B加速度虛軸分量ddB=x(1)*x(4)*cos(x(2)+pi/2)+x(1)*x(3)2*cos(x(2)+pi); x(1)*x(4)*sin(x(2)+pi/2)+x(1)*x(3)2*sin(x(2)+pi); y=d

16、dB;此函數(shù)模塊用于計算轉動副B的加速度的水平分量和垂直分量輸入?yún)?shù)為曲柄的長度、角位移、角速度和角加速度；輸出參數(shù)為曲柄端部（轉動副B）的加速度的水平分量和垂直分量。2.2 6R級桿組運動學分析2.2.1 6R級桿組運動學數(shù)學模型的建立圖2 6R級桿組的位置參數(shù)如圖2所示，在復數(shù)坐標系中，由3個外轉動副(B,C,D)和3個內(nèi)轉動副(E,F,G),4個構件(BE,CF,DG和EFG)構成1個6R級桿組，構件BE,CF,DG的長度分別為ri，rj，rk構件EFG的3個邊為e，f，g，方向如圖所示，規(guī)定所有復向量與實軸正方向逆時針夾角為，并用相應的下標來區(qū)別，用 B,C,D,E,F和G分別表示該轉

17、動副的復數(shù)坐標，則各個構件的運動參數(shù)推導如下： （2.4）（2.5）將式（2.4）代入式（2.5）并整理得： （2.6）將式(2.4)至式(2.6)三式合并成矩陣得： （2.7）將式(2.7)展開整理得：（2.8）點E，F(xiàn)，G的加速度分別為 （2.9）2.2.2 6R級桿組MATLAB運動學仿真模塊M函數(shù)根據(jù)式(2.9)編寫6R級桿組MATLAB的M函數(shù)如下：function y=R6ki(x)%x(1)=riBE桿長%x(2)=rjCF桿長%x(3)=rkDG桿長%x(4)=eFG桿長%x(5)=fGE桿長%x(6)=gEF桿長%x(7)=theta-iBE桿與水平方向夾角%x(8)=the

18、ta-jDF桿與水平方向夾角%x(9)=theta-kDG桿與水平方向夾角%x(10)=theta-eFG桿與水平方向夾角%x(11)=theta-fGE桿與水平方向夾角%x(12)=theta-gEF桿與水平方向夾角%x(13)=dtheta-iBE桿角速度%x(14)=dtheta-jEF桿角速度%x(15)=dtheta-kDG桿角速度%x(16)=dtheta-eGE桿角速度%x(17)=dtheta-fGE桿角速度%x(18)=dtheta-gEF桿角速度%x(19)=ReddB轉動副B加速度實軸分量%x(20)=ImddB轉動副B加速度虛軸分量%x(21)=ReddC轉動副C加速度

19、實軸分量%x(22)=ImddC轉動副C加速度虛軸分量%x(23)=ReddD轉動副D加速度實軸分量%x(24)=ImddD轉動副D加速度虛軸分量%y(1)=ddtheta-iBE桿角加速度%y(2)=ddtheta-jCF桿角加速度%y(3)=ddtheta-kDG桿角加速度%y(4)=ddtheta-eFG桿角加速度%y(5)=ddtheta-fGE桿角加速度%y(6)=ddtheta-fEF桿角加速度%y(7)=ReddE轉動副E加速度實軸分量%y(8)=ImddE轉動副E加速度虛軸分量%y(9)=ReddF轉動副F加速度實軸分量%y(10)=ImddF轉動副F加速度虛軸分量%y(11)

20、=ReddG轉動副G加速度實軸分量%y(12)=ImddG轉動副G加速度虛軸分量a=-x(1)*cos(x(7)+pi/2) x(2)*cos(x(8)+pi/2) 0 0 0 -x(6)*cos(x(12)+pi/2); -x(1)*sin(x(7)+pi/2) x(2)*sin(x(8)+pi/2) 0 0 0 -x(6)*sin(x(12)+pi/2); 0 -x(2)*cos(x(8)+pi/2) x(3)*cos(x(9)+pi/2) -x(4)*cos(x(10)+pi/2) 0 0; 0 -x(2)*sin(x(8)+pi/2) x(3)*sin(x(9)+pi/2) -x(4)

21、*sin(x(10)+pi/2) 0 0; x(1)*cos(x(7)+pi/2) 0 -x(3)*cos(x(9)+pi/2) 0 -x(5)*cos(x(11)+pi/2) 0; x(1)*sin(x(7)+pi/2) 0 -x(3)*sin(x(9)+pi/2) 0 -x(5)*sin(x(11)+pi/2) 0;c=-x(1)*cos(x(7)+pi) x(2)*cos(x(8)+pi) 0 0 0 -x(6)*cos(x(12)+pi); -x(2)*sin(x(7)+pi) x(2)*sin(x(8)+pi) 0 0 0 -x(6)*sin(x(12)+pi); 0 -x(2)*c

22、os(x(8)+pi) x(3)*cos(x(9)+pi) -x(4)*cos(x(10)+pi) 0 0; 0 -x(2)*sin(x(8)+pi) x(3)*sin(x(9)+pi) -x(4)*sin(x(10)+pi) 0 0; x(1)*cos(x(7)+pi) 0 -x(3)*cos(x(9)+pi) 0 -x(5)*cos(x(11)+pi) 0; x(1)*sin(x(7)+pi) 0 -x(3)*sin(x(9)+pi) 0 -x(5)*sin(x(11)+pi) 0;b1=c*x(13)2;x(14)2;x(15)2;x(16)2;x(17)2;x(18)2;b2=x(19

23、)-x(21);x(20)-x(22);x(21)-x(23);x(22)-x(24);x(23)-x(19);x(24)-x(20);b=b1+b2;ddtheta=inv(a)*b;y(1)=ddtheta(1);y(2)=ddtheta(2);y(3)=ddtheta(3);y(4)=ddtheta(4);y(5)=ddtheta(5);y(6)=ddtheta(6);y(7)=x(19)+x(1)*ddtheta(1)*cos(x(7)+pi/2)+x(1)*x(13)2*cos(x(7)+pi);y(8)=x(20)+x(1)*ddtheta(1)*sin(x(7)+pi/2)+x(

24、1)*x(13)2*sin(x(7)+pi);y(9)=x(21)+x(2)*ddtheta(2)*cos(x(8)+pi/2)+x(2)*x(14)2*cos(x(8)+pi);y(10)=x(22)+x(2)*ddtheta(2)*sin(x(8)+pi/2)+x(2)*x(14)2*sin(x(8)+pi);y(11)=x(23)+x(3)*ddtheta(3)*cos(x(9)+pi/2)+x(3)*x(15)2*cos(x(9)+pi);y(12)=x(24)+x(3)*ddtheta(3)*sin(x(9)+pi/2)+x(3)*x(15)2*sin(x(9)+pi);這個模塊用于

25、求級桿組中各桿的加速度的水平及垂直分量。輸入?yún)?shù)為構件2、構件3的角位移和角速度，構件2、構件3和構件4的桿長，構件5的3個邊長，構件2、構件3、構件4的角位移和角速度，構件5的3個邊向量的角位移和3個轉動副B,C,D的加速度；輸出參數(shù)為構件2、構件3和構件4的角加速度，構件5的3個邊向量的角加速度和轉動副E,F,G的加速度。2.3 7R六桿級機構MATLAB仿真積分模塊初值的確定2.3.1 運用牛頓辛普森法進行角位移分析 圖3 7R六桿級機構圖3所示是由原動件（曲柄1）和一個6R級桿組所組成的7R六桿級機構，復數(shù)向量坐標亦如圖所示，各構件的尺寸為r1=120mm，r2=400mm，r3=30

26、0mm，r4=300mm，ReD=250mm，F(xiàn)G=450mm，GE=180mm，EF=350mm，ImD=350mm，ReC=700mm，ImC=350mm，構件1以等角速度10rad/s逆時針方向回轉，試求構件2和構件3的位移、速度和加速度。由圖3可列以下三個方程：AB+BEAD-DG EG=0，即r1+r2-AD-r4-f=0 (2.10)AB+BE+EF-AC-CF=0即r1+r2+g-AC-r3=0 (2.11)EF+FG+GE=0,即g+e+f=0 (2.12)由復向量坐標，可寫出式(2.10)、式(2.11)及式(2.12)的角位移方程為： （2.13） (2.14) (2.15

27、)將式(2.13)、式(2.14)、式(2.15)展開，整理得: (2.16)由式(2.16)求出雅可比矩陣為： (2.17)根據(jù)式(2.16)、式(2.17)，由牛頓辛普森求解方法得編制M函數(shù)如下：function y=r6posi(x)%x(1)=theta-1桿1與水平方向夾角%x(2)=theta-2 桿2與水平方向夾角（估計量）%x(3)=theta-3 桿3與水平方向夾角（估計量）%x(4)=theta-4 桿4與水平方向夾角（估計量）%x(5)=theta-e 桿e與水平方向夾角（估計量）%x(6)=theta-f 桿f與水平方向夾角（估計量）%x(7)=theta-g 桿g與水

28、平方向夾角（估計量）%x(8)=theta-ACAC與水平方向夾角%x(9)=theta-ADAD與水平方向夾角%x(10)=r1桿1長度%x(11)=r2桿2長度%x(12)=r3桿3長度%x(13)=r4桿4長度%x(14)=e桿e長度%x(15)=f桿f長度%x(16)=g桿g長度%x(17)=AC桿AC長度%x(18)=AD桿AD長度%y(1)=theta-22桿與水平方向夾角%y(2)=theta-33桿與水平方向夾角%y(3)=theta-44桿與水平方向夾角%y(4)=theta-ee桿與水平方向夾角%y(5)=theta-ff桿與水平方向夾角%y(6)=theta-gg桿與水平

29、方向夾角%theta2=x(2);theta3=x(3);theta4=x(4);theta5=x(5);theta6=x(6);theta7=x(7);%epsilon=1.0E-6;%f=x(10)*cos(x(1)+x(11)*cos(theta2)-x(18)*cos(x(9)-x(13)*cos(theta4)-x(15)*cos(theta6);x(10)*sin(x(1)+x(11)*sin(theta2)-x(18)*sin(x(9)-x(13)*sin(theta4)-x(15)*sin(theta6);x(10)*cos(x(1)+x(11)*cos(theta2)+x(1

30、6)*cos(theta7)-x(17)*cos(x(8)-x(12)*cos(theta3);x(10)*sin(x(1)+x(11)*sin(theta2)+x(16)*sin(theta7)-x(17)*sin(x(8)-x(12)*sin(theta3);x(14)*cos(theta5)+x(15)*cos(theta6)+x(16)*cos(theta7);x(14)*sin(theta5)+x(15)*sin(theta6)+x(16)*sin(theta7);%while norm(f) epsilon J= -x(11)*sin(theta2) 0 x(13)*sin(the

31、ta4) 0 x(15)*sin(theta6) 0; x(11)*cos(theta2) 0 -x(13)*cos(theta4) 0 -x(15)*cos(theta6) 0; -x(11)*sin(theta2) x(12)*sin(theta3) 0 0 0 -x(16)*sin(theta7); x(11)*cos(theta2) -x(12)*cos(theta3) 0 0 0 x(16)*cos(theta7); 0 0 0 -x(14)*sin(theta5) -x(15)*sin(theta6) -x(16)*sin(theta7); 0 0 0 x(14)*cos(thet

32、a5) x(15)*cos(theta6) x(16)*cos(theta7); dth=inv(J)*(-1.0*f); theta2=theta2+dth(1); theta3=theta3+dth(2); theta4=theta4+dth(3); theta5=theta5+dth(4); theta6=theta6+dth(5); theta7=theta7+dth(6);f=x(10)*cos(x(1)+x(11)*cos(theta2)-x(18)*cos(x(9)-x(13)*cos(theta4)-x(15)*cos(theta6);x(10)*sin(x(1)+x(11)*

33、sin(theta2)-x(18)*sin(x(9)-x(13)*sin(theta4)-x(15)*sin(theta6);x(10)*cos(x(1)+x(11)*cos(theta2)+x(16)*cos(theta7)-x(17)*cos(x(8)-x(12)*cos(theta3);x(10)*sin(x(1)+x(11)*sin(theta2)+x(16)*sin(theta7)-x(17)*sin(x(8)-x(12)*sin(theta3);x(14)*cos(theta5)+x(15)*cos(theta6)+x(16)*cos(theta7);x(14)*sin(theta

34、5)+x(15)*sin(theta6)+x(16)*sin(theta7);norm(f)end;y(1)=theta2;y(2)=theta3;y(3)=theta4;y(4)=theta5;y(5)=theta6;y(6)=theta7;此函數(shù)模塊的輸入量為桿1與水平方向的夾角、其余各桿與水平方向夾角的估計值以及桿的桿長及AC、AD兩個收入?yún)⒘康臈U長和角度的估計值，輸出參量為2，3，4，E,F,G桿與水平方向的夾角。7R六桿級機構在圖1.3.1所示位置，估計構件2，3，4，e，f，g的角位移為2=5.9341rad，3=1.9199rad，4=1.9199rad，e=0，f=2.3562

35、rad，g=0.3491。輸入?yún)?shù)x=40*pi/180 340*pi/180 110*pi/180 110*pi/180 0 135*pi/180 20*pi/180 atan(350/700) atan(350/250) 120 400 300 300 450 180 350 sqrt(3502+7002) sqrt(3502+2502)，解得2=-0.3725rad，3=-1.2735rad，4=-1.2735rad，e=3.1416rad，f=-0.8040rad，g=0.3794rad。2.3.2 運用牛頓辛普森法進行角速度分析對式(2.10)-式(2.12)求導并展開成矩陣形式為：

36、 (2.18)根據(jù)式(2.18)編寫M函數(shù)如下：function y=r6vel(x)%x(1)=theta-1桿1與水平方向夾角%x(2)=theta-2桿2與水平方向夾角%x(3)=theta-3桿3與水平方向夾角%x(4)=theta-4桿4與水平方向夾角%x(5)=theta-e桿e與水平方向夾角%x(6)=theta-f桿f與水平方向夾角%x(7)=theta-g桿g與水平方向夾角%x(8)=dtheta-1桿1角速度%x(9)=r1桿1長度%x(10)=r2桿2長度%x(11)=r3桿3長度%x(12)=r4桿4長度%x(13)=e桿e長度%x(14)=f桿f長度%x(15)=g桿

37、g長度%y(1)=dtheta-2桿2角加速度%y(2)=dtheta-3桿3角加速度%y(3)=dtheta-4桿4角加速度%y(4)=dtheta-e桿e角加速度%y(5)=dtheta-f桿f角加速度%y(6)=dtheta-g桿g角加速度A= -x(10)*sin(x(2) 0 x(12)*sin(x(4) 0 x(14)*sin(x(6) 0; x(10)*cos(x(2) 0 x(12)*cos(x(4) 0 x(14)*sin(x(6) 0; -x(10)*sin(x(2) x(11)*sin(x(3) 0 0 0 -x(15)*sin(x(7); x(10)*cos(x(2)

38、-x(11)*cos(x(3) 0 0 0 x(15)*cos(x(7); 0 0 0 -x(13)*sin(x(5) -x(14)*sin(x(6) -x(15)*sin(x(7); 0 0 0 x(13)*cos(x(5) x(15)*cos(x(6) x(15)*cos(x(7);B=x(10)*sin(x(1);-x(10)*cos(x(1);x(10)*sin(x(1);-x(10)*cos(x(1);0;0*x(8);y=inv(A)*B;此函數(shù)模塊輸入量為各桿與水平方向的角度以及桿1的角速度，輸出參數(shù)為2，3，4，E,F,G桿的角加速度。圖3所示機構，由位移分析計算出的各桿角度和

39、曲柄1的角速度為10rad/s，則輸入?yún)?shù)為：x= 40*pi/180 -0.3725 -1.2735 -1.2735 3.1416 -0.8040 0.3794 10 120 400 300 300 450 180 350，代入上面的M文件，求得2，3，4，e，f，g桿的角速度依次分別為 -3.49rad/s，-4.5298rad/s，-4.5298rad/s，0rad/s，0rad/s，0rad/s。2.4 7R六桿級機構運動學仿真模型及結果2.4.1 7R六桿級機構圖3所示是由原動件（曲柄1）和1個6R級桿組所組成的7R六桿級機構，復數(shù)向量坐標如圖1.4.1所示，各構件的尺寸為的r1=1

40、20mm，r2=400mm，r3=300mm，r4=300mm，ReD=250mm，F(xiàn)G=450mm，GE=180mm，EF=350mm，ImD=350mm，ReC=700mm，ImC=350mm，構件1以等角速度10rad/s逆時針方向回轉，試求構件2和構件3的位移、速度、加速度。 2.4.2 7R六桿級機構MATLAB運動學仿真模型7R六桿級機構MATLAB運動學仿真模型如圖2.4.2所示，在圖2.4.2中各積分模塊的初值是以曲柄1的幅角為0.72rad和角速度等于10rad/s逆時針方向回轉時，相應各個構件的位移、速度的瞬時值，2個MATLAB函數(shù)模塊分別為crank.m和r6ki.m，

41、其中crank.m函數(shù)模塊的輸入?yún)?shù)為曲柄的長度、角位移、角速度和角加速度；輸出參數(shù)為曲柄端部（轉動副B）的加速度的水平分量和垂直分量;r6ki.m函數(shù)模塊的輸入?yún)?shù)為構件2、構件3的角位移和角速度，構件2、構件3和構件4的桿長，構件5的3個邊長，構件2、構件3、構件4的角位移和角速度，構件5的3個邊向量的角位移和3個轉動副B,C,D的加速度；輸出參數(shù)為構件2、構件3和構件4的角加速度，構件5的3個邊向量的角加速度和轉動副E,F,G的加速度。每個數(shù)據(jù)線上標注了相應變量，常量模塊放置了各個構件的尺寸，長度分別為m，角度單位為rad。設置仿真時間為1s，仿真結果輸出到工作空間變量simout中，輸

42、出格式為array,求解器選用ode45,步長選用變步長。2.4.3 7R六桿級機構MATLAB運動學仿真結果由于曲柄轉速為10rad/s,因此每轉動1周的時間是0.628s,用繪圖命令plot(tout,simout(:,1)，plot(tout,simout(:,2)，plot(tout,simout(:,5)，plot(tout,simout(:,6)，plot(tout,simout(:,9)，plot(tout,simout(:,10)，plot(tout,simout(:,3)，plot(tout,simout(:,4)繪制出構件2和構件3的位移、速度、加速度，構件4的位移、向量e

43、的位移，如圖2.4.3所示。從該圖中可以看出這些參數(shù)也都是周期變化的。從圖3所給出的各個構件的尺寸可以看出，由構件3，4，5，6構成平行四邊形，因此構件5就作平動，復數(shù)向量e,f,g的角位移為常量、角速度為零。從圖6（g）中看出復數(shù)向量e仿真結果確實如此 。同時構件3，4的角位移、角速度、角加速度應該相同，比較圖6（b）和（h）曲線也相同，這兩組也說明了所推導的6R級桿組的運動學仿真公式及相應的M函數(shù)和仿真模型是正確的。 圖6（b） 構件3的角位移 (縱坐標表示角位移的大小，單位為rad；橫坐標表示時間，單位為s。) 圖6（a） 構件2的角位移(縱坐標表示角位移的大小，單位為rad；橫坐標表示

44、時間，單位為s。) 圖6（c） 構件3的角速度(縱坐標表示角速度的大小，單位為rad/s；橫坐標表示時間，單位為s。) 圖6（d） 構件2的角速度(縱坐標表示角位移的大小，單位為rad/s；橫坐標表示時間，單位為s。) 圖6（e） 構件2的角加速度(縱坐標表示角加速度的大小，單位為rad/s2；橫坐標表示時間，單位為s。) 圖6（f） 構件3的角加速度(縱坐標表示角加速度的大小，單位為rad/s2；橫坐標表示時間，單位為s。) 圖6（g） 向量e的角位移(縱坐標表示角位移的大小，單位為rad；橫坐標表示時間，單位為s。) 圖6（h） 構件4的角位移(縱坐標表示角位移的大小，單位為rad；橫坐標

45、表示時間，單位為s。)3 7R六桿級機構動力學仿真3.1 曲柄原動件動力學數(shù)學模型的建立3.1.1 曲柄原動件動力學數(shù)學分析 圖7 曲柄的受力模型如圖7所示，已知曲柄AB向量的模ri為常數(shù)，幅角i為變量，質(zhì)心到轉動副A的距離為rci,質(zhì)量為mi，繞質(zhì)心轉動慣量為Ji，作用于質(zhì)心上的外力為Fxi和Fyi、外力矩為Mi，曲柄與機架聯(lián)接，轉動副A的約束反力為Rxa和Rya，驅動力矩為Ml。 由理論力學可得： (3.1) (3.2) (3.3) 由運動學知識可推得： (3.4) (3.5)將式（3.4）、式（3.5）代入式（3.1）、式（3.2），并與式（3.3）合并得: (3.6)3.1.2 曲柄M

46、ATLAB動力學仿真模塊M函數(shù)根據(jù)式(3.6)編寫曲柄原動件MATLAB的M函數(shù)如下：function y=crankdy(x)%x(1)=theta-i曲柄與水平方向夾角%x(2)=dtheta-i曲柄角速度%x(3)=ddtheta-i曲柄角加速度%x(4)=RxB轉動副B約束反力水平分量%x(5)=RyB轉動副B約束反力垂直分量%y(1)=RxA轉動副A約束反力水平分量%y(2)=RyA轉動副A約束反力垂直分量%y(3)=M1轉動副A的驅動力矩%g=9.8ri=0.4;rci=0.2;mi=1.2;Ji=0.016;Fxi=0;Fyi=0;Mi=0;ReddA=0;ImddA=0;y(1

47、)=mi*ReddA+mi*rci*x(3)*cos(x(1)+pi/2)+mi*rci*x(2)2*cos(x(1)+pi)-Fxi+x(4);y(2)=mi*ImddA+mi*rci*x(3)*sin(x(1)+pi/2)+mi*rci*x(2)2*sin(x(1)+pi)-Fyi+x(5)+mi*g;y(3)=Ji*x(3)-y(1)*rci*sin(x(1)+y(2)*rci*cos(x(1)-x(4)*(ri-rci)*sin(x(1)+x(5)*(ri-rci)*cos(x(1)-Mi;此函數(shù)模塊轉用于求動副A的約束反力和曲柄上作用的驅動力矩。輸入?yún)?shù)是曲柄原動件的角位移、角速度和

48、角加速度 以及轉動副B的反作用力；輸出參數(shù)是轉動副A的約束反力和曲柄上作用的驅動力矩。3.2 6R級桿組動力學數(shù)學模型的建立3.2.1 6R級桿組動力學數(shù)學分析由3個外轉動副（B,C,D）和3個內(nèi)轉動副（E,F,G）組成的6R級桿組如圖3.2.1(a)所示，其中具有2個轉動副的3個構件：BE構件、FC構件和GD構件的受力情況如圖3.2.1(b)、(c)、(d)所示。除不受驅動力矩外，其受力分析同曲柄的受力分析是一樣的，直接給出公式為： (3.7) (3.8) 具有3個轉動副EFG構件受力如圖8(e)所示，可推導出： (3.9) (3.10) (3.11) 圖8（d）GD桿件受力圖 圖8（c）F

49、C桿件受力圖 圖8 (b)BE桿件受力圖 圖8（a）6R級桿組 圖8（e）EFG桿件受力圖 圖8 各桿受力圖將式(3.7)與式（3.11）整理并寫成矩陣形式為： (3.12) 3.2.2 6R級桿組MATLAB動力學仿真模塊M函數(shù)根據(jù)式(3.12)編寫6R級桿組MATLAB的M函數(shù)如下：function y=R6dy(x)%x(1)=theta-iBE桿與水平方向夾角%x(2)=theta-jCF桿與水平方向夾角%x(3)=theta-kDG桿與水平方向夾角%x(4)=theta-eFG桿與水平方向夾角%x(5)=theta-fGE桿與水平方向夾角%x(6)=theta-gEF桿與水平方向夾角

50、%x(7)=dtheta-iBE桿角速度%x(8)=dtheta-jCF桿角速度%x(9)=dtheta-kDG桿角速度%x(10)=dtheta-gFG桿角速度%x(11)=ddtheta-iBE桿角加速度%x(12)=ddtheta-jCF桿角加速度%x(13)=ddtheta-kDG桿角加速度%x(14)=ddtheta-gFG桿角加速度%x(15)=ReddB轉動副B加速度實軸分量%x(16)=ImddB轉動副B加速度虛軸分量%x(17)=ReddE轉動副E加速度實軸分量%x(18)=ImddE轉動副E加速度虛軸分量%y(1)=RxB轉動副B約束反力水平分量%y(2)=RyB轉動副B約束反力垂直分量%y(3)=RxC轉動副C約束反力水平分量%y(4)=RyC轉動副C約束反力垂直分量%y(5)=RxD轉動副D約束反力水平分量%y(6)=RyD轉動副D約束反力垂直分量%y(7)=RxE轉動副E約束反力水平分量%y(8)=RyE轉動副E約束反力垂直分量%y(9)=RxF轉動副F約束反力水平分量%y(10)=RyF轉動副F約束反力垂直分量%y(11)=RxG轉動副