机械原理编程

function data=analy()
%杆件上某点的轨迹，速度，加速度图像
%用解析法求解

l1=input('请输入 l1(mm): '); %输入各项参数
l2=input('请输入 l2(mm): ');
l3=input('请输入 l3(mm): ');
lbp=input('请输入 lbp(mm): ');
ang=input('请输入 ang(。): ');
dx=input('请输入 dx(mm): ');
dy=input('请输入 dy(mm): ');
n=input('请输入 n(。/s): ');
a=input('请输入 a(mm/s2): ');

ang=ang*pi/180;
ang1=0:pi/18:2*pi;
ag=0:10:360;
l1=l1.*ones(size(ang1));
l2=l2.*ones(size(ang1));
l3=l3.*ones(size(ang1));
lbp=lbp.*ones(size(ang1));
dx=dx.*ones(size(ang1));
dy=dy.*ones(size(ang1));
ang=ang.*ones(size(ang1));
k1=dx-l1.*cos(ang1);
k2=dy-l1.*sin(ang1);
k3=(k1.^2+k2.^2+l3.^2-l2.^2)./(2.*l3);
y3=k2-sqrt(k1.^2+k2.^2-k3.^2);
x3=k1-k3;
[m3,n3]=size(y3);
for i=1:1:n3
if (x3(m3,i)>0&&y3(m3,i)>0)
ang3(m3,i)=atan(y3(m3,i)/x3(m3,i));
elseif (x3(m3,i)>0&&y3(m3,i)<0)
ang3(m3,i)=atan(y3(m3,i)/x3(m3,i))+2*pi;
elseif x3(m3,i)<0
ang3(m3,i)=atan(y3(m3,i)/x3(m3,i))+pi;
elseif (x3(m3,i)==0&&y3(m3,i)>0)
ang3(m3,i)=pi/2;
else ang3(m3,i)=3*pi/2;
end
end
ang3=2.*ang3;
for i=1:1:n3
if ang3(m3,i)>2*pi;
ang3(m3,i)=ang3(m3,i)-2*pi;
elseif ang3(m3,i)<0;
ang3(m3,i)=ang3(m3,i)+2*pi;
end
end

x2=k1+l3.*cos(ang3);
y2=k2+l3.*sin(ang3);
[m2,n2]=size(y2);
for i=1:1:n2
if (x2(m2,i)>0&&y2(m2,i)>0)
ang2(m2,i)=atan(y2(m2,i)/x2(m2,i));
elseif (x2(m2,i)>0&&y2(m2,i)<0)
ang2(m2,i)=atan(y2(m2,i)/x2(m2,i))+2*pi;
elseif x2(m2,i)<0
ang2(m2,i)=atan(y2(m2,i)/x2(m2,i))+pi;
elseif (x2(m2,i)==0&&y2(m2,i)>0)
ang2(m2,i)=pi/2;
else ang2(m2,i)=3*pi/2;
end
end

px=l1.*cos(ang1)+lbp.*cos(ang2+ang); %P点的轨迹图
py=l1.*sin(ang1)+lbp.*sin(ang2+ang);
figure(1);
plot(px,py,'k');
grid on;
axis tight;
xlabel('X /（mm）');
ylabel('Y /（mm）');
title('P点的轨迹');

w1=n*pi/180; %P点的速度图像
w1=w1*ones(size(ang2));
w3=w1.*(l1.*sin(ang1-ang2))./(l3.*sin(ang3-ang2));
w2=-w1.*(l1.*sin(ang1-ang3))./(l2.*sin(ang2-ang3));
vpx=-l1.*w1.*sin(ang1)-lbp.*w2.*sin(ang2+ang);
vpy=l1.*w1.*cos(ang1)+lbp.*w2.*cos(ang2+ang);
vp=sqrt(vpx.^2+vpy.^2);
vpx=vpx./1000;
vpy=vpy./1000;
vp=vp./1000;
figure(2);
plot(ag,vpx,'k:',ag,vpy,'k-',ag,vp,'k-.');
grid on;
axis tight;
xlabel('anger /（。）');
ylabel('v /（m/s）');
title('P点的速度');

a1=a*ones(size(ang2)); %P点的加速度图像
a2=(l2.*w2.^2+l1.*a1.*sin(ang1-ang2)+l1.*w1.^2.*cos(ang1-ang2)-l3.*w3.^2.*cos(ang3-ang2))./(l3.*sin(ang3-ang2));
a3=(l3.*w3.^2-l1.*a1.*sin(ang1-ang3)-l1.*w1.^2.*cos(ang1-ang3)-l2.*w2.^2.*cos(ang2-ang3))./(l2.*sin(ang2-ang3));
apx=-l1.*(a1.*sin(ang1)+w1.^2.*cos(ang1))-lbp.*(a2.*sin(ang2+ang)+w2.^2.*cos(ang2+ang));
apy=l1.*(a1.*cos(ang1)-w1.^2.*sin(ang1))+lbp.*(a2.*cos(ang2+ang)-w2.^2.*sin(ang2+ang));
ap=sqrt(apx.^2+apy.^2);
apx

=apx/1000;
apy=apy/1000;
ap=ap/1000;
figure(3);
plot(ag,apx,'k-.',ag,apy,'k:',ag,ap,'k-');
grid on;
axis tight;
xlabel('anger / (。)');
ylabel('a / (m/s2)');
title('P点的加速度 ');

ang1=ang1.*180./pi; %输出各项参数的计算值
ang2=ang2.*180./pi;
ang3=ang3.*180./pi;
data=[ang1' ang2' ang3' w2' w3' a2' a3' px' py' vpx' vpy' apx' apy'];
end