機動學第十三次作業

1.
開始時，假設要設計一組複式齒輪使其轉速比為125。其考慮的因素及過程如下：
i 決定組合數：複式齒列可用串聯多個組合之方式達到整體轉速比。設計總是由最簡單的結構開始，由於每組之轉速比以維持在10以內為佳，超過此值時則需考慮增加齒輪組數。而轉速比125很明顯可以看出開平方的結果比10大，而剛好125又是5的三次方，所以直接開立方。

ii 設驅動之小齒輪數最小為12齒(N2>=12)*，則依序可以得到對應大齒輪之齒數如下：

N3=5x12=60
N3=5x13=65
N3=5x14=70
N3=5x15=75

與網頁上例子不同的地方在於，網頁上的複式齒輪的齒數並非整數，所以要選最接近整數的，而作業中的齒輪剛好都為整數
故不用再加以討論

可選用60:12
65:13
70:14
75:15等

依此類推，選出60:12;60:12;60:12三組(或改用其他數字)


2.

我覺得讓我覺得做的作好的一次作業就是第五次作業吧!
因為當時為了這個作業把整個假日，兩天的時間都花在
這個作業上，雖然做出來的手並不是非常的漂亮，但是
其在實際功能的說明對於沒有學過機動學的人也是相當
清楚，我想寫出來程式的那一瞬間應該是那幾天最快樂
的事情吧！

機動學第十二次作業

B94611029 林軍耀
1.

5/31日曾全程來上課。

2.1~2一組標準全齒輪齒輪之徑節為8（亦可使用自設值），
齒數分別為30T與48T，其工作壓力角為20度（可為
14.5或25度，自選）。

設其徑節為8，
齒數為30T和48T
工作壓力角20度
使用程式function [c_ratio,c_length,ad,pc,pb,d2,d3,ag]=contact_ratio(pd,n2,n3, phi)

其中pd=8
n2=30
n3=48
phi=20
可得

接觸比c_ratio=1.7005
接觸長度c_length =0.6275
齒冠ad =0.1250
周節pc =0.3927
基周節pb = 0.3690
齒輪一節圓直徑d2 = 3.7500
齒輪二節圓直徑d3 = 6
齒輪一之接近角、遠退角及作用角ag = 10.4850 9.9211 20.4061
齒輪二之接近角、遠退角及作用角 6.5532 6.2007 12.7538

其中

齒輪一的節徑為30/8=3.75
齒輪二的節徑為48/8=6

齒輪一的基圓直徑為3.75*cos(20)=3.523
齒輪二的基圓直徑為 6*cos(20)=5.638

單位皆為吋


2.3此組齒輪是否會產生干涉現象？試列式證明之。

(N2²+2N2 x N3)sin²>= 4 + 4N3
滿足此式就不會產生干涉

利用講義中公式可以知道

(20^2+2*20*48)(sin(20))^2 >= 4*(1+48)
271.1 > 196

故可知，不會產生干涉現象

2.4

使用draw_gear和move2_gear
利用move2_gear(8,30,48,20,5)
產生動畫


附上draw_gear和move2_gear程式

function move2_gear(Dpitch,nn1,nn2,phi,omega1)
% move2_gear(Dpitch,nn1,nn2,phi,omega1)
% To draw a whole gear
% Inputs:
% Inputs:
% Dpitch:dimetral pitch
% nn1,nn2: no. of teeth for both gears
% phi:pressure angle, degrees
% omega1: angular velocity of gear 1
% Example move2_gear(10,15,20,20,10)
clf;
d2r=pi/180;delt=0.01;
[coord1,r1,rb1]=one_tooth(Dpitch,nn1,phi,360,0,0);
[coord2,r2,rb2]=one_tooth(Dpitch,nn2,phi,360,0,0);
st=180/nn2;if nn1+nn2>2*fix((nn1+nn2)/2),st=0;end
coord2=rotate2D(coord2,180+st,0,0);
xc1=coord1(:,1);yc1=coord1(:,2);
xc2=coord2(:,1);yc2=coord2(:,2);
height=max(r1,r2)*1.2;
ar=min(abs(r1),abs(r2));
coord=bushing(ar/5,0,0); % Get the coordinates of 1st bushing
xb1=coord(:,1)-r1;yb1=coord(:,2);
xb2=coord(:,1)+r2;yb2=coord(:,2);
coord=bushing(-r1,-r1,0);%Get the 1st pitch circle
xp1=coord(:,1);yp1=coord(:,2);
coord=bushing(-r2,r2,0);% Get the 2nd pitch circle
xp2=coord(:,1);yp2=coord(:,2);
plot(xb1,yb1,'r-');hold on;
plot(xb2,yb2,'k-');
plot(xp1,yp1,'r:');
plot(xp2,yp2,'k:');
plot([-r1,r2]',[0,0]','r:');
xx1=min([r1,r2])/2;phir=(90-phi)*d2r;
plot([0,0]',[-xx1*2,xx1*2]','b:');
plot([-xx1 xx1]',[-xx1*tan(phir), xx1*tan(phir)]','b:');

