matlab曲柄连杆机构分析(共5页)

1、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(

2、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;%下面黄金分割法求滑块的速度与加速度最大值epsi

3、lon=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)=epsilon if y1=epsilon if y1=epsilo

4、n if y1=y2 b=a2;a2=a1;y2=y1;a1=b-0.618*(b-a); y1=-(-R.*omiga.2.*cos(omiga.*a1)-(R.*omiga.*cos(omiga.*a1).2./ (l.*(cos(alpha1).3)-R.*omiga.2.*sin(omiga.*a1).*tan(alpha1); else a=a1;a1=a2;y1=y2;a2=a+0.618*(b-a); y2=-(-R.*omiga.2.*cos(omiga.*a2)-(R.*omiga.*cos(omiga.*a2).2./ (l.*(cos(alpha1).3)-R.*omig

5、a.2.*sin(omiga.*a2).*tan(alpha1); end n1=n1+1;enddisp(经过,num2str(n1),次计算,找到的加速度最大值对应的w是:,num2str(omiga*(a+b)/2)disp(速度最大值点与加速度最大值相差的wt是:,num2str(omiga*(a+b)/2)-vbm1),弧度。)曲柄连杆机构动力学分析图clear;clc;close all;n=750;h=0.2;f=0.15;g=9.8;l=0.975;m2=3.8;R=0.0381;m1=2.1;m3=2.1*2.3;omega=n*pi/30;tmax=2*pi/omega;t

6、=linspace(0,tmax,100);%计算曲柄转一周的总t值for n=1:100 alpha1(n)=atan(h+R.*sin(omega.*t(n)./sqrt(l.*l-(h+R.*sin(omega.*t(n)+pi; alpha1p(n)=-(R.*omega.*cos(omega*t(n)./(l.*cos(alpha1(n);%alpha1的一阶导数alpha1pp(n)=(R.*omega.2.*sin(omega.*t(n)+l.*sin(alpha1(n).*alpha1p(n).2)./ (l.*cos(alpha1(n);%alpha1的二阶导数 xA=R.*

7、cos(omega*t(n); yA=R.*sin(omega*t(n); xApp=-R.*omega.2.*cos(omega.*t(n); yApp=-R.*omega.2.*sin(omega.*t(n); %曲柄末端A的加速度 x1=xA/2; y1=yA/2; x1pp=xApp./2;y1pp=yApp./2; %曲柄质心的加速度 xBpp=xApp-l*(alpha1pp(n).*sin(alpha1(n)+alpha1p(n).2.*cos(alpha1(n); yBpp=0; %刀头(滑块)B的加速度 x2=xA+l.*cos(alpha1(n)./2;y2=yA+l.*s

8、in(alpha1(n)./2; x2pp=xApp-l*(alpha1pp(n).*sin(alpha1(n)+alpha1p(n).2.*cos(alpha1(n)./2; y2pp=yApp+l*(alpha1pp(n).*cos(alpha1(n)-alpha1p(n).2.*sin(alpha1(n)./2; %连杆质心的加速度%对于曲柄m1有方程% Ox-m1.*x1pp+Ax=0% Oy-m1.*y1pp+Ay-m1.*g=0% M-Ay.*xA-Ax.*yA+m1.*g.*x1+m1.*y1pp.*x1+m1.*x1pp*y1=0%对于连杆m2有方程% -Ax+Bn.*cos(

9、alpha1)+Bt.*cos(alpha1+pi/2)-m2.*x2pp=0% -Ay+Bn.*sin(alpha1)+Bt.*sin(alpha1+pi/2)-m2.*y2pp-m2.*g=0% Bt.*l+m2.*(g+y2pp).*(x2-xA)+m2.*x2pp.*(y2-yA)=0%对于刀头(滑块)m3有方程% -f.*N+Bn.*cos(alpha1-pi)-Bt.*cos(alpha1-pi/2)-m3.*xBpp=0% N-m3.*g+Bn.*sin(alpha1-pi)+Bt.*sin(alpha1-pi/2)=0% 8个方程，未知数也有Ox,Oy,Ax,Ay,Bn,Bt,

10、N,M共8个A=1,0,1,0,0,0,0,0; 0,1,0,1,0,0,0,0; 0,0,-yA,-xA,0,0,0,1; 0,0,-1,0,cos(alpha1(n),cos(alpha1(n)+pi/2),0,0; 0,0,0,-1,sin(alpha1(n),sin(alpha1(n)+pi/2),0,0; 0,0,0,0,0,l,0,0; 0,0,0,0,cos(alpha1(n)-pi),-cos(alpha1(n)-pi/2),-f,0; 0,0,0,0,sin(alpha1(n)-pi),sin(alpha1(n)-pi/2),1,0;B=m1.*x1pp; m1.*y1pp+

11、m1.*g; -m1.*g.*x1-m1.*y1pp.*x1-m1.*x1pp*y1; m2.*x2pp; m2.*y2pp+m2.*g; -m2.*(g+y2pp).*(x2-xA)-m2.*x2pp.*(y2-yA); m3.*xBpp; m3.*g; X=AB; %X向量的顺序为Ox,Oy,Ax,Ay,Bn,Bt,N,M XT(n,:)=X; %将N转置后以行向量形式存入矩阵XT中 Bn(n)=X(5);Bt(n)=X(6);N(n)=X(7); %为求滑块最大惯性力Qm从XT中取出的Bn、Bt及N值 Q(n)=Bn(n)*cos(alpha1(n)-pi)-Bt(n)*cos(alpha1(n)-pi/2)-f*N(n);endfor k=1:8 subplot(3,3,k);plot(t,XT(:,k);grid on;endsubpl