[工学]机械原理课程设计——牛头刨床.doc

《[工学]机械原理课程设计——牛头刨床.doc》由会员分享，可在线阅读，更多相关《[工学]机械原理课程设计——牛头刨床.doc（33页珍藏版）》请在三一办公上搜索。

1、牛头刨床运动方案改进和对比分析 专 业： 机械设计制造及其自动化 学 号： 2009100XXXX 姓 名： XXXXXXX 指导教师： XXXXXX 日 期： 2011.7.8 目录一、 概述31.1、课程设计的题目31.2、课程设计的任务和目的31.3、课程设计的要求3二、 方案一运动分析（原方案）2.1、方案一的原始数据32.2、机械运动简图42.3、运动综合设计42.4、运动分析的数学模型5三、方案二运动分析（新方案）3.1、方案二的原始数据 73.2、机械运动简图 93.3、运动综合设计 93.4、运动分析的数学模型 93.5、受力分析的数学模型11四、 公共部分设计4.1、行星轮系

2、设计144.2、变位齿轮设计16五、两种方案的比较18六、设计心得 20七、参考文献 21八、附录 21一、 概述1.1、课程设计的题目牛头刨床运动方案设计和分析1.2、课程设计的任务和目的1） 任务：1、 设计更完善的机构和运动方案2、 对倒杆机构进行运动分析3、 导杆机构的动态静立分析根据要求发挥自己的创新能力，设计牛头刨床的主传动机构，使其可以满足牛头刨床的传动需要2） 目的：机械原理课程设计是培养学生掌握机械系统运动方案设计能力的技术基础课程，它是机械原理课程学习过程中的一个重要实践环节。其目的是以机械原理课程的学习为基础，进一步巩固和加深所学的基本理论、基本概念和基本知识，培养学生分

3、析和解决与本课程有关的具体机械所涉及的实际问题的能力，使学生熟悉机械系统设计的步骤和方法，其中包括选型、运动方案的确定、运动学和动力学的分析和整体设计等，并进一步提高计算、分析、计算机辅助设计、绘图以及查阅和使用文献的综合能力。1.3、课程设计的要求牛头刨床的主传动的从动机构是刨床，在设计主传动机构时，要满足所涉及的机构要能使牛头刨床正常的运转，同时设计的主传动机构的行程要有急回运动的特性，以及很好的动力特性，尽量使设计的结构简单，实用，能很好的实现传动功能。二、 方案一运动分析2.1、方案一的原始数据H：刨头行程 ； K：行程速比系数； Fc切削阻力 ； 分别为导杆、连杆及刨头的质量；分别为

4、导杆4及导杆5绕各自质心的转动惯量；分别为行星减速器中心轮及齿轮4、5的模数；为齿轮4及5的齿数；：电机转速；：曲柄2及齿轮5的转速；k：行星轮个数。符 号方案运动分析和运动综合动力分析HmmKmmkgkgkgkg56001.83700.50.30.5223520.90.0151400符 号方案行星轮系设计变位齿轮rpmrpmK类型mm0.520mm510008032K-H5144916202.2、机构运动简图2.3运动综合设计设 1、导杆的摆角K=1.82、导杆的长度L43、连杆的长度L54、刨头导路中心线xx至O3点的垂直距离L6根据已知xx被认为通过圆弧BB的绕度ME的中点D知5、曲柄的

5、长度6、切削越程长度0.05H，如图所示则切削越程长度为0.05H=0.05500=25mm7、计算机构的自由度F=3527=12.4、运动分析的数学模型如图所示，先建立一直角坐标系，并标出各杆矢量及其方位角。其中共有四个未知量、。为求解需建立两个封闭的矢量方程，为此需利用两个封闭的图形及，由此可得：并写成投影方程为：由上述各式可解得：由以上各式即可求得、四个运动变量，而滑块的方位角=。然后，分别将上式对时间取一次、二次导数，并写成矩阵形式，及得一下速度和加速度方程式。+根据以上各式，将已知参数代入，即可应用计算机计算，应用MATLAB编程并求得的数列与表1中(见附录)。并根据所得数值绘出机构