cir1=line('xdata',[],'ydata',[],'erasemode','xor','linewidth',1,'color','r');
cir2=line('xdata',[],'ydata',[],'erasemode','xor','linewidth',1,'color','k');
line1=line('xdata',[],'ydata',[],'erasemode','xor','linewidth',2,'color','r');
line2=line('xdata',[],'ydata',[],'erasemode','xor','linewidth',2,'color','k');
lx1=[0 -r1]';ly1=[0,0]';
lx2=[r2,0]';ly2=[0,0]';
axis([-2.5*r1 2.5*r2 -height height]);
axis equal;
title('Press Ctl-C to stop');
theta1=180;theta2=180;s1=omega1*delt/d2r;
while 1,
z1=rotate2D([xc1,yc1],theta1,-r1,0);
z2=rotate2D([xc2,yc2],theta2,r2,0);
L1=rotate2D([lx1,ly1],theta1,-r1,0);
L2=rotate2D([lx2,ly2],theta2,r2,0);
set(cir1,'xdata',z1(:,1),'ydata',z1(:,2)); % For 1st circle moving
set(cir2,'xdata',z2(:,1),'ydata',z2(:,2)); % For 2nd circle moving
set(line1,'xdata',L1(:,1),'ydata',L1(:,2)); % For 1st line
set(line2,'xdata',L2(:,1),'ydata',L2(:,2)); % For 2nd line
drawnow;
pause(1/s1); %Stop for a while so we can see the graph
theta1=theta1+s1;
theta2=theta2-s1*r1/r2;
if theta1>360, theta1=theta1-360;end; %Reverse the direction at bondary line
if theta2>360,theta2=theta2-360;end;
end




%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function [coords,rp,rb]=one_tooth(Dp,N,phi,range,x0,y0)
% [coords]=draw_gear(Dp,N,phi,range,x0,y0)
% To draw a whole gear
% Inputs:
% Dp: Diametrial pitch
% N: no of teeth in a gear
% phi: pressure angle, degrees
% range: the section range to be drawn
% Example one_tooth(10,10,20,360)
nn=10;
d2r=pi/180;
phir=phi*d2r;
rp=N/Dp/2;
pc=pi/Dp;
ra=rp+1/Dp;
rb=rp*cos(phir);
rd=rp-1.25/Dp;
thpb=pc/rp;% angle respect to one pitch
tp=pc/2;
rr=linspace(rb,ra,nn)';
invphi=tan(phir)-phir;
for i=1:nn
beta=acos(rp/rr(i)*cos(phir));
tt(i)=(tp/rp/2+invphi-tan(beta)+beta);
end
coord1=[rr.*cos(tt') rr.*sin(tt')];
coord3=reverse(coord1);
th1=linspace(thpb/2,thpb/2-tt(nn),nn)';
coord0=[cos(th1) sin(th1)]*rd;
th2=linspace(tt(nn),-tt(nn),nn)';
coord2=[cos(th2) sin(th2)]*ra;
coord4=reverse(coord0);
coord=[coord0;coord1;coord2;coord3;coord4];
theta=thpb/d2r;
coords=[];i=0;
while i coord1=rotate2D(coord,-i,x0,y0);
coords=[coords;coord1];
i=i+theta;
end

function [B]=reverse(A)
nn=length(A); B=A;
for i=1:nn,B(i,:)=A(nn-i+1,:);end
B=[B(:,1) -B(:,2)];

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function [coords] = bushing(rr,x0,y0)
d2r=pi/180;
theta=[360:-10:0]*d2r;
r=abs(rr);
rx=r*cos(theta);ry=r*sin(theta);
if rr<0, rx=rx+x0; ry=ry+y0; coords=[rx' ry']; return; end; rx1=rx/2;rx2=-rx1; ry1=ry/2;ry2=-ry1; r4=r+r/4;r3=r/3; bx=[ 0 0 0 -r -r -r4 -r4 r4 r4 -r r r]; by=[r3 -r3 0 0 -r -r -r4 -r4 -r -r -r 0]; coords(:,1)=[bx rx rx1 rx2]'+x0; coords(:,2)=[by ry ry1 ry2]'+y0;
%%%%%%%%%%%%%%%%%%%
function [coords]=rotate2D(coord,theta,x0,y0)
th=theta*pi/180; c=cos(th);
s=sin(th);fact=[c s;-s c];
coords=coord*fact;
coords(:,1)=coords(:,1)+x0; coords(:,2)=coords(:,2)+y0;
function [coords]=draw_gear(Dp,N,phi,range,x0,y0)
% [coords]=draw_gear(Dp,N,phi,range,x0,y0)
%To draw a whole gear % Inputs: % Dp: Diametrical pitch
% N: no of teeth in a gear % phi: pressure angle, degrees
% range: the section range to be drawn % x0,y0: the location of the gear center
% Example [coords]=draw_gear(10,15,20,360,0,0)
[coord,theta,rp,rb]=tooth(Dp,N,phi); coords=[];
i=0; while i coord1=rotate2D(coord,-i,x0,y0);
coords=[coords;coord1];
i=i+theta;
end
plot(coords(:,1),coords(:,2));hold on;
[coord]=bushing(rp/8,x0,y0);
plot(coord(:,1),coord(:,2),'b-');
[coord]=bushing(-rp,x0,y0);
plot(coord(:,1),coord(:,2),'r:');
[coord]=bushing(-rb,x0,y0);
plot(coord(:,1),coord(:,2),'b:');
axis equal;




