MATLAB的曲柄滑块和四杆机构的综合设计解析
- 格式:docx
- 大小:237.74 KB
- 文档页数:10
基于matlab的四杆机构运动分析一、四杆机构基本概念四杆机构是一种通过变换连杆长度,改变机构运动形态的机械系统。
四杆机构通常由固定连杆、推动连杆、连接杆和工作连杆四个连杆组成,其中固定连杆和推动连杆固定不动,连接杆和工作连杆则沿固定轴线的方向做平动或旋转运动。
四杆机构的基本构造如下图所示:四杆机构的四个连杆的长度和构造参数,以及驱动机构的运动决定了机构的运动特性。
在进行四杆机构运动分析时,需要通过求解运动学关系式和动力学方程,得到连杆的运动规律和力学特性。
二、四杆机构运动学分析1.运动学基本方程四杆机构的运动学分析基本方程是连杆长度变化的定理,即:l₁²+l₂²-2l₁l₂cosθ₂=l₃²+l₄²-2l₃l₄cosθ₄其中,l₁,l₂分别为固定连杆和推动连杆长度;l₃,l₄分别为连接杆和工作连杆长度;θ₂,θ₄分别为推动连杆和工作连杆的夹角。
2.运动学求解方法根据四杆机构运动学基本方程,可以求解机构中任意连杆的角度和位置,从而分析机构运动规律。
在matlab程序中,运动分析可以采用分析法或图解法。
分析法通常采用向量法或坐标法,即将四杆机构中各连杆和运动副的运动量表示为向量或坐标,然后根据连杆长度变化的定理,求解四个未知角度θ₁、θ₂、θ₃、θ₄。
图解法则先通过画图确定机构的运动规律,在图上求解连杆的角度。
比如可以采用伯格(Bourgeois)图法或恰普利恩(Chaplygin)图法等。
四杆机构动力学分析基本方程包括平衡方程和力平衡方程。
平衡方程:当四杆机构处于平衡状态时,连杆的受力关系可以表示为:ΣF=0其中ΣF为各连杆受力的合力。
ΣF=m×a其中,m为每个连杆的质量,a为连杆的加速度。
四杆机构动力学求解方法以matlab为工具,可借助matlab的求解器完成求解。
具体可以利用matlab的优化工具箱、控制工具箱和系统动态学工具箱等,来实现机构模型的动态模拟、仿真和优化设计。
石河子大学毕业设计(论文)题目:基于MATLAB的四杆机构运动分析与动画模拟系统院(系):机械电气工程学院专业:机械设计制造及其自动化学号: 2002071189姓名: 娄元建指导教师:葛建兵完成日期:二零零六年五月基于MATLAB的四杆机构运动分析与动画模拟系统[摘要] 本文介绍MATLAB开发机构运动分析和动画模拟系统的方法,并且利用MATLAB软件实现平面四杆机构的运动仿真。
以MATLAB程序设计语言为平台,将参数化设计与交互式相结合,设计出四杆机构仿真系统,能够实现四杆机构的参数化设计,并且能够进行机构的速度和加速度分析。
系统具有方便用户的良好界面,并给出界面设计程序,从而使机构分析更加方便、快捷、直观和形象,设计者只需输几参数就可得到仿真结果,为平面四杆机构的设计与分析提供一条便捷的途径。
[关键词] 机构;运动分析;动画模拟;仿真;参数化;MATLABAbstract:The kinematical analysis and animation method of the mechanism using MATLAB was discussed in the paper , and the kinematic simulation of planar four-bar mechanism with software MATLAB . And emulational system was developed , the system adopted Matlab as a design , It combined parametic design with interactive design and had good interface for user , that can realize parametic design of four-bar mechanism , also to make real speed and acceleration of mechanism 。
m a t l a b曲柄连杆机构分析clear;clc;n=750;l=0.975;R=0.0381;h=0.2;omiga=n.*pi/30;tmax=2.*pi/omiga;t=0:0.001:tmax; %计算曲柄转一圈的总t值alpha1=atan((h+R.*sin(omiga.*t))./sqrt(l.*l-(h+R.*sin(omiga.*t))))+pi;alpha1p=-(R.*omiga.*cos(omiga.*t))./(l.*cos(alpha1));vb=-R.*omiga.*sin(omiga.*t)+R.*omiga.*cos(omiga.*t).*tan(alpha1);ab=-R.*omiga.^2.*cos(omiga.*t)-(R.*omiga.*cos(omiga.*t)).^2./(l.*(cos(alpha1)).^3)-R.*omiga.^2.*sin(omiga.*t).*tan(alpha1);subplot(1,2,1);plot(t,vb);title('曲柄滑块机构的滑块v-t图');xlabel('时间t(曲柄旋转一周)');ylabel('滑块速度v');grid on;subplot(1,2,2);plot(t,ab);title('曲柄滑块机构的滑块a-t图');xlabel('时间t(曲柄旋转一周)');ylabel('滑块加速度a');grid on;%下面黄金分割法求滑块的速度与加速度最大值epsilon=input('根据曲线初始区间已确定,请输入计算精度epsilon(如输入0.001):');a=0;b=0.04; %初始区间n1=0; %n1用于计算次数a1=b-0.618*(b-a);y1=-R.*omiga.*sin(omiga.*a1)+R.*omiga.*cos(omiga.*a1).*tan(alpha1);a2=a+0.618*(b-a);y2=-R.*omiga.*sin(omiga.*a2)+R.*omiga.*cos(omiga.*a2).*tan(alpha1);while abs(a-b)>=epsilonif y1<=y2b=a2;a2=a1;y2=y1;a1=b-0.618*(b-a);y1=-R.*omiga.*sin(omiga.*a1)+R.*omiga.*cos(omiga.*a1).*tan(alpha1);elsea=a1;a1=a2;y1=y2;a2=a+0.618*(b-a);y2=-R.*omiga.*sin(omiga.*a2)+R.*omiga.*cos(omiga.*a2).*tan(alpha1);endn1=n1+1;endvbm1=omiga*(a+b)/2;disp(['经过',num2str(n1),'次计算,用黄金分割法找到速度最大值对应的wt是:', num2str(vbm1),'弧度。
基于matlab的平面四连杆机构设计以及该机构的运动仿真分析摘要四连杆机构因其结构方便灵活,能够传递动力并实现多种运动形式而被广泛应用于各个领域,因此对其进行运动分析具有重要的意义。
传统的分析方法主要应用几何综合法和解析综合法,几何综合法简单直观,但是精确度较低;解析法精确度较高,但是计算工作量大。
随着计算机辅助数值解法的发展,特别是MATLAB软件的引入,解析法已经得到了广泛的应用。
对于四连杆的运动分析,若应用MATLAB 则需要大量的编程,因此我们引入proe软件,我们不仅可以在此软件中建立实物图,而且还可以对其进行运动仿真并对其运动分析。
在设计四连杆时,我们利用解析综合法建立数学模型,再根据数学模型在MATLAB中编程可以求得其他杆件的长度。
针对范例中所求得的各连杆的长度,我们在proe软件中画出其三维图(如图4)并在proe软件中进行仿真分析得出CB,的角加速度的变化,从而得到CB,两接触处所受到的力是成周期性变化的,可以看出CB,两点处的疲劳断裂,我们提B,两点处极易疲劳断裂,针对C出了在设计四连杆中的一些建议。
关键字:解析法 MATLAB 软件 proe 软件 运动仿真建立用解析法设计平面四杆机构模型对于问题中所给出的连架杆AB 的三个位置与连架杆CD 的三个位置相对应,即三组对应位置为:332211,,,,,ψϕψϕψϕ,其中他们对应的值分别为: 52,45,82,90,112,135,为了便于写代数式,可作出AB 与CD 对应的关系,其图如下:图—2 AB 与CD 三个位置对应的关系通过上图我们可以通过建立平面直角坐标系并利用解析法来求解,其直角坐标系图如下:φααi θi φi图—3 平面机构直角坐标系通过建立直角坐标系OXY ,如上图所示,其中0α与0φ为AB 杆与CD 杆的初始角,各杆件的长度分别用矢量d c b a ,,,,表示,将各矢量分别在X 轴与Y 轴上投影的方程为⎩⎨⎧=++=+)sin(*)sin(*)sin(*)cos(*)cos(*)cos(*φθαφθαc b a c d b a在上述的方程中我们可以消除θ,从而可以得到α与φ之间的关系如下:)cos(2)cos(2)cos(2)(2222αφαφab ac cd b d c a +-=+-++ (1) 为便于化简以及matlab 编程我们可以令:⎪⎪⎪⎩⎪⎪⎪⎨⎧==-++=c d H a d H ac b d c a H 32222212 (2) 通过将(2)式代入(1)式中则可以化简得到如下等式: )cos()cos()cos(321αφαφH H H +-=+ (3)我们可以通过(3)式将两连架杆对应的位置带入(3)式中,我们可以得到如下方程:⎪⎩⎪⎨⎧+-=++-=++-=+)cos()cos()cos()cos()cos()cos()cos()cos()cos(333332123222211311121ϕψϕψϕψϕψϕψϕψH H H H H H H H H (4) 联立(4)方程组我们可以求得321,,H H H ,再根据(2)中的条件以及所给定的机架d 的长度,我们可以求出其它杆件的长度为:⎪⎪⎪⎩⎪⎪⎪⎨⎧-++===1222322acH d c a b H d c H d a (5)四连杆设计范例:在日常生活中,我们经常看到消防门总能自动关上,其实它是利用四连杆机构与弹簧组成的。
MATLAB 解题1.设有如图所示四杆机构,其中→R 4为机架(常矢),→R1为主动杆,→R3为从动杆,→R 2为连杆。
设在某一工作位置时各杆的角速度和角加速度分别取如下值:ω1=20 rad/s, ε1= 0;ω2=8.5 rad/s, ε2=-10 rad /s 2;ω3=13 rad/s, ε3=-160rad /s 2.试根据上述要求确定该机构尺寸比。
根据图(2),回路闭合方程可写为:→R 1 +→R 2 +→R 3=-→R 4 回路闭合方程对时间求导一次,利用(6)式,可得: 图2 ω1→R 1 +ω2→R 2 +ω3→R 3 = 0回路闭合方程对时间求导两次,利用(7)式,可得c 1→R 1 + c 2 →R 2 + c 3→R 3 = 0其中 c 1=ε1+j ω12 , c 2=ε2+j ω22, c 3=ε3+j ω32解关于→R 1 ,→R 2 和→R 3的线性方程组:⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡-=⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎣⎡⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡→→→001111321321321R R R c c c ωωω→R 4 (13) 可得 →R 1=DDx →R 4, →R 2=DDy →R 4 , →R 3=DDz →R 4注意到上述解中含有相同的分母D,它是一个复数,不妨记为D =k<j α|,被它除的效果是把各杆的长度都缩小k 倍,同时方向都顺时针旋转α角,相当于机构不动,坐标轴逆时针旋转α角。
设计机构时,重要的是机构的形状与尺寸比例。
基于这种考虑,可设→R 4 / D =1,则有→R 1=D x =32320111c c ωω-=1230-j497.3 ; →R 2= D y =311030111c c ωω-=-3200-j1820 ; →R 3= D z =001112121c c ωω-=200+j1955 . 于是:→R 4 = -(→R 1 +→R 2+→R 3) = 1770+j362.3在坐标系上作出上述各杆矢量图,根据各杆矢量图作出机构的闭合矢量图,再根据实际需要选定某一杆长度,其它各杆长度按图比例相似放大。
子函数%子函数slider_crank文件function[theta2,s3,omega2,v3,alpha2,a3]=slider_crank(theta1,omega 1,alpha1,l1,l2,e)%计算连杆2的角位移和滑块3的线位移theta2=asin((e-l1*sin(theta1))/l2);s3=l1*cos(theta1)+l2*cos(theta2);%计算连杆2的角为速度和滑块的线速度A=[-l1*sin(theta1),1;-2*cos(theta2),0];B=[-l1*sin(theta1);l1*cos(theta1)];omega=A\(omega1*B);omega2=omega(1);v3=omega(2);%计算连杆2的角加速度和滑块3的线加速度At=[omega2*l2*cos(theta2),0;omega2*l2*sin(theta2),0];Bt=[-omega1*l1*cos(theta1);-omega1*l1*sin(theta1)];alpha=A\(-At*omega+alpha1*B+omega1*Bt);alpha2=alpha(1);a3=alpha(2);主函数%住程序slider_crank_main文件%输入已经知道的数据clear;l1=100;l2=300;e=0;hd=pi/180;du=180/pi;omega1=10;alpha1=0;%调用子函数slider_ank计算曲柄滑块机构位移,速度,加速度for n1=1:720theta1(n1)=(n1-1)*hd;[theta2(n1),s3(n1),omega2(n1),v3(n1),alpha2(n1),a3(n1)]=slider_cr ank...(theta1(n1),omega1,alpha1,l1,l2,e);end%位移,速度,加速度和曲柄滑块机构图形输出figure(l1);n1=1:720;subplot(2,2,1); %绘制位移图[AX,H1,H2]=plotyy(theta1*du,theta2*du,theta1*du,s3);set(get(AX(1),'ylabel'),'String','连杆角位移/\circ')set(get(AX(2),'ylabel'),'String','滑块位移/mm')title('位移线图');xlabel('曲柄转角\theta_1/\circ')grid on;subplot(2,2,2); %绘制速度图[AX,H1,H2]=plotyy(theta1*du,omega2,theta1*du,v3);set(get(AX(2),'ylabel'),'String','滑块速度/mm\cdots^{-1}') title('速度线图');xlabel('曲柄转角\theta_1/\circ')ylabel('连杆角速度/rad\cdots^{-1}')grid on;subplot(2,2,3); %绘制加速度图[AX,H1,H2]=plotyy(theta1*du,alpha2,theta1*du,a3);set(get(AX(2),'ylabel'),'String','滑块加速度/mm\cdots^{-2}') title('加速度线图');xlabel('曲柄转角\theta_1/\circ')ylabel('连杆加速度/rad\cdots^{-2}')grid on;subplot(2,2,4);%绘曲柄滑块机构图x(1)=0;y(1)=0;x(2)=l1*cos(70*hd);y(2)=l1*sin(70*hd);x(3)=s3(70);y(3)=e;x(4)=s3(70);y(4)=0;x(5)=0;y(5)=0;x(6)=x(3)-40;y(6)=y(3)+10;x(7)=x(3)+40;y(7)=y(3)+10;x(8)=x(3)+40;y(8)=y(3)-10;x(9)=x(3)-40;y(9)=y(3)-10;x(10)=x(3)-40;y(10)=y(3)+10;i=1:5;plot(x(i),y(i));grid on;hold on;i=6:10;plot(x(i),y(i));title('曲柄滑块机构');grid on;hold on;xlabel('mm');ylabel('mm')axis([-50 400 -20 130]); plot(x(1),y(1),'o');plot(x(2),y(2),'o');plot(x(3),y(3),'o');%曲柄滑块的仿真运动figure(2)m=moviein(20);j=0;for n1=1:5:360j=j+1;clf;%x(1)=0;y(1)=0;x(2)=l1*cos(n1*hd); y(2)=l1*sin(n1*hd); x(3)=s3(n1);y(3)=e;x(4)=(l1+l2+50);y(4)=0;x(5)=0;y(5)=0;x(6)=x(3)-40;y(6)=y(3)+10;x(7)=x(3)+40;y(7)=y(3)+10;x(8)=x(3)+40;y(8)=y(3)-10;x(9)=x(3)-40;y(9)=y(3)-10;x(10)=x(3)-40;y(10)=y(3)+10;%i=1:3;plot(x(i),y(i));grid on;hold on;i=4:5;plot(x(i),y(i));i=6:10;plot(x(i),y(i));plot(x(1),y(1),'o');plot(x(2),y(2),'o');plot(x(3),y(3),'o');xlabel('mm');ylabel('mm')axis([-50 450 -150 150]); m(j)=getframe;endmovie(m)。
matlab曲柄滑块机构课程设计一、课程目标知识目标:1. 理解并掌握曲柄滑块机构的基本原理与运动特性;2. 学会使用MATLAB软件进行曲柄滑块机构的运动仿真;3. 掌握运用MATLAB分析曲柄滑块机构的运动数据及性能参数。
技能目标:1. 能够运用所学知识,设计简单的曲柄滑块机构;2. 熟练操作MATLAB软件,进行曲柄滑块机构的运动分析与仿真;3. 能够通过MATLAB处理数据,优化曲柄滑块机构的设计。
情感态度价值观目标:1. 培养学生的团队协作精神,提高沟通与表达能力;2. 激发学生探索科学、技术问题的兴趣,培养创新意识和实践能力;3. 增强学生对机械工程领域的认识和热爱,提高专业认同感。
课程性质:本课程为机械工程专业课程设计,旨在通过实践操作,使学生掌握曲柄滑块机构的设计与分析方法。
学生特点:学生已具备一定的机械原理、力学和MATLAB基础,具有较强的动手能力和学习兴趣。
教学要求:结合实际工程案例,以实践为主,注重培养学生的实际操作能力、分析问题和解决问题的能力。
通过课程学习,使学生能够独立完成曲柄滑块机构的设计与分析任务。
二、教学内容1. 曲柄滑块机构原理及运动特性分析:- 曲柄滑块机构基本组成与工作原理- 曲柄滑块机构的运动学分析- 运动方程的建立及求解2. MATLAB软件在曲柄滑块机构中的应用:- MATLAB软件的基本操作与常用函数- MATLAB曲线拟合、数值计算等功能在曲柄滑块机构分析中的应用- MATLAB/Simulink环境下曲柄滑块机构的运动仿真3. 曲柄滑块机构设计及优化:- 设计原则与步骤- 参数化设计方法- 基于MATLAB的曲柄滑块机构设计优化4. 实践操作与案例分析:- 实际工程案例介绍与分析- 曲柄滑块机构设计及运动分析的实践操作- 数据处理与结果分析教学内容安排与进度:1. 第一周:曲柄滑块机构原理及运动特性分析2. 第二周:MATLAB软件在曲柄滑块机构中的应用3. 第三周:曲柄滑块机构设计及优化4. 第四周:实践操作与案例分析教材章节:1. 《机械原理》中曲柄滑块机构相关章节2. 《MATLAB基础与应用》中相关章节3. 《机械设计》中机构设计及优化相关章节教学内容注重理论与实践相结合,通过系统性的教学,使学生掌握曲柄滑块机构的设计与分析方法,并能够运用MATLAB软件进行实际操作。
KUNMING UNIVERSITY OF SCIENCE AND TECHNOLOGY《计算机仿真技术》课程设计报告冯叶/ 浦合旳201410302544/ 201410302547刘孝保2015年6月姓名: 学号: 专业班级: 指导教师:机械卓目录©区肌理乂殳申KUMWBG sngn OF SCIENCE MO TCCWlOGr目录1 •仿真问题描述.........................................................................2•仿真问题数学模型......................................................................3. Mat lab实现方法 .....................................................................4・Mat lab代码..........................................................................5•仿真结论..............................................................................6.遇到的问题和解决的方式.................................................................7 •课程学习意见与建议...................................................................《计算机仿真技术》课程设计报告 艮咽疗N 乂孝 ItnVH ; WmJSTY :f SCOtCE MP TOCtXCCf 1 •仿真问题描述已知机架AD 长为L1,曲柄AB 长为L2,连杆BC 长L3,另一机架长CD 长为L4,与AB 杆相 连的是一滑块E 。
《计算机仿真技术》课程设计报告姓名:冯叶 / 浦合昀学号: 201410302544/ 201410302547专业班级:机械卓越141 指导教师:刘孝保2015年 6月目录目录1.仿真问题描述 ...................................................................................................................................................2.仿真问题数学模型 ...........................................................................................................................................3.Matlab实现方法 ..............................................................................................................................................4.Matlab代码 ......................................................................................................................................................5.仿真结论 ...........................................................................................................................................................6.遇到的问题和解决的方式 ...............................................................................................................................7.课程学习意见与建议 .......................................................................................................................................1.仿真问题描述已知机架AD 长为L1,曲柄AB 长为L2,连杆BC 长L3,另一机架长CD 长为L4,与AB 杆相连的是一滑块E 。
BE 杆长为L5,设计一个四杆加滑块的机构,其中L1-L5杆长可变。
并且可以通过输入的杆长,来判别,该机构到底可不可行。
L3L4L2 L5L12.仿真问题数学模型(1)四杆机构的设计:在用矢量法建立机构的位置方程时,需将构件用矢量来表示,并作出机构的封闭矢量多边形。
如图1所示,先建立一直角坐标系。
设各构件的长度分别为1L 、2L 、3L 、4L ,其方位角为1θ、2θ、3θ、 4θ。
以各杆矢量组成一个封闭矢量多边形,即ABCDA 。
其个矢量之和必等于零。
易知:2314L L L L +=+角位移方程的分量形式为:2233114422331144cos cos cos cos sin sin sin sin L L L L L L L L θθθθθθθθ+=+⎧⎫⎨⎬+=+⎩⎭ 要求th3,那么 ()'''sin cos cos cos sin sin d d dt dt d d dt dt uv vu uv θθθωθθθθωθ⎧⎫==⎪⎪⎪⎪⎪⎪=-=-⎨⎬⎪⎪⎪⎪=+⎪⎪⎩⎭在角位移方程分量形式中,由于假定机架为参考系,矢量1与x 轴重合,1θ=0,则有非线性超越方程组:1342233144234223344(,)cos cos cos 0(,)sin sin sin 0f L L L L f L L L θθθθθθθθθθ=+--=⎧⎫⎨⎬=+-=⎩⎭可以借助牛顿-辛普森数值解法或Matlab 自带的fsolve 函数求出连杆3的角位移3θ和摇杆4的角位移4θ。
A B C E D求解具有n 个未知量i x (i=1,2,…,n )的线性方程组:111122112111221211221n n n n n n nn n n a x a x a x b a x a x a x b a x a x a x b +++=⎧⎫⎪⎪+++=⎪⎪⎨⎬⎪⎪⎪⎪+++=⎩⎭式中,系列矩阵A 是一个*n n 阶方阵:1111n m mn a a A a a ⎛⎫ ⎪= ⎪ ⎪⎝⎭A 的逆矩阵为1A -;常数项b 是一个n 维矢量:12(,,,)T n b b b b =因此,线性方程组解的矢量为:12(,,,)T T n x x x x A b == 非线性超越方程组是求解连杆3和摇杆4角速度和角加速度的依据。
(2)曲柄滑块的设计:由图可知,C 滑块的位移总是与AB,BC 和他们之间的角度存在着一定的关系,关系如下:LAC=AB ×cos (∅)+√BC 2—AB 2×cos∅2通过以上这个式子,我们就可以来求C 点的位移,速度,加速度。
3.Matlab 实现方法(1)怎么设计四杆机构:创建函数FoutBarPosition ,函数fsolve 通过他确定34,θθ,然后知道34,θθ后,来求取各个点的坐标,通过plot 命令在指定的区域内连线,取点,画图。
(2)怎么设计曲柄滑块机构:通过解方程的方法,用solve 来求取C 滑块的坐标,用diff 函数求取C 滑块的速度,加速度曲线,通过plot 命令在指定的区域内连线,取点,画图。
4.Matlab 代码(1)建新的函数在点m文件中:function t=fourbarposition(th,th2,L2,L3,L4,L1)t=[L2*cos(th2)+L3*cos(th(1))-L4*cos(th(2))-L1;…L2*sin(th2)+L3*sin(th(1))-L4*sin(th(2))];(2)主程序如下:%获取杆长l1=str2double(get(handles.edit1,'string'));l2=str2double(get(handles.edit2,'string'));l3=str2double(get(handles.edit3,'string'));l4=str2double(get(handles.edit4,'string'));l5=str2double(get(handles.edit8,'string'));%滑块和四杆机构的设计syms t s; %定义变量f=l5^2-l2^2-s^2+2*l2*s*cos(t);ff=solve(f,s);vv=diff(ff,1);aa=diff(ff,2);th2=0:pi/15:6*pi;times=length(th2);for i=1:91wyy(1,i)=eval(subs(ff(1),t,th2(i)));wyy(2,i)=eval(subs(ff(2),t,th2(i)));vyy(1,i)=eval(subs(vv(1),t,th2(i)));vyy(2,i)=eval(subs(vv(2),t,th2(i)));ayy(1,i)=eval(subs(aa(1),t,th2(i)));ayy(2,i)=eval(subs(aa(2),t,th2(i)));endfor i=1:timesif wyy(1,i)>0wy(i)=wyy(1,i);elsewy(i)=wyy(2,i);endif vyy(1,i)>0vy(i)=vyy(1,i);elsevy(i)=vyy(2,i);endif ayy(1,i)>0ay(i)=ayy(1,i);elseay(i)=ayy(2,i);endendth34=zeros(length(th2),2); %%建立一个N行2列的零矩阵options=optimset('display','off');for m=1:length(th2) %用fsove函数求解关于th3,th4的非线性超越方程,结果保存在th34中th34(m,:)=fsolve('fourbarposition',[1 1],...options,th2(m),l2,l3,l4,l1);end%求各个的坐标Ex=wy;Ey=zeros(size(th2));Cy=l2*sin(th2)+l3*sin(th34(:,1)');Cx=l2*cos(th2)+l3*cos(th34(:,1)');Bx=[l2*cos(th2)];By=[l2*sin(th2)];Ax=zeros(size(th2));Ay=zeros(size(th2));Dx=l1+zeros(size(th2));Dy=zeros(size(th2));Ev=vy;Ew=zeros(size(th2));Ea=ay;En=zeros(size(th2));%求位移,速度,加速度的范围:g=[Ax Bx Cx Dx Ex];m=[Ay By Cy Dy Ey];maxX=max(g);minX=min(g);maxY=max(m);minY=min(m);maxwy=max(Ex);minwy=min(Ex);maxvy=max(Ev);minvy=min(Ev);maxay=max(Ea+50);minay=min(Ea-50);%画动画图for i=1:timesaxes(handles.axes1);plot([Ax(i),Bx(i)],[Ay(i),By(i)],'lineWidth',3); hold onplot([Bx(i),Cx(i)],[By(i),Cy(i)],'lineWidth',3); plot([Ax(i),Dx(i)],[Ay(i),Dy(i)],'lineWidth',3); plot([Cx(i),Dx(i)],[Cy(i),Dy(i)],'lineWidth',3); plot([Bx(i),Ex(i)],[By(i),Ey(i)],'lineWidth',3); plot([-10000,Ax(i)],[0,Ay(i)],'lineWidth',3);plot([10000,Ax(i)],[0,Ay(i)],'lineWidth',3);plot([Ex(i)+10,Ex(i)-10],[Ey(i)+10,Ey(i)+10]);plot([Ex(i)+10,Ex(i)-10],[Ey(i)-10,Ey(i)-10]);plot([Ex(i)-10,Ex(i)-10],[Ey(i)+10,Ey(i)-10]);plot([Ex(i)+10,Ex(i)+10],[Ey(i)+10,Ey(i)-10]);plot(0,0,'or','lineWidth',3)plot(Bx(i),By(i),'or','lineWidth',3)plot(Cx(i),Cy(i),'or','lineWidth',3)plot(Dx(i),Dy(i),'or','lineWidth',3)plot(Ex(i),Ey(i),'or','lineWidth',3)axis equal;axis([minX,maxX,minY,maxY]);axis off;hold off;pause(0.1)%画滑块位移图axes(handles.axes3);plot(th2(1:i),wy(1:i),'lineWidth',3);axis([0,2*pi,minwy,maxwy])pause(0.1)%画滑块速度图axes(handles.axes4);plot(th2(1:i),vy(1:i),'lineWidth',3);axis([0,2*pi,minvy,maxvy])pause(0.1)%画滑块加速度图axes(handles.axes5);plot(th2(1:i),ay(1:i),'lineWidth',3);axis([0,2*pi,minay,maxay])pause(0.1)end%判断杆长的代码:l1=str2double(get(handles.edit1,'string'));l2=str2double(get(handles.edit2,'string'));l3=str2double(get(handles.edit3,'string'));l4=str2double(get(handles.edit4,'string'));lall=l1+l2+l3+l4;lmax=max([l1 l2 l3 l4]);lmin=min([l1 l2 l3 l4]);if (lmax+lmin)>(lall-lmax-lmin)set(handles.edit5,'string','不符合杆长条件,请重新输入');set(handles.edit6,'string','');set(handles.edit7,'string','');elseif l2>l4set(handles.edit6,'string','不符合A点的周转条件,请重新输入'); set(handles.edit5,'string','');set(handles.edit7,'string','');elseset(handles.edit7,'string','·符合条件,可以运动'); set(handles.edit5,'string','');set(handles.edit6,'string','');endend%关闭代码:close5.仿真结论(1)如图所示为界面。