6、的位置线图、速度线图、加速度线图（见方案比较）。三、方案二运动分析3.1、课程设计的原始数据H：刨头行程 ； K：行程速比系数； Fc切削阻力 ； m4 m5 m6分别为导杆、连杆及刨头的质量；为扇形齿轮质心到的距离；J4、J5分别分别为导杆4及导杆5绕各自质心的转动惯量；m1、mH分别为行星减速器中心轮及齿轮4、5的模数；Z4，Z5为齿轮4及5的齿数；n1：电机转速；n2：曲柄2及齿轮5的转速；k：行星轮个数。运动分析和运动综合动力分析变位齿轮HKFc单位mmmmkg kg数据6001.83706723652814001620 行星轮设计扇形齿轮及齿条n1n2K类型m1单位rpm2K-Hmm

7、mmmm数据1000803 5144943365620285.23.2、机械运动简图3.3、运动综合设计以下计算所涉及的尺寸： = 1、 导杆的摆角:2、 =0.8929rad3、 行程速比系数K: K=4、 曲柄长度 =370mm =159.8mm5、 切削越程长度0.05H，如图所示H=600mm 则切削越程长度为0.05H=0.05600=30mm6、计算机构的自由度F=3.4、运动分析数学模型如图所示，先建立一个直角坐标系，并标出各杆矢量及其方位角，如下图所示。其中共有2个未知量。为求解需建立一个封闭的矢量方程，为此需利用这个封闭的图形，由此可得：并写成投影方程为由上述各式可解得： 由

8、以上各式即可求得两个运动变量然后，分别将上式对时间取一次、二次导数，并写成矩阵形式，及得一下速度和加速度方程式。 另外， 根据以上各式，将已知参数代入，即可应用计算机计算，应用MATLAB编程并求得的数列与表2中。并根据所得数值作出机构的位置线图、速度线图、加速度线图。这些线图称为机构的运动线图。通过这些线图可以一目了然的看出机构在一个运动循环中位移、速度、加速度的变化情况，有利于进一步掌握机构的性能。3.5、动态静力分析数学模型用矩阵法对机构进行受力分析，须先建立一个直角坐标系，如下图将各力都分解为沿两坐标轴的两个分立，再分别就各个构件列出它们的平衡力方程，具体过程如下：对于曲柄1 受力分析

9、如图可分别依据，以及 列出三个平衡方式，并将含有未知要素的项写在等号左边，故有: 同理，对于构件2受力分析如下图 其中 对构件3受力分析如下图 以上十一个方程式，可按构件1、2、3尚待定的位置离的次序整理成矩阵式，如下：=上式即为该机构的静态动力分析的矩阵方程。矩阵方程可以简化为下列形式:式中，为未知力的列阵，而D和C分别为已知力和未知力的矩阵系数。应用计算机计算，用MATLAB编程并求得的数列于表2、表3中（见附录）。并根据所得数值作出机构的惯性力及惯性力矩图、原动件上的平衡力矩Mp图。求得各运动副中的约束反力如表4（见附录），并根据所得数据绘制各运动副中的约束反力图，如下：1、O1点对曲柄

10、的约束反力2、O2点对扇形齿轮的约束反力3、工作台对刨头的约束反力及曲柄和扇形齿轮间的约束力四、公共部分设计4.1、行星轮系设计已知， ，在行星轮系中，各齿轮数的选配必须满足以下四个条件：（1） 尽可能近似的实现给定的传动比因为 故有 （2） 满足同心要求要行星轮系能正常运转，其基本构件的回转轴线必须在同一直线上，此此即同心条件。为此，对于本机构的轮系来说，必须满足：当采用标准齿轮传动或等变位齿轮传动时，上式变为：（3) 满足均布要求为使各星轮能均布的安装，行星轮的个数和齿数之间必须满足一定的关系，否则将会因行星轮和太阳轮的干涉而不能安装。由公式如这时太阳轮1恰好转过整个齿数N,即式中，N为整

11、数；为太阳轮1的齿距角。两式结合后，得由以上公式可得：为使结构架凑，此处k=3，且均为整数取（4） 验证是否满足满足邻接条件为了保证相邻两个行星轮不会相互碰撞，需使中心距大于两轮齿顶圆的半径之和，对于标准齿轮传动有得到64.95135，满足上述条件，故能顺利装配。4.2、变位齿轮的设计已知：(1) 确定传动类型 所以传动类型为等变位齿轮传动。(2) 确定变位系数对于变位齿轮，为有利于强度的提高，小齿轮4采用正变位，大齿轮5采用负变位，使大小齿轮的强度趋于接近，从而使齿轮承载能力提高。则取满足(3) 计算几何尺寸名称符号公式结果变位系数节圆直径啮合角齿顶高齿根高齿顶圆直径齿根圆直径中心距中心距变