%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function [coords,theta,rp,rb]=tooth(Dp,N,phi)
% Example tooth(10,10,20)
nn=10;
d2r=pi/180;
phir=phi*d2r;
rp=N/Dp/2;
pc=pi/Dp;
ra=rp+1/Dp;
rb=rp*cos(phir);
rd=rp-1.25/Dp;
thpb=pc/rp;% angle respect to one pitch
tp=pc/2;
rr=linspace(rb,ra,nn)';
invphi=tan(phir)-phir;
for i=1:nn
beta=acos(rp/rr(i)*cos(phir));
tt(i)=(tp/rp/2+invphi-tan(beta)+beta);
end
coord1=[rr.*cos(tt') rr.*sin(tt')];
coord3=reverse(coord1);
th1=linspace(thpb/2,thpb/2-tt(nn),nn)';
coord0=[cos(th1) sin(th1)]*rd;
th2=linspace(tt(nn),-tt(nn),nn)';
coord2=[cos(th2) sin(th2)]*ra;
coord4=reverse(coord0);
coords=[coord0;coord1;coord2;coord3;coord4];
theta=thpb/d2r;

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function [B]=reverse(A)
nn=length(A); B=A;
for i=1:nn,B(i,:)=A(nn-i+1,:);end
B=[B(:,1) -B(:,2)];

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function [coords] = bushing(rr,x0,y0)
d2r=pi/180;
theta=[360:-10:0]*d2r;
r=abs(rr);
rx=r*cos(theta);ry=r*sin(theta);
if rr<0,
rx=rx+x0;
ry=ry+y0;
coords=[rx' ry'];
return;
end;
rx1=rx/2;rx2=-rx1;
ry1=ry/2;ry2=-ry1;
r4=r+r/4;r3=r/3;
bx=[ 0 0 0 -r -r -r4 -r4 r4 r4 -r r r];
by=[r3 -r3 0 0 -r -r -r4 -r4 -r -r -r 0];
coords(:,1)=[bx rx rx1 rx2]'+x0;
coords(:,2)=[by ry ry1 ry2]'+y0;

%%%%%%%%%%%%%%%%%%%%%%%%

function [coords]=rotate2D(coord,theta,x0,y0)
th=theta*pi/180;
c=cos(th);s=sin(th);fact=[c s;-s c];
coords=coord*fact;
coords(:,1)=coords(:,1)+x0;
coords(:,2)=coords(:,2)+y0;

機動學第十一次作業

b94611029林軍耀

1.本人本週(5/24)有來上課。

2.
使用程式
function plot_dwell(ctheta,s,pattern,range)
%ctheta = cam angle (deg)--can be a matrix
%pattern = denote the type of motion used(a 3 element-row matrix)
% 1:uniform 2:parabolic 3:simple harmonic 4: cycloidal
% 5:polynomial motion
% example [4 3]%range =the degrees the specific motion starts
% Output: y is for displacement, yy is the derivative of the displacement with
% respect to theta, and yyy the second derivative with respect
% to theta.
% Example plot_dwell(0:10:360,2,[4 3],[90 180 240]);figure(1);
clf;
[y,yy,yyy]=dwell(ctheta,range,pattern);
h1=plot(ctheta,y*s,'b-',ctheta,yy*s,'k-',ctheta,yyy*s,'r-');
legend('Displacement','Velocity','Acceleration',3);
xlabel('Elapsed Angle, degrees');
grid;



因為等加速度的位移軌跡
就是拋物線
所以在pattern中輸入2


返程等速
plot_dwell(0:10:360,5,[2 1],[100 200 260]);


返程等加速度
plot_dwell(0:10:360,5,[2 2],[100 200 260]);


返程簡諧
plot_dwell(0:10:360,5,[2 3],[100 200 260]);


返程擺線
plot_dwell(0:10:360,5,[2 4],[100 200 260]);


返程多項式
plot_dwell(0:10:360,5,[2 5],[100 200 260]);



求出各種返程型態之圖形

其中
1:等速運動uniform 2:抛物線parabolic 3:簡諧simple harmonic
4:擺線cycloidal　5:多項式polynomial motion

使用function [y,yy,yyy]=dwell(ctheta,range,pattern)

