1712_RPP平面连杆机构的动态仿真
1712_RPP平面连杆机构的动态仿真,_rpp,平面,连杆机构,动态,仿真
南昌航空大学科技学院学士学位1RPP 平面连杆机构的运动仿真1 引言在大学四年学习,通过老师的讲解和自己的学习,收获了很多,深深的喜欢上了机械这个行业,对机械设计方面很是喜欢,我所研究的课题就是通过 MATLAB 做平面连杆 RPP 仿真,进一步加深了对平面连杆的研究。此次毕业设计,要求我对MATLAB 软件有一定的认识,在加上我学完了机械原理,进行了生产实习之后,进行的一个重要的实践性环节。这要求我们把所学的机械原理和实践知识相结合,在实际的设计中综合地加以运用,这有助与提高了我们分析和解决实际问题的能力,为以后从事相关的技术工作奠定的基础。1.1 平面连杆机构概述平面连杆机构是将各构件用转动副或移动副联接而成的平面机构。最简单的平面连杆机构是由四个构件组成的,简称平面四杆机构。它的应用非常广泛,而且是组成多杆机构的基础。全部用回转副组成的平面四杆机构称为铰链四杆机构。仅能在某一角度摆动的连架杆,称为摇杆。对于铰链四杆机构来说,机架和连杆总是存在的,因此可按照连架杆是曲柄还是摇杆,将铰链四杆机构分为三种基本型式:曲柄摇杆机构、双曲柄机构和双摇杆机构。在实际机械中,平面连杆机构的型式是多种多样的,但其中绝大多数是在铰链四杆机构的基础上发展和演化而成。如曲柄滑块机构、导杆机构等。1.2 杆组机构的从动件系统一般还可以进一步分解成若干个不可再分的自由度为零的构件组合,这种组合称为基本杆组,简称杆组。对于只含低副的平面连杆机构,若杆组中有 N 个活动构件、 个低副,因杆组自由LP度为零,故有: 023Lpn为保证 n 和 均为整数,n 只能取偶数。根据 n 的取值不同,杆组可以分为以下几LP种情况。(1)n=2, 的 II 级杆组3LII 级杆组为最简单,也是应用最多的基本杆组。根据 3 个运动副(转动副用 R 表示,运动副用 P 表示)的不同组合,II 级杆组分为 5 种,RRRII 级杆组、 RRPII 级杆组、 RPR II 级杆组、PRPII 级杆组和 RPPII 级杆组 。南昌航空大学科技学院学士学位2(2)n=4, 的 III 级杆组6LPIII 级杆组特别是 III 级以上杆组早实际应用中较少,故在这里不再介绍。1.3 机构的组成原理任何机构都可以看做是由若干个基本杆组依次联接于原动件和机构而构成的,者就是所谓机构的组成原理。把由最高级别为 II 级杆组的基本杆组构成的机构称为II 级机构,把最高级别为 III 级杆组的基本杆组构成的机构称为 III 级机构。其余类推。2 曲柄原动件、RPP 平面连杆运动学数学模型的建立2.1 曲柄原动件运动学分析2.11 曲柄原动件运动学数学模型的建立如图 2.1 所示,在复数坐标系中,曲柄 AB 复向量的模 rj 为常数、幅角 j 为变量,通过转动副 A 与机架连接,转动副 A 的复向量的模 ri 为常量、幅角 i 为常量,曲柄 AB 端点 B 的位移、速度和加速度的推导如下:图 2.1 曲柄的复数坐标系(1.1)jijijjji errABe,将方程(1.1)两边对时间 t 求两次导数得:(1.2))2/()2/( jj jje由式(1.2)写成矩阵形式有:(1.3) )sin()2/sin(cocoImR2jjjj jjjj rrB2.12 曲柄 MATLAB 运动学仿真模块 M 函数根据式编写曲柄原动件 MATLAB 的 M 函数如下:南昌航空大学科技学院学士学位3function y=crank(x)%Function to compute the acceleration of crank%Input parameters%x(1)=r1 (r1 的杆长)%x(2)=theta-1 (r1 的角位移)%x(3)=dtheta-1 (r1 的角速度)%x(4)=ddtheta-1 (r1 的角加速度)%Output parameters%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=ddB;各构件的初值为:r1=0.4, theta-1=0, dtheta-1=10,ddtheta-1=0 。2.2 RPP 四杆运动学分析2.21 RPP 四杆运动学数学模型的建立如图 2.2 所示,在复数坐标系中,由 1 个转动副(B)、2 个移动副(C,D)和 2 个滑块(C,D)组成 RPPII 级杆组,滑块 C 的滑动方向与滑块 D 的滑动方向的夹角 j 为常量,滑块 D 的幅度 j 也为常量,滑块 C 相对滑块 D 位移 S i 为变量,滑块 D 相对固定点 K 的位移也为变量,则滑块 C 相对滑块 D 的加速度和滑块 D 相对固定点 K 的加速度推导如下页图:C=B=K+ (1.4)整理(1.4)得: (1.5)KBesjijjij )(南昌航空大学科技学院学士学位4图 2.2 RPPII 级杆组的位置参数式(1.5)对 t 求导并整理得到:(1.6)KBSeejijjji )(式(1.6)对时间 t 求导并整理得:(1.7)BSeejijjji )(由(1.7)写成矩阵形式有(1.8)KSjijjjiji sincosncoRPPII 级杆组 MATLAB 运动学仿真 M 函数function y =RPPki(x)% compute veloccity of RPP bar group% Input parameters% x(1) = theta - i (构件 2 的移动方向)% x(2) = theta - j (构件 3 的移动方向)南昌航空大学科技学院学士学位5% x(3) = ReddB (转动副 B 的水平加速度)% x(4) = ImddB (转动副 B 的垂直加速度)% x(5) = ReddK (转动副 3 的水平加速度)% x(6) = ImddK (转动副 B 的垂直加速度)% Output parameters% y(1) = dds - 2 (构件 2 的加速度)% y(2) = dds - 3 (构件 3 的加速度)%a = cos(x(1)+x(2) cos(x(2);sin(x(1)+x(2) sin(x(2);b = x(3)-x(5);x(4)-x(6);ddsij = inv(a) * b; y(1) = ddsij(1);y(2) = ddsij(2);各构件的初值为:theta-i=1.5708,theta-j=0,ReddK=0,ImddK=0 。2.3 RPP 四杆机构 MATLAB 运动仿真 1、 如图 2.3 所示为 RPP 四杆机构,它由原动件(曲柄 1)和一个 RPP 杆组构成。构件的尺寸为 r1=400mm,转动副 A 到移动副 DE 的距离 r4=800mm,复数向量坐标见图上,构件 1 以等角速度 10rad/s 逆时针方向回转,试求构件 2 和构件 3 的速度和加速度?图 2.3 RPP 四杆机构2.31RPP 四杆机构 MATLAB 仿真模型图中的各积分模块的初值是以曲柄 1 的幅度为 0和角速度等于 10rad/s 逆时针方向回转,相应各个构件的位移和速度的瞬时值。theta-1=0;dtheta-1=10;南昌航空大学科技学院学士学位6theta-i=1.5708;s-2=0;ds-2=0.8;s-3=0.89;ds-3=0.4。二个 MATLAB 函数模块分别为 crank.m 和 rppki.m图 2.4 RPP 四杆机构 Simunlink南昌航空大学科技学院学士学位72.32 用 MATLAB 实现牛顿-辛普森求解方法图 2.5 RPPII 级杆组的位置参数如图 2.5 所示,以复数形式表示为 )(jijjjessAKBC整理上式得: AKBejijj )(按欧拉公式展开得: AKBjSjS jijiijj )sn()cos(sinco 上式展开整理得: 0sinsin)sin(si co43 AKBSf jijjj将上式求出雅克比矩阵为: jiijj jiijjjij jii SJ coscsicosfunction y = RPPposi(x)% Script used to implement Newton - Raphson mechod for% solving nonlinear position of RPP bar group % Input parameters% x(1) = theta-1 guess value (构件 1 的角位移)% x(2) = theta-2 guess value (构件 2 的角位移)% x(3) = theta-3 guess value (构件 3 的角位移)南昌航空大学科技学院学士学位8% x(4) = theta-4 guess value (构件 4 的角位移)% x(5) = s1 guess value (构件 1 的相对位移)% x(6) = s2 ( 构件 2 的相对位移)% x(7) = s3 guess value (构件 3 的相对位移)% x(8) = s4 guess value (构件 4 的相对位移)% Output parameters% y(1) = theta - 1 (构件 1 的角位移)% y(2) = theta - 2 (构件 2 的角位移)%theta1 = x(1);theta2 = x(2);%epsilon = 1.0E-6;%f=x(6)*cos(x(2)+x(5)*cos(x(1)+x(2)+x(7)*cos(x(3)-x(8)*cos(x(4);x(6)*sin(x(2)+x(5)*sin(x(1)+x(2)+x(7)*sin(x(3)-x(8)*sin(x(4);%while norm(f)epsilonJ=-x(5)*sin(theta1+theta2) -x(6)*sin(theta2)-x(5)*sin(theta1+theta2);x(5)*cos(theta1+theta2) x(6)*cos(theta2)+x(5)*cos(theta1+theta2);dth = inv(J)*(-1.0*f);theta1=theta1+dth(1);theta2=theta2+dth(2);f=x(6)*cos(theta2)+x(5)*cos(theta1+theta2)+x(7)*cos(x(3)-x(8)*cos(x(4);x(6)*sin(theta2)+x(5)*sin(theta1+theta2)+x(7)*sin(x(3)-x(8)*sin(x(4);norm(f);end南昌航空大学科技学院学士学位9y(1)=theta1;y(2)=theta2;估计杆 1 和杆 2 的角位移为 ,则输入参数 x=90*pi/180 0 0 0,92163.435*pi/180 0.8 0.4 0 0.89,带入上面的函数,得到构件 1 和构件 2 的角位移分别为 。.,583.21function y = RPPvel(x)% Input parameters% x(1)= theta-1 (构件 1 的角位移)% x(2) =theta-2 (构件 2 的角位移)% x(3)= theta-3 (构件 3 的角位移)% x(4)= dtheta-1 (构件 1 的角速度)% x(5)= r1 (构件 1 的杆长)% x(6)= r2 (构件 2 的杆长)% x(7)= r3 (构件 3 的杆长)% Output parameters% y(1)=dtheta-2% y(2)=dtheta-3%A = cos(x(1)+x(2) cos(x(2);sin(x(1)+x(2) sin(x(2);B = x(7)*cos(x(3);x(7)*sin(x(3);y = inv(A)*B;由位移分析计算出的 和曲柄 1 的角速度09.,5831.21及各个构件长度,则输入参数为srad/10x=1.5831 -0.0099 0 10 0.8 0.4 0带入上面的函数得到构件 1 和构件 2 的角速度分别为 。sradsrad/4.0,/8.0322.4 RPP 四杆机构 MATLAB 仿真结果南昌航空大学科技学院学士学位10曲柄 1 的幅度为 0和角速度等于 10rad/s 逆时针方向回转,相应各个构件的位移和速度的瞬时值。theta-1=0;dtheta-1=10;theta-i=1.5708;s-2=0;ds-2=0.8;s-3=0.89;ds-3=0.4。由于曲柄转速为 10rad/s,因此每转动 1 周的时间是 0.628s,用绘画命令plot(tout,simout(:,3),plot(tout,simout(:,4) ,plot(tout,simout(:,5)和 plot(tout,simout(:,6)分别绘画出构件 2 和构件 3 的速度和加速度。见下图(a)构件 2 的速度(纵坐标表示构件 2 的速度,单位为 ;横坐标表示时间,单位为 s)sm/南昌航空大学科技学院学士学位11(b)构件 3 的速度(纵坐标表示构件 3 的速度,单位为 ;横坐标表示时间,单位为 s)sm/(c)构件 2 的加速度(纵坐标表示构件 2 的加速度,单位为 ;横坐标表示时间,单位为 s)/sm南昌航空大学科技学院学士学位12(d)构件 3 的加速度(纵坐标表示构件 3 的加速度,单位为 ;横坐标表示时间,单位为 s)2/sm3 RPRRPP 六杆机构 MATLAB 运动学仿真3.1 RPRII 级杆组 MATLAB 运动学仿真模块如下图 RPRII 级杆组由 2 个转动副、一个移动副和导杆、滑块组成 RPRII 级杆组,导杆幅角 为变量,滑块 C 相对 D 的移动 S 也为变量,则导杆角速度 、滑块 C 相j j对 D 的加速度 和导杆上点 E 的加速度推导如下 sRPRII 级杆组的位置参数BseDCj移项整理得 sj南昌航空大学科技学院学士学位13上式对时间 t 求导并整理得: DBsesejjj )2/(上式对时间 t 求导并整理得: DBseeesese jjjjjj )2/()2/()2/()2/(把上式写成矩阵的形式: Bss jijjjjjjjjjjj ImRe)2/sin(co)2/in(cosinco)2/in(co 2 同样可以推出点 E 的加速度矩阵方程形式: 2)sin(co)2/sin(coImReI jjjjjjjjjj rrD RPRII 级杆组 MATLAB 运动学仿真模块 M 函数function y = RPRki(x)% function to compute the acceleration for RPR bar group% Input parameters% x(1) = r-3 (构件 3 的杆长)% x(2) = theta-3 (构件 3 的角位移)% x(3) = s-2 (构件 2 的相对位移)% x(4) = dtheta-3 (构件 3 的角速度)% x(5) = ds-2 (构件 2 的相对速度)% x(6) = ReddB (转动副 B 的水平加速度)% x(7) = ImddB (转动副 B 的垂直加速度)% x(8) = ReddD (转动副 D 的水平加速度)% x(9) = ImddD (转动副 D 的垂直加速度)% Output parameters%南昌航空大学科技学院学士学位14% y(1) = ddtheta -3 (构件 3 的角加速度) % y(2) = dds-2 (构件 2 的相对加速度)% y(3) = ReddE (转动副 E 的水平加速度)% y(4) = ImddE (转动副 E 的垂直加速度)%a=x(3)*cos(x(2)+pi/2) cos(x(2);x(3)*sin(x(2)+pi/2) sin(x(2);b=-x(5)*cos(x(2)+pi/2)+x(3)*x(4)*cos(x(2)+pi) x(4)*cos(x(2)+pi/2);x(5)*sin(x(2)+pi/2)+x(3)*x(4)*sin(x(2)+pi) x(4)*sin(x(2)+pi/2);b = b*x(4);x(5)+x(6)-x(8);x(7)-x(9);ddths = inv(a)*b;dde = x(8);x(9)+x(1)*ddths(1)*cos(x(2)+pi/2);x(1)*ddths(1)*sin(x(2)+pi/2)+x(1)*x(4)2*cos(x(2)+pi);x(1)*x(4)2*sin(x(2)+pi);y(1) = ddths(1);y(2) = ddths(2);y(3) = dde(1);y(4) = dde(2);各构件的初值为:r-3=1.6,theta-3 =1.5708,s-2 =1.4,dtheta-3 =2.8571,ds-2 =0,ReddD =0,ImddD=0 。3.2 RPRRPP 六杆机构图右所示是由原动件(曲柄1)和 一个 RPR 杆组、RPP 杆组所组成的 RPRRPP 六杆机构,各构件的尺寸为 r1=400mm,r3=1600mm,AD=1000,复数向量坐标见图所示,转动副 A 到滑块5的滑道的垂直距离为800mm,构件1以等角速度10rad/s 逆时针方向回转,试求构件4相对构件5的位移、速度、加速度和构件5的位移、速度、加速度。南昌航空大学科技学院学士学位15图 RPRRPP 六杆机构3.21 RPRRPP 六杆机构 MATLAB 仿真模型RPRRPP 六杆机构的 MATLAB 仿真模型如下图所示,在图中各积分模块的初值是以曲柄1的幅角为90和角速度等于10rad/s 逆时针方向回转,相应各构件的位移、速度的瞬时值分别为 dtheta-1=10;theta-1=1.5708;ds-2=0;s-2=1.4;dtheta-3=2.8571;theta-3=1.5708;ds-4=0;s-4=0.4。ds-5=-4.5714;s-5=0 。3个 MATLAB 函数模块分别为 crank.m、rprki.m 和 rppki.m。南昌航空大学科技学院学士学位16RPRRPP 六机构 Simulink 仿真模型南昌航空大学科技学院学士学位173.22 RPRRPP 杆组 M 函数为:function y = crank(x)% Function to compute the acceleration of crank% Input parameters% x(1) = r-1 (构件 1 的杆长)% x(2) = theta-1 (构件 1 的角位移)% x(3) = dtheta-1 (构件 1 的角速度)% x(3) = ddtheta-1 (构件 1 的角加速度)% Output parameters% 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 = ddB;function y =RPPki(x)% compute veloccity of RPP bar group% Input parameters% x(1) = theta - s4 (构件 4 的角位移)% x(2) = theta - s5 (构件 5 的角位移)% x(3) = ReddE (转动副 E 加速度的水平分量)% x(4) = ImddE (转动副 E 加速度的垂直分量)南昌航空大学科技学院学士学位18% x(5) = ReddK (转动副 K 加速度的水平分量)% x(6) = ImddK (转动副 K 加速度的水平分量)% Output parameters% y(1) = dds - 4 (构件 4 的加速度)% y(2) = dds - 5 (构件 5 的加速度)%a = cos(x(1)+x(2) cos(x(2);sin(x(1)+x(2) sin(x(2);b = x(3)-x(5);x(4)-x(6);ddsij = inv(a) * b; y(1) = ddsij(1);y(2) = ddsij(2);3.23 RPRRPP 六杆机构 MATLAB 仿真结果曲柄1的幅角为 和角速度等于10rad/s 逆时针方向回转,相应各构件的位移、速度的瞬时值分别为 dtheta-1=10;theta-1=1.5708;ds-2=0;s-2=1.4;dtheta-3=2.8571;theta-3=1.5708;ds-4=0;s-4=0.4。ds-5=-4.5714;s-5=0 由于曲柄转速为10rad/s,因此每转动1周的时间是0.628s,用绘画命令plot(tout,simout(:,1),plot(tout,simout(:,2),plot(tout,simout(:,3),plot(tout,simout(:,4)plot(tout,simout(:,5),plot(tout,simout(:6)绘制构件4相对构件5的位移、速度、加速度和构件5的位移、速度、加速度,如下图:南昌航空大学科技学院学士学位19(a)构件 4 相对构件 5 的位移(纵坐标表示构件 4 相对构件 5 的位移,单位为 m;横坐标表示时间,单位为 s)(b)构件 5 的位移(纵坐标表示构件 5 的位移,单位为 m;横坐标表示时间,单位为 s)(c)构件 4 相对构件 5 的速度(纵坐标表示构件 4 相对构件 5 的速度,单位为 m/s;横坐标表示时间,单位为 s)南昌航空大学科技学院学士学位20(d)构件 5 的速度(纵坐标表示构件 5 的速度,单位为 m/s;横坐标表示时间,单位为 s)(e)构件 4 相对构件 5 的加速度(纵坐标表示构件 4 相对构件 5 的加速度,单位为 ;横坐标表示时间,单位为 s)2/sm南昌航空大学科技学院学士学位21(f)构件 5 的加速度(纵坐标表示构件 5 的加速度,单位为 m/s2;横坐标表示时间,单位为 s)4 曲柄、RPP 杆组的 MATLAB 动力学分析4.1 曲柄 MATLAB 动力学仿真模块4.11 曲柄的动力学矩阵表达式如图4.1,已知曲柄 AB 向量的模 为常量,幅角 为常量,质点到转动副 A 的距离为iri,质量为 ,绕质点的转动惯量为 ,作用于质点上的外力为 和 ,外力矩cirimiJxiFyi为 ,曲柄与机架联接,转动副 A 的约束反力为 和 ,驱动力矩为 。iMxiRyi 1M图 4.1由理论力学可得: iiciiyBiciixBiciyBicixAi iiyiBxAx JrRrRrrRMsmF os)(sn)(onIe1由运动学知识可推出得:南昌航空大学科技学院学士学位22)sin()2/sin(ImcocoRe2 iiiii iiii rrAs 将上面的式子合并整理得: iiciiyBiciixBiciyicixAi iiiiiii xBiiiciiciyAx MrRrRrrRJ gmFmM os)(sn)(ossn)2/(I )o(e 21 曲柄 MATLAB 动力学仿真模块 M 函数:根据上式编写曲柄原动件 MATLAB 的 M 函数如下:function = crankdy_3(x)% Dynamic analysis for arank% Input parameters% x(1) = theta - 1 (构件1的角位移)% x(2) = dtheta - 1 (构件1的角速度)% x(3) = ddtheta - 1 (构件1的角加速度)% x(4) = -RxB (转动副 B 反作用力的水平分量)% x(5) = -RyB (转动副 B 反作用力的垂直分量)% % Output parameters% y(1) = RxA (转动副 A 反作用力的水平分量) % y(2) = RxB (转动副 A 反作用力的垂直分量)% y(3) = M1 (曲柄1上的驱动力矩)%g = 9.8;ri = 0.4; rci = 0.2;mi = 1.2; Ji = 0.016;Fxi = 0; Fyi = 0; Mi = 0;Redda = 0; Imdda = 0;南昌航空大学科技学院学士学位23y(1) = mi * ReddA + mi * rci * x(3) * cos(x(1) + pi/2) + mi * rci * (2)2*cos(x(1) + pi)-Fxi + x(4);y(2) = mi * ImddA + mi * rci * x(3) * sin(x(1) + pi/2) + mi * rci * (2)2*sin(x(1) + pi)-Fxi + 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;各构件的初值为:theta - 1=1.5708,dtheta - 1=10, ddtheta - 1=0 。4.2 RPRII 级杆组的动力学仿真模块4.21 RPRII 级杆组动力学矩阵表达式RPRII 级杆组由2个构件滑块和导杆组成。滑块的质量为 ,导杆的质量为 ,imjm转动惯量为 。滑块和导杆的受力分析图如下所示,则转动副 B 和 E 及移动副 CjJ的约束反力推导如下。由图 b 受力分析得: BmgRFiijcyiyBijxix Iosein由图 c 受力分析得: jjyDjcjxDciCj jjjyjj jjxjj JRrRrsMsgs coin)( Ioein南昌航空大学科技学院学士学位24由动力学可推出得: )sin()2/sin(ImcocoRe2 jjjjj jjjjj rrDs 将上式分别代入得: )sin()2/sin(Imcos cocoRein 2 jjjjjyDjjC jjjjjxjj rrFR 将上式整理成矩阵为: jj jyjjjcjjc xj yiiixyDxCyBxjcjjcjcjjj MJ gmFrrmFgRrrs )sin()2/sn(I ooeosin010oics1in0 24.22 RPRII 级杆组 MATLAB 动力学仿真模块 M 函数根据上式编写 RPRII 级杆组 MATLAB 的 M 函数如下:function y =RPRdy_2(x)% Dynamic analysis of RPR bar group% Input parameters% x(1) = theta-3 (构件3 的角位移)% x(2) = s-2 (构件2的相对位移)% x(3) = dtheta-3 (构件3 的角速度)% x(4) = ddtheta-3 (构件3 的角加速度)% x(5) = ReddB (转动副 B 加速度的水平分量)% x(6) = ImddB (转动副 B 加速度的垂直分量)% x(7) = -RxE (转动副 E 的约束反力的水平分量)% x(8) = -RyE (转动副 E 的约束反力的垂直分量)% x(9) = M3 (构件3上的力矩)南昌航空大学科技学院学士学位25% % Output parameters% y(1) = RxB (运动副 B 的约束反力的水平分量)% y(2) = RyB (运动副 B 的约束反力的垂直分量)% y(3) = RC (移动副 C 的约束反力)% y(4) = RxD (移动副 D 的约束反力的水平分量)% y(5) = RyD (移动副 D 的约束反力的垂直分量)%g = 9.8;rcj = 0.8;mi = 1.5; mj = 10;Jj = 2.2;ReddD=0;ImddD=0;Fxi = 0; Fyi = 0;a = zeros(5);a(1,1) = 1; a(1,3) = -sin(x(1);a(2,2) = 1; a(2,3) = cos(x(1);a(3,3) = sin(x(1); a(3,4) = 1;a(4,3) = -cos(x(1); a(4,5) = 1;a(5,3) = - (x(2) - rcj);a(5,4) = -rcj*sin(x(1);a(5,5) = -rcj*cos(x(1);b = zeros(5,1);b(1,1) = mi*x(5) - Fxi;b(2,1) = mi*x(6) + mi*g-Fyi;b(3,1) =mj*(ReddD+rcj*x(4)*cos(x(1)+pi/2)+rcj*x(3)2*cos(x(1)+pi)+x(7);b(4,1) =mj*(ImddD+rcj*x(4)*sin(x(1)+pi/2)+rcj*x(3)2*sin(x(1)+pi)+x(8)+mj*gb(5,1) = Jj*x(4)-x(9);y = inv(a) * b;各构件的初值为:theta-3=1.5708,s-2 =1.4,dtheta-3=2.8571 。 4.3 RPPII 级杆组的动力学仿真模块南昌航空大学科技学院学士学位264.31 RPPII 级杆组动力学矩阵表达式RPPII 级杆组由2个滑块 i 和 j 组成,2个滑块的滑动方向成固定角度 ,滑块 j 的方向分别为受力分析对象做出受力图分别如图所示,有关参数标注于图中,下面将推导转动副 B 的约束力、移动副 C 的约束反力和移动副 D 的约束力和约束反力矩。由图(b)受力分析得: BmRFjjCyiyB jjxix I)cos(Rein由图(c)受力分析得: jCcjcjjDj jjjjjyj jjjjxj JRrsMsgo)( Ioe)in(i南昌航空大学科技学院学士学位27上式联立整理成矩阵形式为: jj jyjxjj iyiixDCyBxcjjcj jjj MJgmFsmRrsr ncoIe1os00)o(inincs1)i(0 4.32 RPPII 级杆组 MATLAB 动力学仿真模块 M 函数根据上式编写 RPPII 级杆组 MATLAB 的函数如下:function y = RPPdy(x)% Dynamic analysis of RPP bar group% Input parameters% x(1) = s-4 (构件4的相对位移)% x(2) = dds-5 (构件5的加速度)% x(3) = ReddE (转动副 E 的加速度的水平分量)% x(4) = ImddE (转动副 E 的加速度的垂直分量)% % Output parameters% y(1) = RxE (转动副 E 的约束反力的水平分量)% y(2) = RyE (转动副 E 的约束反力的垂直分量)% y(3) = RF (转动副 F 的约束反力)% y(4) = RG (转动副 G 的约束反力)% y(5) = MG (转动副 G 的力矩)%g = 9.8;rcj = 0.2;thj = 0; delta = -pi/2; ddthj = 0;mi = 1.5; mj = 20;Jj = 2;南昌航空大学科技学院学士学位28Fxi = 0; Fyi = 0;Fxj = 1000; Fyj = 0;Mj = 0;a = zeros(5);a(1,1) = 1; a(1,3) = -sin(thj + delta);a(2,2) = 1; a(2,3) = cos(thj + delta);a(3,3) = sin(thj + delta); a(3,4) = -sin(thj);a(4,3) = -cos(thj + delta); a(4,4) = cos(thj);a(5,3) = -(x(1) - rcj);a(5,4) = -rcj*cos(delta);a(5,5) = 1;b = zeros(5,1);b(1,1) = mi*x(3)-Fxi;b(2,1) = mi*x(4)-Fyi + mi*g;b(3,1) = mj*x(2)*cos(thj)-Fxj;b(4,1) = mj*x(2)*sin(thj)-Fyj + mj*g;b(5,1) = Jj *ddthj - Mj;y = inv(a)*b;4.4 RPRRPP 六杆机构 MATLAB 动力学仿真由原动件(曲柄 1)和一个 RPR 杆组、RPP 杆组所组成的 RPRRPP 六杆 II 级机构。各机构的尺寸为 转动副 A 到滑块的滑道的垂mADrmr 10,160,403直距离为 800mm,各构件质心在构件的中心处,各构件的质量为;转动惯量kgkgkgm2,5.7,5.1,.24321构件 1 以等角速度 10rad/s 逆时针方29004JJJ 向回转,执行构件 5 的工作阻力 F=1000N,试求在不计摩擦力时,转动副 A 的约束反力、驱动力矩 移动副 G 的约束反力和力矩以及驱动力矩 所作的功。1M1M南昌航空大学科技学院学士学位29RPRRPP 六杆机构南昌航空大学科技学院学士学位30RPRRPP 六杆机构仿真模型南昌航空大学科技学院学士学位314.41 RPRRPP 六杆机构 MATLAB 仿真模型用 MATLAB/Simlink 对该曲柄滑块机构动力学仿真求解转动副 A 的约束反力、驱动力矩 ,移动副 G 的约束反力和力矩以及驱动力矩 所作的功。用到曲柄原动件、1M1MRPR 杆组和 RPP 杆组的 MATLAB 3个运动仿真模块和曲柄原动件、RPR 杆组和 RPP 杆组的 MATLAB 3个动力学仿真模块。function y = crank_3(x)% Function to compute the acceleration of crank% Input parameters% x(1) = theta-1 (构件1的角位移)% x(2) = dtheta-1 (构件1的角速度)% x(3) = ddtheta-1 (构件1的角加速度)% Output parameters% y(1) = ReddB (转动副 B 加速度的水平分量)% y(2) = ImddB (转动副 B 加速度的垂直分量)%r1=0.4ddB = r1*x(3)*cos(x(1)+pi/2)+r1*x(2)2*cos(x(1)+pi);r1*x(3)*sin(x(1)+pi/2)+r1*x(2)2*sin(x(1)+pi);y = ddB;各构件的初值为:theta-1=1.5708, dtheta-1=10, ddtheta-1=0 。function y = RPRki_1(x)% function to compute the acceleration for RPR bar group% Input parameters南昌航空大学科技学院学士学位32% x(1) = theta-3 (构件3的角位移 ) % x(2) = s-2 (构件2的相对位移)% x(3) = dtheta-3 (构件3的角速度 )% x(4) = ds-2 (构件2的相对速度)% x(5) = ReddB (转动副 B 加速度的水平分量)% x(6) = ImddB (转动副 B 加速度的垂直分量)% Output parameters% y(1) = ddtheta -3 (构件3的角加速度)% y(2) = dds-2 (构件2的相对加速度)% y(3) = ReddE (转动副 E 加速度的水平分量)% y(4) = ImddE (转动副 E 加速度的垂直分量)%r2=1.6;rd=0;id=-1;a=x(2)*cos(x(1)+pi/2) cos(x(1);x(2)*sin(x(1)+pi/2) sin(x(1);b=-x(4)*cos(x(1)+pi/2+x(2)*x(3)*cos(x(1)+pi) x(3)*cos(x(1)+pi/2);x(4)*sin(x(1)+pi/2+x(2)*x(3)*sin(x(1)+pi) x(3)*sin(x(1)+pi/2);b = b*x(3);x(4)+x(5)-rd;x(6)-id;ddths = inv(a)*bdde = rd;id+r2*ddths(1)*cos(x(1)+pi/2);r2*ddths(1)*sin(x(1)+pi/2)+r2*x(3)2*cos(x(1)+pi);r2*x(3)2*sin(x(1)+pi);y(1) = ddths(1);y(2) = ddths(2);y(3) = dde(2);y(4) = dde(2);各构件的初值为:theta-3=1.5708, s-2 =1.4, dtheta-3 =2.8571,ds-2 =0 。function y =RPPki_1(x)%南昌航空大学科技学院学士学位33% compute veloccity of RPP bar group% Input parameters% x(1) = ReddE (转动副 E 加速度的水平分量)% x(2) = ImddE (转动副 E 加速度的垂直分量)% Output parameters% y(1) = dds-4 (构件4的加速度)% y(2) = dds-5 (构件5的加速度)%s(4)=1.5708;s(5)=0;rk=0;ik=0;a = cos(s(4)+s(5) cos(s(5);sin(s(4)+s(5) sin(s(5);b = x(1)-rk;x(2)-ik;ddsij = inv(a) * b; y(1) = ddsij(1);y(2) = ddsij(2);function y = M3_1(x)% Compute moment of bar 3% % Input parameters% x(1)=theta3 (构件3的角位移)% x(2)=-RxE (转动副 E 约束反力的水平分量)% x(3)=-RyE (转动副 E 约束反力的垂直分量)% Output parameters% y(1)=M3 (构件3的力矩)南昌航空大学科技学院学士学位34r3=1.6;rc3=0.8;%y(1)=-x(3)*(r3-rc3)*cos(x(1)+x(2)*(r3-rc3)*sin(x(1)4.42 RPRRPP 六杆机构 MATLAB 仿真结果由于 MATLAB/Simulink 仿真默认的时间变量 tout,因此在仿真过程中未再设置时间变量。用绘图命令:plot(tout,simout1(:,1),plot(tout,simout1(:,2),plot(tout,simout1(:,3),plot(tout,simout1(:,4),plot(tout,simout3(:,4)和 plot(tout,simout3(:,5)分别绘制出转动副 A 的约束反力,驱动力矩 M1、驱动力矩 M1 所作的功、移动副 G的约束反力和力矩。(a)转动副 A 的水平方向力(纵坐标表示转动副 A 的水平方向力,单位为 N;横坐标表示时间,单位为 s)南昌航空大学科技学院学士学位35(b)转动副 A 的垂直方向力(纵坐标表示转动副 A 的垂直方向力,单位为 N;横坐标表示时间,单位为 s)(c)曲柄上作用的力矩(纵坐标表示曲柄上作用的力矩,单位为 ;横坐标表示时间,单位为 s)mN南昌航空大学科技学院学士学位36(d)曲柄力矩作的功(纵坐标表示曲柄力矩作的功,单位为 J;横坐标表示时间,单位为 s)(e)移动副 G 上作用的约束反力(纵坐标表示移动副 G 上作用的约束反力,单位为 N;横坐标表示时间,单位为 s)南昌航空大学科技学院学士学位37(f)移动副 G 上作用的约束反力矩(纵坐标表示移动副 G 上作用的约束反力矩,单位 ;横坐标表示时间,单位 s)mN总结本文是用 MATLAB 的软件对 RPP 平面连杆的仿真(运动仿真和动力仿真),在做RPP 平面连杆仿真上进行了对连杆的运动和动力分析,并编制了用于MATLAB/simulin 仿真的 M 函数。现总结如下:对运动学分析:要求我仔细的分析连杆的运动,滑块的运动方向。对动力学分析:要求我仔细的分析连杆的受力和运动,滑块的运动方向.转动惯量,约束反力。编制 M 函数:要求对运动学对位移和速度的分析准确,编制 M 函数才能运行。通过努力,我成功地完成了毕设的任务。但其中充满了艰辛,刚开始拿到课题,我对 MATLAN 很陌生,完全从零开始。在老师的指导帮助和自己的不懈努力下,我一步一步地解决问题,享受着学习的快乐,沉浸在知识的海洋。这期间我学到了好多东西,收获了丰硕的果实。其实,好多东西我们不懂,但只要认真地行动起来,没有什么可以阻碍我们。南昌航空大学科技学院学士学位38由于时间和能力有限,该设计还存在不足,望各位领导、各位老师批评指正,我在此谢谢大家,也祝愿我国的 MARLAB 运用发展更上一层楼。参考文献1 孙桓,陈作模主编.机械原理.第七版.北京:高等教育出版社,2006.122 曲秀全主编.基于 MATLAB/Simulink 平面连杆机构的动态仿真.哈尔滨:哈尔滨工业大学出版社,2007.43 邱晓林主编. 基于 MATLAB 的动态模型与系统仿真工具. 西安:西安交通大学出版社,2003.104 张策主编. 机械动力学. 北京:高等教育出版社,20005 Ye Zhonghe, Lan Zhaohui. Mechanisms and Machine Theory. Higher Education Press, 2001.76陈桂明等.应用 MATLAB 建模与仿真.北京:科学出版社,2001南昌航空大学科技学院学士学位39致谢本课题及论文是在我的指导老师朱保利博士的悉心指导下完成的。课题从实施到
收藏