12、动系数0齿顶高降低系数0（4） 变位齿轮的重合度的验证重合度的计算：其中 因为1，所以这对变位齿轮能连续传动(5) 齿顶圆齿厚验证由公式得其中所以6.5mm因为齿轮4为正变位齿轮，所以其齿顶圆齿厚将减小，所承载的压力将减小，所以齿顶圆齿厚不能过小，及即齿顶圆不能很尖，6.5mm，可满足正常的传动需要，故可以正常使用。五、 两种方案的比较由对两种方案的运动分析并通过上机应用MATLAB绘制出其运动曲线图，如下：1、 该图为方案一刨头位移和方案二刨头位移（其中计算时初始位置不同）结论：有上图可以看出两种方案，两刨头的位移以相同的斜率增加，并同时达到最高点，并以相同的斜率下降达到最低点，说明两种方案

13、中刨头的位移规律是一样的。两种刨头有相似的效率。2、方案一刨头速度和方案二刨头速度结论：由上图可以看出两种总方案的包头速度趋于重合，也就是说包头的速度一样，而方案二的最低速度大于方案一的最低速度，所以方案二较方案一更平稳。3、方案一刨头速度(绿色)和方案二刨头速度(红色)结论：由上图不难看出两曲线近乎一致，说明方案一和方案二的刨头有相同的加速度，但从细节中不难看出，方案二的最低加速度和最大加速度的绝对值均小于方案已相对应加速度的绝对值，说明方案二更为稳定。综合分析：（1）运动是否具有确定的运动 新方案的机构由齿条3、扇形齿轮2、滑块4和曲柄1组成，其中曲柄1为主动件。机构中有一个原动件，原动件

14、的个数等于该机构的自由度。所以，该机构具有确定的运动。且原动件较原始方案构件少，传动效率高，机械效率损失少等特点。（2）机构功能的实现 根据新方案机构简图可知，整个机构的运转是由原动件1带动的。杆1通过滑块4带动扇形齿轮2的运动。扇形齿轮2和与刨头连接的齿条啮合。从而实现刨刀的往复运动。（3）机构的工作性能 该机构中原动件1对滑块4的压力角一直在改变。但是原动件1的长度较小，扇形齿轮的半径较大，即原动件1的变化速度对于扇形齿轮2的影响不是很大，同时机构是在转速不大的情况下运转的，也就是说，在扇形齿轮作用下的齿条的速度在切削过程中变化不大。趋于匀速运行。 原动件1在滑块4上的速度始终不变，但是随

15、着原动件1的运转，在一个周期里，的长度由小到大，再变小。而的长度是扇形齿轮3的回转半径，也就是说，在机构的运行过程中，推程的速度趋于稳定，在刨头回程时，由于扇形齿轮受到齿条的反作用力减小。还有扇形齿轮3的回转半径减小，使扇形齿轮的回程速度远大于推程时的速度。即可以达到刨床在切削时速度较低，在回程时有速度较高的急回运动。（4）机构的传递性能 该机构中除了有扇形齿轮和齿条接触的两个高副外，所有的运动副都是低副，齿轮接触的运动副对于载荷的承受能力较强，所以，该机构对于载荷的承受能力较强，适于加工一定硬度的工件。同时。扇形齿轮是比较大的工件，强度比较高，不需要担心因为载荷的过大而出现机构的断裂。 在整

16、个机构的运转过程中，原动件1是一个曲柄，扇形齿轮3只是在一定的范围内活动，机构的是设计上不存在运转的死角，机构可以正常的往复运行。（5）机构的动力性能分析 该机构的主传动机构采用导杆机构和扇形齿轮，齿条机构。齿条固结于刨头的下方。扇形齿轮的重量较大，运转时产生的惯量也比较大，会对机构产生一定的冲击，使机构震动，不过在低速运转情况下，影响不会很大。（6）机构的合理性 该机构的设计简单，尺寸可以根据机器的需要而进行选择，不宜过高或过低。同时，扇形齿轮的重量有助于保持整个机构的平衡。使其重心稳定。由于该机构的设计较为简单。所以维修方便，除了齿轮的啮合需要很高的精确度外没有什么需要特别设计的工件，具有