function [y,yy,yyy]=dwell(ctheta,range,pattern)
%
% This function determines the follower displacement and derivatives
% for a full rotation cam.  The routine is set up for the displacement
% schedule in Examples 6.7 and 6.8% The input values are:
%ctheta   = individual cam angles (deg)--can be a matrix to be processed
%pattern  = denote the types of motion used(a 2 element-row matrix). types
%         are:
%         1:uniform 2:parabolic 3:simple harmonic 4: cycloidal
%         5:polynomial motion%         e.g. [4 3]:cycloidal motion for rise and harmonic motion for
%         return%range    =the degrees the specific motion starts(array of 3 elements), e.g.
%          [90 180 240] starts to rise at 90 deg. and ends at 180 deg.
%           starts to return at 240 deg. and ends at 360 deg.
% Output: y is for displacement, yy is the derivative of the displacement with
% respect to theta, and yyy the second derivative with respect  to theta.
% Example dwell(60,[90 180 240],[4 3]);
% Author:DSFon, BIME, NTU Revise Date:May 18, 2007
d2r=pi/180;
theta=ctheta*d2r;range=range*d2r;
dim=length(ctheta);y=zeros(size(ctheta));
yy=y;yyy=y;
for i=1:dim
    if theta(i)>=range(3)    %for the last motion(downward)
        mode=pattern(2);
        betax=2*pi-range(3);
        switch mode, 
            case 1, [y(i),yy(i),yyy(i)]=uniform(theta(i), range(3),betax,-1);
            case 2, [y(i),yy(i),yyy(i)]=parabolicm(theta(i), range(3),betax,-1);
            case 3, [y(i),yy(i),yyy(i)]=harmonicm(theta(i), range(3),betax,-1);
            case 4, [y(i),yy(i),yyy(i)]=cycloidm(theta(i), range(3),betax,-1);
            case 5, [y(i),yy(i),yyy(i)]=polynorm(theta(i), range(3),betax,-1);
        end;
    elseif theta(i)>=range(2) % dewell on the top
        y(i)=1;elseif theta(i)>=range(1) % for the 1st motion(upward)
        mode=pattern(1);betax=range(2)-range(1);
        switch mode,
            case 1, [y(i), yy(i), yyy(i)]=uniform(theta(i), range(1),betax,+1);
            case 2, [y(i), yy(i), yyy(i)]=parabolicm(theta(i), range(1),betax,+1);
            case 3, [y(i), yy(i), yyy(i)]=harmonicm(theta(i), range(1),betax,+1);
            case 4, [y(i), yy(i), yyy(i)]=cycloidm(theta(i), range(1),betax,+1);
            case 5, [y(i), yy(i), yyy(i)]=polynorm(theta(i), range(1),betax,+1);end
    end
end
%*********************************************
function [t1, t2, t3]=uniform(th, thinit,beta,direct)
% code = 1 for uniform motion
%th=cam angle, radians%beta=motion range, radians
%thinit=starting cam angle, radians%beta=motion range, radians
%direct=motion type; +1 for upward, -1 for downward
theta=th-thinit;t1=theta/beta;
if direct==-1,t1=1-t1;end;
t2=direct*1/beta;t3=0;
%*********************************************
function [t1, t2, t3]=parabolicm(th,thinit,beta,direct)
% code = 2 for parabolic motion%th=cam angle, radians
%beta=motion range, radians%thinit=starting cam angle, radians
%beta=motion range, radians%direct=motion type; +1 for upward, -1 for downward
theta=th-thinit;thmed=thinit+beta/2;thx=theta/beta;
if direct==1,
    if th＜thmed
        t1=2*thx^2;
        t3=4/beta/beta;
        t2=t3*th;
    else
            t1=1-2*(1-thx)^2;
            t2=4/beta*(1-thx);
        t3=-4/beta/beta;
    end
else
    if th        t1=1-2*thx^2;
        t3=-4/beta/beta;
        t2=t3*th;
    else
            t1=2*(1-thx)^2;
            t2=-4/beta*(1-th/beta);
            t3=4/beta/beta;
    end
end
            %*********************************************
            function [t1, t2, t3]=harmonicm(th,thinit,beta,direct)
                % code = 3 for harmonic motion%th=cam angle, radians
                %beta=motion range, radians%thinit=starting cam angle, radians
                %beta=motion range, radians%direct=motion type; +1 for upward, -1 for downward
                theta=th-thinit;t1=0.5*(1-cos(pi*theta/beta));
                if direct==-1,  t1=1-t1;end;
                t2=direct*(0.5*pi/beta)*sin(pi*theta/beta);
                t3=direct*0.5*(pi/beta)^2*cos(pi*theta/beta);
                %*********************************************
                function [t1, t2, t3]=cycloidm(th,thinit,beta,direct)
                    %% code = 4 for cycloidal motion%th=cam angle, radians
                    %thinit=starting cam angle, radians%beta=motion range, radians
                    %direct=motion type; +1 for upward, -1 for downward
                    theta=th-thinit;
                    t1=theta/beta-(0.5/pi)*sin(2*pi*theta/beta);
                    if direct==-1,t1=1-t1;end;
                    t2=direct*(1-cos(2*pi*theta/beta));
                    t3=direct*2*pi/beta/beta*sin(2*pi*theta/beta);
                    %*********************************************
                    function [t1, t2, t3]=polynorm(th,thinit,beta,direct)
                        %% code = 5 for polynormial motion%th=cam angle, radians
                        %thinit=starting cam angle, radians
                        %beta=motion range, radians%direct=motion type; +1 for upward, -1 for downward
                        theta=th-thinit;
                        thx=theta/beta;t1=thx*thx*thx*(10+thx*(-15+thx*6));
                        if direct==-1, t1=1-t1;end;
                        t2=direct*(30/beta)*thx*thx*(1+thx*(-2+thx));
                        t3=direct*(60/beta/beta)*thx*(1+thx*(-3 +2*thx));

3.

利用[x y]=pincam([0:10:360],15,s,0,10,[100 200 260],[2 1],-1)
以及[s]=drawcam(15,y,1)
作出圖形
與第四小題的動畫一同呈現

4.

利用上面的程式
[s]=drawcam(15,y,1)

再用
for i=1:1:360
x2=s(:,1);
y2=s(:,2);
x3=x2*cosd(i)-y2*sind(i);
y3=x2*sind(i)+y2*cosd(i);
axis([-30 30 -30 30]);
plot(x3,y3);
pause(0.01);
axis([-30 30 -30 30]);
axis equal;
grid on;
end;

作出動畫


附件1
function [x,y]=pincam(cth,r0,s,e,L,range,pattern,cw)
%Find the pin type cam with an offsect e%Inputs:
% cth:angle of cam, degrees% r0:radius of base circle
% e:offset% s:stroke
% L:length of pin
% cw:rotation direction of cam(-counterclockwise,+clockwise
%pattern = denote the type of motion used(a 3 element-row matrix)
% 1:uniform 2:parabolic 3:simple harmonic 4: cycloidal
% 5:polynomial motion% example [4 3]
%range =the degrees the specific motion starts, eg.[90 180 240]
% Example: [x y]=pincam([10 60],5,2,1,10,[90 180 240],[4 3],-1)
figure(1);
clf;
th=cth*pi/180;
s0=sqrt(r0*r0-e*e);
for i=1:length(cth)
t=th(i)*cw;
A=[cos(t) -sin(t);sin(t) cos(t)];
[ym,yy,yyy]=dwell(cth(i),range,pattern);
x0=s0+ym*s;
Sx=[0 x0 x0+L;e e e];
X=A\Sx;
x(i)=X(1,2)
y(i)=X(2,2)
line(X(1,1:2),X(2,1:2));
line(X(1,2:3),X(2,2:3),'linewidth',3,'color','red')
end
hold on;
plot([0 x],[0 y],'ro',x,y,'k-')

附件2
function [rp,rb]=drawcam(r0,y,direct)
% To draw a cam profile with base radius of r0 and rise in y
% The program may work with P8_11% Input: direct: +1 for clockwise; -1 for counterclockwise
%
% Example: [s,T]=drawcam(35,y)
d2r=pi/180;nn=length(y);
LM=max(y+r0);figure(2);
line([-LM LM]',[0 0]');
line([0 0]',[-LM LM]');
if direct==1,
theta=linspace(0,360,nn)'*d2r;
else
theta=linspace(360,0,nn)'*d2r;
end
theta=theta+pi/2;
rb=[r0*cos(theta) r0*sin(theta)];
rp=[(r0+y').*cos(theta) (r0+y').*sin(theta)];
line(rb(:,1),rb(:,2),'color','r');
line(rp(:,1),rp(:,2));
axis equal;
grid on;

機動學第十次作業

本週(5/17)有來上課。

10.2



當一桿以桿上一點M做等角速度迴轉運動時，端點之速度為端點至M點之距離乘以角速度w
若M點為桿外一點，V=仍然為端點至M點之距離*ω
a=端點至M之距離*ω*ω

假若M以V等速水平移動
θ為桿與水平面之夾角
圖請見網頁
則P點此時之速度為

V+r*ω*sinθ - r*ω*cosθ

加速度維持不變



假若M又具有加速度a)時
則P點之速度假設與上述相同，
則加速度為

V+r*ω^2*cosθ + r*ω^2*sinθ


10.3

function sldwork(R,L,e)
th1=slider_limit1(R,L,e);
th2=90;
angle=linspace(th1,th2,100);
d=slider_solved(angle,R,L,e,1);
x=R*cosd(angle);
y=R*sind(angle);

for i=1:100
hold on
plot([0,x(i),d(i)],[0,y(i),e],'linewidth',4);
plot([d(i)-2,d(i)+2,d(i)+2,d(i)-2,d(i)-2],[e-2,e-2,e+2,e+2,e-2]);
plot([x(i),0],[y(i),e-d(i)*(y(i)-e)/(x(i)-d(i))],'b.','linewidth',8)
plot([d(i),d(i)],[0,y(i)*d(i)/x(i)],'b.')
axis equal;
axis ([-10 50 -10 50]);
pause(0.03);
end

程式如上

動畫中出現的兩個點為瞬心
還有三個瞬心分別出現在旋轉結以及滑塊與桿之連結點上
第六個瞬心在無限遠的地方
為滑塊的瞬心