17、较好的合理性。（7）机构的经济性能 该机构中扇形齿轮与齿条的加工的精度要求很高，在工艺上需要比较麻烦的工艺过程，制作起来不是很容易。此方案经济成本较高。综上：方案二具有构件少，使用齿轮传动更精确，效率高，能量损失少，更适合加工硬度较高的工件，总体效果更稳定等特点六、心得体会 这次的机械原理课程设计，让我提高了综合运用机械原理课程理论的能力，养了分析和解决一般机械运动实际问题的能力，并使所学知识得到进一步巩故、深化和扩展。掌握了一些常用执行机构、传动机构或简单机器的设计方法过程 。 这这次机会让我利用一整段时间巩固旧的知识和学习新的知识，把所学用到现实生活当中。在所学理论知识的基础上也充分的发挥

18、了创造性。各类资查也熟练了很多。自己的计算机绘图水平也有了一定的提高，并对所学的专知识有了进一步的理解。 不过，作为自己的第一次课程设计，其中有很多的不足之处，希望在今后的学和生活中，能够得到改正，使自己日益臻于成熟，专业知识水平日益深厚。“君与善其事必先利其器”只有通过不锻炼,调高自己知识水平这个“利器”，能迎接更大的挑战和机遇,我相信自己一定能够在未来的道路上得到锻炼成长.七、 参考文献1、机械原理（第七版） 孙恒，陈作模，葛文杰，高等教育出版社，20002、机械原理课程设计手册 邹慧君主编，高等教育出版社，19983、MATLAB从入门到精通 周建兴、岂兴明、矫津毅、 等 人民邮电出版社

19、 (2008-11出版)4、机械原理与机械设计(上册)(第2版) 张策 机械工业出版社 (2011-01出版)八、 附录8.1、应用MATLAB求解表1 各构件的位置速度和加速度(m/s)12345678910111201.4832.16842.46672.53162.39341.98611.0881-0.8526-4.4811-6.2197-2.669633.79515.5997.16472.6745-0.54151-4.0052-9.4584-20.434-43.778-63.92227.4460.5840-1.9153-1.9863-1.08120.238991.47072.10841.

20、5072-1.2241-4.11951.88653.3128-49.55-13.0088.374718.96821.58816.1872.4634-24.374-60.4035.0171121.07-52.2540-0.90904-1.4568-1.7216-1.7416-1.5937-1.3055-0.726510.572052.94764.34211.6635-19.357-11.043-6.4099-2.08461.21893.36626.141213.40629.20242.83-18.229-42.199表2 各构件的位置速度和加速度(m/s)12345678910111201.47

21、782.16182.45962.52372.38411.97441.0725-0.87068-4.4761-6.1657-2.647333.615.5567.1542.6659-0.5572-4.0342-9.5085-20.503-43.735-63.09327.52259.99600.99311.45271.65281.69591.60211.32680.72075-0.5851-3.0079-4.1433-1.77922.57910.4544.80751.7915-0.3744-2.711-6.3897-13.778-29.39-42.39918.49540.317表3 各机构的惯性力及

22、惯性力矩、原动件上的平衡力矩Mp：F2M2F3Mp123456789101112-344.98-159.72-73.451-27.3715.720441.4297.626210.51449.03647.79-282.57-615.99-268.8-124.45-57.232-21.3274.457232.27476.068164.02349.88504.75-220.18-479.97-812.86-376.33-173.07-64.49313.47997.596230.03496.0110581526.4-665.81-1451.4-110.28-1106.41009.12403.52761

23、.41998.4307.0235.096185.29444.42-223.77-345.348表4 求得各运动副中的约束反力:123456789101112805.465640.25566.35526.75498.45467.75419.5329.07124.51-45.948751.941037.9-799.03-20794-18587-17634-17406-17777-18904395.381195.42698.9-1395.1-2083.6-382.62-8646.7-5277.1-2296.9497.463315.36338.4-178.05-558.68-764.46-127.49