機動學第九次作業

本人本週(5/3)有來上課

b94611029 林軍耀


由於偏置量會影響滑桿之活動範圍
所以都曲柄如果要轉360度時
其連桿與曲桿及滑塊三者必須構成
某種三角形的關係

連桿長-曲桿長 大於 偏置量

由課本4-35頁的slider_limit程式可知道，
由於偏置量e大於連桿與曲桿之差，
所以會產生極限點
在theta = 6.9199度時，有右極限角
　
令R=29+10=39
　L=39+5=44
　e=10
參考下圖


以及利用move_sldpaths([10,R,L,10],0,0,3,10,10,0,1,0,4,100)以及
slider_draw(39,44,10)做出動畫

機動學第八次作業

本人（4／26）曾來上課。
8-1
利用f4bar程式

function [data,form] = f4bar(r,theta1,theta2,td2,tdd2,mode,linkdrive)
%
%function [data,form] = f4bar(r,theta1,theta2,td2,tdd2,mode,linkdrive)
% This function analyzes a four-bar linkage when the crank is the
% driving link. The input data are:
% theta1,theta2 are angles in degrees
% r=[r1 r2 r3 r4]= lengths of links(1=frame)
%td2 = crank or coupler angular velocity (rad/sec)
%tdd2 = crank or coupler angular acceleration (rad/sec^2)
%mode = +1 or -1. Identifies assembly mode
%linkdrive = 0 for crank as driver; 1 for coupler as driver
%data (1:4,1) = link positions for 4 links
%data (1:4,2) = link angles in degrees
%data (1:4,3) = link angular velocities
%data (1:4,4) = link angular accelerations
%data (1,5) = velocity of point Q
%data (2,5) = velocity of point P
%data (3,5) = acceleration of point Q
%data (4,5) = acceleration of point P
%data (1,6) = position of Q
%data (2,6) = position of P
%form = assembly status. form = 0, mechanism fails to
% assemble.
% program revised from Waldron's Textbook
% Revised:DSFON, BIME, NTU. Date: Feb. 7, 2007
if nargin<7,linkdrive=0;end
if nargin<6,mode=1;end
data=zeros(4,6);
% if coupler is the driver, interchange the vetor 3 & 2
if linkdrive==1,r=[r(1) r(3) r(2) r(4)];end
rr=r.*r;d2g=pi/180;
[theta,td,tdd]=deal(zeros(4,1));
theta(1:2)=[theta1 theta2]*d2g;
t1=theta(1);tx=theta(2);
s1=sin(t1);c1=cos(t1);
sx=sin(tx);cx=cos(tx);
% position calculations
A=2*r(1)*r(4)*c1-2*r(2)*r(4)*cx;
C=rr(1)+rr(2)+rr(4)-rr(3)-2*r(1)*r(2)*(c1*cx+s1*sx);
B=2*r(1)*r(4)*s1-2*r(2)*r(4)*sx;
pos=B*B-C*C+A*A;
if pos>=0,
form=1;
% Check for the denominator equal to zero
if abs(C-A)>=1e-5
t4=2*atan((-B+mode*sqrt(pos))/(C-A));
s4=sin(t4);c4=cos(t4);
t3=atan2((r(1)*s1+r(4)*s4-r(2)*sx),(r(1)*c1+r(4)*c4-r(2)*cx));
s3=sin(t3);c3=cos(t3);
else
% If the denominator is zero, compute theta(3) first
A=-2*r(1)*r(3)*c1+2*r(2)*r(3)*cx;
B=-2*r(1)*r(3)*s1+2*r(2)*r(3)*sx;
C=rr(1)+rr(2)+rr(3)-rr(4)-2*r(1)*r(2)*(c1*cx+s1*sx);
pos=B*B-C*C+A*A;
if pos>=0,
t3=2*atan((-B-mode*sqrt(pos))/(C-A));
s3=sin(t3); c3=cos(t3);
t4=atan2((-r(1)*s1+r(3)*s3+r(2)*sx),...
(-r(1)*c1+r(3)*c3+r(2)*cx));
s4=sin(t4);c4=cos(t4);
end
end
theta(3)=t3;theta(4)=t4;
%velocity calculation
td(2)=td2;
AM=[-r(3)*s3, r(4)*s4; -r(3)*c3, r(4)*c4];
BM=[r(2)*td(2)*sx;r(2)*td(2)*cx];
CM=AM\BM;
td(3)=CM(1);td(4)=CM(2);

%acceleration calculation

tdd(2)=tdd2;
BM=[r(2)*tdd(2)*sx+r(2)*td(2)*td(2)*cx+r(3)*td(3)*td(3)*c3-...
r(4)*td(4)*td(4)*c4;r(2)*tdd(2)*cx-r(2)*td(2)*td(2)*sx-...
r(3)*td(3)*td(3)*s3+r(4)*td(4)*td(4)*s4];
CM=AM\BM;
tdd(3)=CM(1);tdd(4)=CM(2);
%store results in array data
% coordinates of P and Q
if linkdrive==1,
c2=c3;c3=cx;s2=s3;s3=sx;
r(2:3)=[r(3) r(2)];theta(2:3)=[theta(3) theta(2)];
td(2:3)=[td(3) td(2)];tdd(2:3)=[tdd(3) tdd(2)];
else
c2=cx;s2=sx;
end
for j=1:4,
data(j,1:4)=[r(j)*exp(i*theta(j)) theta(j)/d2g td(j) tdd(j)] ;
end % position vectors
data(1,5)=r(2)*td(2)*exp(i*theta(2));%velocity for point Q
data(2,5)=r(4)*td(4)*exp(i*theta(4));%velocity for point P
data(3,5)=r(2)*(i*tdd(2)-td(2)*td(2))*exp(i*theta(2));%acc of Q
data(4,5)=r(4)*(i*tdd(4)-td(4)*td(4))*exp(i*theta(4));%acc of P
data(1,6)=data(2,1);%position of Q, again
data(2,6)=data(1,1)+data(4,1);% position of P

%find the accelerations
else
form=0;
if linkdrive==1,
r=[r(1) r(3) r(2) r(4)];
for j=1:4, data(j,1)=r(j).*exp(i*theta(j));end % positions
end
end


輸入[val,form]=f4bar([4 3 3 5],0,45,10,0,-1,0)



位置　　　　　　　角速度　　　　　角加速度

I (0.0,0.0)　　　　　　　0　　　　　　　　0
II (2.1,2.1)　　　　10.0000　　　　　　　　0
III (3.2,4.9)　　　　16.2681　　　　 491.4428
IV (4.0,0.0)　　　　4.9677　　　　　383.6120

8-2
利用drawlinks程式為之
此程式有利用到f4bar的值
所以在程式中只要設定好桿長以及角度
再加上模式與驅動桿的設定，就可以完成

function [values]=drawlinks(r,th1,th2,mode,linkdrive)
if nargin<5,linkdrive=0;end
if nargin<4,mode=1;end
[values b]=f4bar(r,th1,th2,0,0,mode,linkdrive);
rr=values(:,1);
rr(3,1)=rr(1,1)+rr(4,1);
rx=real(rr(:,1));rx(4)=0;
ry=imag(rr(:,1));ry(4)=0;
if b==1
plot([0 rx(1)],[0 ry(1)],'k-','LineWidth',4);
hold on;
if linkdrive==0
plot([0 rx(2)],[0 ry(2)],'b-','LineWidth',1.5);
plot([rx(2) rx(3)],[ry(2) ry(3)],'r-','LineWidth',2);
else
plot([0 rx(2)],[0 ry(2)],'r-','LineWidth',2);
plot([rx(2) rx(3)],[ry(2) ry(3)],'b-','LineWidth',1.5);
end
plot([rx(1) rx(3)],[ry(1) ry(3)],'g-','LineWidth',1.5);
plot(rx,ry,'bo');
text(0,0,' O');text(rx(1),ry(1),' R');
text(rx(2),ry(2),' P');text(rx(3),ry(3),' Q');
else
fprintf('Combination of links fail at degrees %6.1f\n',th2);
end
axis equal
grid on




8-3
利用drawlimits程式
function drawlimits(r,th1,sigma,driver)
[Qstart, Qstop]=fb_angle_limits(r,th1,driver)
clf;
[values b]=f4bar(r,th1,Qstart,0,0,sigma,driver);
rr=values(:,1);
rr(3)=rr(1)+rr(4);
rx=real(rr);rx(4)=0;
ry=imag(rr);ry(4)=0;
if b==1
plot([0 rx(1)],[0 ry(1)],'k-','LineWidth',4);
hold on;
if driver==0
plot([0 rx(2)],[0 ry(2)],'b-','LineWidth',1.5);
plot([rx(2) rx(3)],[ry(2) ry(3)],'r-','LineWidth',2);
else
plot([0 rx(2)],[0 ry(2)],'r-','LineWidth',2);
plot([rx(2) rx(3)],[ry(2) ry(3)],'b-','LineWidth',1.5);
end
plot([rx(1) rx(3)],[ry(1) ry(3)],'-g');
plot(rx,ry,'bo');
text(0,0,' O');text(rx(1),ry(1),' R');
text(rx(2),ry(2),' P');text(rx(3),ry(3),' Q');
text(rx(2)/2,ry(2)/2,['s1=' num2str(Qstart,'%6.1f')]);
else
fprintf('Combination of links fails at degrees %6.1f\n',Qstart);
end
[values b]=f4bar(r,th1,Qstop,0,0,sigma,driver);
rr=values(:,1);
rr(3)=rr(1)+rr(4);
rx=real(rr);rx(4)=0;
ry=imag(rr);ry(4)=0;
if b==1
if driver==0
plot([0 rx(2)],[0 ry(2)],'b-','LineWidth',1);
plot([rx(2) rx(3)],[ry(2) ry(3)],'r-','LineWidth',1.5);
else
plot([0 rx(2)],[0 ry(2)],'r-','LineWidth',1.5);
plot([rx(2) rx(3)],[ry(2) ry(3)],'b-','LineWidth',1);
end
plot([rx(1) rx(3)],[ry(1) ry(3)],'g-');
plot(rx,ry,'bo');
text(rx(2),ry(2),' p');text(rx(3),ry(3),' q');
text(rx(2)/2,ry(2)/2,[' s2=' num2str(Qstop,'%6.1f')]);
else
fprintf('Combination of links fail at degrees %6.1f\n',Qstop);
end
axis equal
grid on