24、-829.31324.97-6550.4-4623.7-3822.7-3704.7-4195.3-5506.4-292.57-269.4549.28-447.89-59.823588.563930.2593.27-2363.9-5133.3-7913-10865797.831486.71849.3263.67881.648.2、源程序clear all;clc;%初始条件theta0=0.446429;theta1=linspace(-theta0,2*pi-theta0,120); %单位度H=0.6; %行程 单位mL1=0.1598; %曲柄L1的长度 单位mL2=0.672; %L2的

25、长度 单位mL4=0.37; %O1O2的长度 单位mL5=0.672; %DO2的长度 单位mW1=80*pi/30; %角速度 单位rad/s%初始条件J2=8;m2=36;m3=52;yc=0.2852;%求解S2、theta2,S3S2=(L4)2+(L1)2-2*L4*L1*cos(theta1+pi/2).0.5; %求出S3的值for i=1:120 %求解theta2 theta2(i)=acos(L1*cos(theta1(i)/S2(i); S3(i)=(theta2(i)-theta2(1)*L2;end %求解完成%求解VS2,2,V3for i=1:120 J=inv

26、(cos(theta2(i),-S2(i)*sin(theta2(i);%矩阵逆运算 sin(theta2(i),S2(i)*cos(theta2(i); K=J*W1*-L1*sin(theta1(i);L1*cos(theta1(i); VS2(i)=K(1); % S2一次求导 即S2速度 W2(i)=K(2); % 2角速度2 V3(i)=W2(i)*L2; %刨刀移动速度end %求解完成%求解aS2，a2,a3刨刀加速度for i=1:120 J=inv(cos(theta2(i),-S2(i)*sin(theta2(i);%矩阵逆运算 sin(theta2(i),S2(i)*co

27、s(theta2(i); P=W1*W1*-L1*cos(theta1(i);-L1*sin(theta1(i); M=-W2(i)*sin(theta2(i),-VS2(i)*sin(theta2(i)-S2(i)*W2(i)*cos(theta2(i); W2(i)*cos(theta2(i),VS2(i)*cos(theta2(i)-S2(i)*W2(i)*sin(theta2(i); N=VS2(i);W2(i); K=J*(-M*N+P); aS2(i)=K(1); %L2线加速度 a2(i)=K(2); %L2角加速度 a3(i)=a2(i)*L2;end %求解结束k=1;for

28、 i=1:12 qW2(i)=W2(k); qa2(i)=a2(k); qV3(i)=V3(k); qa3(i)=a3(k); k=k+10;endBBB=qW2;qa2;qV3;qa3;%动态静力分析s%切削阻力赋值for i=1:120 if(abs(S3(1)-S3(i)0.05*H & abs(S3(1)-S3(i)0.95*H) & theta1(i)(-theta0) Fc(i)=-1400*9.8; else Fc(i)=0; endend%完成求解for i=1:120 F3(i)=-m2*a3(i); M2(i)=-J2*a2(i); F23x(i)=-F3(i)-Fc(i)

29、; F23y(i)=F23x(i)*tan(pi/9); Fr3(i)=m3*9.8+F23y(i); F2x(i)=-a2(i)*yc*sin(theta2(i)*m2+m2*W2(i)*W2(i)*yc*sin(theta2(i); F2y(i)=-a2(i)*yc*cos(theta2(i)*m2+m2*W2(i)*W2(i)*yc*cos(theta2(i); F2(i)=(F2x(i)2+F2y(i)2).0.5; J=inv(tan(theta2(i)*L1*sin(theta1(i)-L1*cos(theta1(i),1,0,0; -tan(theta2(i),0,1,0; 1,

30、0,0,1; S2(i)*cos(theta2(i)+tan(theta2(i)*S2(i)*sin(theta2(i),0,0,0); K=0;F23x(i)-F2x(i);m2*9.8-F23y(i)-F2y(i);F2x(i)*yc*sin(theta2(i)+(m2*9.8-F2y(i)*yc*cos(theta2(i)-F23x(i)*L5-M2(i); GG=J*K; Foy(i)=GG(1); Fox(i)=Foy(i)*tan(theta2(i); Mp(i)=GG(2); Fo2x(i)=GG(3); Fo2y(i)=GG(4); F12x(i)=Fox(i); F12y(i

31、)=Foy(i);endk=1;for i=1:12 uF2(i)=F2(k); uM2(i)=M2(k); uF3(i)=F3(k); uFr3(i)=Fr3(k); uFox(i)=Fox(k); uFoy(i)=Foy(k); uFo2x(i)=Fo2x(k); uFo2y(i)=Fo2y(k); uF23x(i)=F23x(k); uF23y(i)=F23y(k); uMp(i)=Mp(k); k=k+10;endFF=uF2;uM2;uF3;uMp;uFr3;uFox;uFoy;uFo2x;uFo2y;%方案二%初始条件stheya0=0.448811;stheta1=linspac

32、e(-stheya0,2*pi-stheya0,120);%单位度sW1=80*pi/30;%角速度 单位rad/ssH=0.6;%行程 单位msL1=0.1605;%O2A的长度 单位msL4=0.6914;%O3B的长度 单位msL5=0.2074;%BF的长度 单位msL6=0.370;%O2O3的长度 单位msL6u=0.6572;%O3D的长度 单位m%求解S3、Theta3、Theta4和SE四个变量sS4=(sL6)2+(sL1)2-2*sL6*sL1*cos(stheta1+pi/2).0.5; %求出O3A的值for i=1:120 %求解角度theta3、Theta4和SE

33、的长度 stheta4(i)=acos(sL1*cos(stheta1(i)/sS4(i); stheta5(i)=asin(sL6u-sL4*sin(stheta4(i)/sL5); sSE(i)=sL4*cos(stheta4(i)+sL5*cos(stheta5(i);end%求解完成%求解完成%求解VS3、W3、W4和VE四个变量for i=1:120 J= inv(cos(stheta4(i),-sS4(i)*sin(stheta4(i),0,0; sin(stheta4(i),sS4(i)*cos(stheta4(i),0,0; 0,-sL4*sin(stheta4(i),-sL5

34、*sin(stheta5(i),-1; 0,sL4*cos(stheta4(i),sL5*cos(stheta5(i),0); K=J*W1*-sL1*sin(stheta1(i);sL1*cos(stheta1(i);0;0; sVS4(i)=K(1); sW4(i)=K(2); sW5(i)=K(3); sVE(i)=K(4);end%求解完成%求解aS3、a3、a4、aE四个变量for i=1:120 J= inv(cos(stheta4(i),-sS4(i)*sin(stheta4(i),0,0; sin(stheta4(i),sS4(i)*cos(stheta4(i),0,0; 0,

35、-sL4*sin(stheta4(i),-sL5*sin(stheta5(i),-1; 0,sL4*cos(stheta4(i),sL5*cos(stheta5(i),0); P=W1*W1*-sL1*cos(stheta1(i);-sL1*sin(stheta1(i);0;0; M=-sW4(i)*sin(stheta4(i),-sVS4(i)*sin(stheta4(i)-sS4(i)*sW4(i)*cos(stheta4(i),0,0; sW4(i)*cos(stheta4(i),sVS4(i)*cos(stheta4(i)-sS4(i)*sW4(i)*sin(stheta4(i),0,

36、0; 0,-sL4*sW4(i)*cos(stheta4(i),-sL5*sW5(i)*cos(stheta5(i),0; 0,-sL4*sW4(i)*sin(stheta4(i),-sL5*sW5(i)*sin(stheta5(i),0; N=sVS4(i);sW4(i);sW5(i);sVE(i); K=J*(-M*N+P); saS3(i)=K(1); sa4(i)=K(2); sa5(i)=K(3); saE(i)=K(4);end%求解完成%画图%画运动图figure(1);plot(theta1,S3,g*,stheta1,-sSE,r+);grid on;xlabel(1(rad

37、);ylabel(2(rad);title(方案一刨头位移(红色)和方案二刨头位移(绿色);axis(theta1(1) ,theta1(120),-0.6,0.7);figure(2);plot(theta1,V3,go,stheta1,-sVE,r.);hold on;grid on;xlabel(1(rad);ylabel(2(rad/s);title(方案一刨头速度(绿色)和方案二刨头速度(红色);axis(theta1(1) ,theta1(120),-7,3);figure(3);plot(theta1,a3,g*,stheta1,-saE,r.);hold on;grid on;xlabel(1(rad);ylabel(a2(rad/s)、a3(m/s);title(方案一刨头速度(红色)和方案二刨头速度(绿色);axis(theta1(1) ,theta1(120),-80,80);figure(4);plot(theta1,Mp,g*);xlabel(