輸入drawlimits([4 3 3 5],0,1,0)

Qstart = 28.9550


Qstop = 331.0450


如圖所示



8-4
for i=0:20:360
drawlinks([4 3 3 5],0,i,-1,0);
pause(1)
end
Combination of links fail at degrees 0.0
Combination of links fail at degrees 20.0
Combination of links fail at degrees 340.0
Combination of links fail at degrees 360.0

我們會發現在某些角度譬如說20，340，360等沒有值
因為在旋轉四連桿的過程中，必須符合葛拉索定理
其實也就是連桿之間要構成封閉三角形的樣子
如果沒有構成封閉三角刑，那麼四連桿當然無法運作

r1 + r2 ＜ r3 +r4
|r1 -r2|＞|r3 - r4|

如下圖




8-5
利用網頁上講義的程式move_4paths
此程式必須用到相當多其他的子程式
例如f4bar.m funcion, f4limits.m, fb_angle_limits.m
等等，可以說是一個整合的結果
function move_4paths(r,r6,th6,nlink,th1,td2,tdd2,sigma,driver,ntimes,npts)
%
%function move_4paths(r,r6,th6,nlink,th1,td2,tdd2,sigma,driver,ntimes,npts)
%
%draw the positions of four-bar links
%call f4bar.m funcion, f4limits.m, fb_angle_limits.m, body.m
%
%Inputs:
% r: row vector for four links
% th1: frame angle
% th2: crank angle or couple angle
% td2,tdd2:angular velocity and acceleration of the driving link.
% sigma: assembly mode
% driver: 0 for crank, 1 for coupler
% ntimes: no. of cycles
% npts: number of points divided
% r6,rh6,nlink:additional length and angle for nlink link.
%example:
% move_4paths([4 2 3 4],2,-30,3,0,10,0,1,0,4,100)
%
%clf;
if nargin<10, ntimes=3;npts=100;end;
figure(1);
[Qstart, Qstop]=fb_angle_limits(r,th1,driver);
npoint=abs(npts);
th2=linspace(Qstart,Qstop,npoint);
val=zeros(6,npoint);
for i=1:npoint,
[vr b]=f4bar(r,th1,th2(i),td2,tdd2,sigma,driver);
[para]=body(r6,th6,vr,nlink);
val(1:3,i)=[vr(1,1);vr(2,1);vr(1,1)+vr(4,1)];
val(4:6,i)=[para(1);para(3);para(2)];
end
x=real(val);y=imag(val);
h=f4limits(r,th1,sigma,driver);
line(x(5,:)',y(5,:)','color','r','linestyle',':');
line(x(4,:)',y(4,:)','color','b','linestyle','-.');
line(x(6,:)',y(6,:)','color','k','linestyle',':');
range=1.2*([min(min(x)) max(max(x)) min(min(y)) max(max(y))]);
axis(range);axis equal;grid off;
for i=2:4,set(h(i),'erasemode','xor');end
h0=patch('xdata',[],'ydata',[],'erasemode','xor','facecolor','r',...
'marker','o');
i=0;s=-1;axis off;
for m=1:ntimes
s=-s;
if abs(Qstop-Qstart-360)<1,i=0;s=1;end;
while 1,
i=i+s;
if i>npoint|i==0,break;end;
set(h(2),'xdata',[0 x(2,i)], 'ydata',[0 y(2,i)]);%crank
set(h(3),'xdata',[x(2,i) x(3,i)], 'ydata',[y(2,i) y(3,i)]);%coupler
set(h(4),'xdata',[x(1,i) x(3,i)], 'ydata',[y(1,i) y(3,i)]);%Rocker
set(h0,'xdata',[x(4:6,i)], 'ydata',[y(4:6,i)]);
drawnow; %flush the draw buffer
pause(0.1);
end
end % for m loop
function h=f4limits(r,th1,sigma,driver)
[Qstart, Qstop]=fb_angle_limits(r,th1,driver)
[values b]=f4bar(r,th1,Qstart,0,0,sigma,driver);
if b==1,
h=draw4link(values,driver);
else
fprintf('Combination of links fails at degrees %6.1f\n',Qstart);
end
[values b]=f4bar(r,th1,Qstop,0,0,sigma,driver);
if b==1,
h=draw4link(values,driver);
else
fprintf('Combination of links fails at degrees %6.1f\n',Qstart);
end
axis equal
grid on


輸入
move_4paths([4 3 3 5],0,45,3,0,0,0,1,0,4,100)

可得下面的影片