2/6/07

7.3 建立滑塊程式

曲桿滑塊之分析仍可當做一般四連桿之型式進行分析,故其所需之程式僅需針對一般四連桿加以修改即可。曲桿滑塊雖仍有四連桿之結構,但是第一桿則可代表滑塊滑動的方向,故其角度會影響滑塊之滑動路徑,而其長短亦代表滑塊在該路徑移動的距離。故此時之第一桿之方位雖仍為固定,但長度將會變動,此外第四桿其長度雖不變動,但角度需與第一桿相互垂直,其長度因此可以表示滑塊與水平線之偏置距離。在這種情況下,第四桿之長度可為零,代表滑塊完全在第一桿之方向上滑動,這也可能是滑塊之最常用之位置。第四章中之曲桿滑塊之四連桿分析可供參考。本節則針對其動作為MATLAB動作程式加以說明。目前所設計之程式有sldlink.m、drawsldlinks.m、 sld_angle_limits.m、drawsldlimits.m等四個程式,茲分別說明如下:

一、sldlink函數



sldlink 函數之呼叫格式如下:

function [values,form] = sldlink(r,theta1,theta2,td2,tdd2,sigma,driver)

輸入變數:


  • r(1:4) = 各桿之長度,r(1)為固定桿,其餘分別為曲桿、結合桿及被動桿,但r(4)值可有可無,因為會被更改。
  • theta1 = 第一桿之水平角,或為四連桿之架構角,以角度表示。
  • theta2 = 驅動桿之水平夾角,以角度表示。一般為曲桿角,但若為結合桿驅動,則為結合桿之水平夾角。
  • td2 = 驅動桿(第二桿或第三桿)之角速度(rad/sec)。
  • tdd2 = 驅動桿(第二桿或第三桿)之角加速度(rad/sec^2)。
  • sigma = +1 or -1. 組合模式,負值表示閉合型,正值為分支型,但有時需視實際情況而定。
  • driver = 0 (驅動桿為第二桿); 1 (驅動桿為第三桿) ;2(滑塊驅動)

輸出變數:



  • form = 組合狀態, 0 :表示無法組合; 1:可以正確組合
  • values = 輸出矩陣,其大小為 4 X 8,各行之資料分配如下:








桿序
12(deg)3(rad/s)4(rad/s2)5678
桿1位置r1θ1ω1α1r1'Vq|Vq|∠Vq
桿2位置r2θ2ω2α2r1"Vp|Vp|∠Vp
桿3位置r3θ3ω3α3Q2AQ|AQ|∠AQ
桿4位置r4θ4ω4α4P3AP|AP|∠AP


其中第一行之連桿位置向量,屬於單桿的位置向量。第二行為各桿之水平夾角,第三及第四行為各桿之角度速度及角加速度。第五行為第一桿之軸向速度與加速度、P及Q點之位置向量。第六至八行則為P點與Q點之速度與加速度量,第五行為向量,第六行為絕對量,第七行為夾角。
內容上,值得一提的是第一行、三行、四行及五行之向量表示法屬於複數之型式。故若要得到其絕對值僅需在MATLAB指令檔中,以abs()這一個函數指令即可求得,而以函數angle()則可求得其夾角,雖然第二行與第七行之輸出亦有相對應之夾角。

程式內容




function [data,form] = sldlink(r,theta1,theta2,td2,tdd2,sigma,driver)
%
% function [data,form] = sldlink(r,theta1,theta2,td2,tdd2,sigma,driver)
% rivesed from Waldron's
% This function analyzes a slider crank mechanism
% driver =0,1,2 for link 2,3 & slide as crank,with td2 & tdd2 being their
% angular velocity & acceleration.
%Input:
% theta1,theta2= Angles of frame & crank in degrees
% r = [r1 r2 r3 r4], r4:slider offset.
% td2 = crank or coupler angular velocity (rad/sec)
% tdd2 = crank or coupler angular acceleration (rad/sec^2)
% sigma = +1 or -1. Identifies assembly mode
% data = data stored accordingly,
% | 1 | 2 | 3 | 4 | 5 | 6 | 7 | 8 |
% I | r1 | t1 | td1 | tdd1 | r1' | Vq | |Vq|| ThVq |
% II | r2 | t2 | td2 | tdd2 | r1" | Vp | |Vp|| ThVp |
% III| r3 | t3 | td3 | tdd3 | Q2 | Aq | |Aq|| ThAq |
% IV | r4 | t4 | td4 | tdd4 | P3 | Ap | |Ap|| ThAp |
%form = assembly flag. form = 0 for failure of assembly.
%Example: [val,form] = sldlink([1 2 3 4],0,60,10,0,1,0)
f=@(num,ndg) round(num*10^ndg)/10^ndg;
data=zeros(4,8);
if driver==1, r(2:3)=[r(3) r(2)];end
rr=r.*r;d2g=pi/180;
theta=[theta1 theta2 0 theta1+90]*d2g;
t1=theta(1);t4=theta(4);
[td,tdd,rd,rdd]=deal(zeros(4,1));
s1=sin(t1);c1=cos(t1);
s4=sin(t4);c4=cos(t4);
switch driver
case 2 % for the sliding block as driver
rd(1)=td2;
rdd(1)=tdd2;
A=-2*r(1)*r(2)*c1-2*r(2)*r(4)*c4;
B=-2*r(1)*r(2)*s1-2*r(2)*r(4)*s4;
C=rr(1)+rr(2)+rr(4)-rr(3)+2*r(1)*r(4)*(c1*c4+s1*s4);
arg=B*B-C*C+A*A;
arg=f(arg,5);
if (arg>=0)
form=1; %assembly flag
t2=2*atan((-B+sigma*sqrt(arg))/(C-A));
s2=sin(t2); c2=cos(t2);
t3=atan2((r(1)*s1+r(4)*s4-r(2)*s2),(r(1)*c1+r(4)*c4-r(2)*c2));
theta(2)=t2; theta(3)=t3;
s3=sin(t3); c3=cos(t3);
%velocity calculation
AM=[-r(2)*s2, -r(3)*s3; r(2)*c2, r(3)*c3];
BM=[rd(1)*c1;rd(1)*s1];
CM=AM\BM;
td(2)=CM(1); td(3)=CM(2);
%acceleration calculation
BM=[r(2)*td(2)*td(2)*c2+r(3)*td(3)*td(3)*c3+rdd(1)*c1;...
r(2)*td(2)*td(2)*s2+r(3)*td(3)*td(3)*s3+rdd(1)*s1];
CM=AM\BM;
tdd(2)=CM(1); tdd(3)=CM(2);
else
form=0; return;
end
case {0,1} % for link 2 or 3 as a crank
td(2)=td2; tdd(2)=tdd2; tx=theta(2);
sx=sin(tx); cx=cos(tx);
% position calculations
A=2*r(4)*(c1*c4+s1*s4)-2*r(2)*(c1*cx+s1*sx);
B=rr(2)+rr(4)-rr(3)-2*r(2)*r(4)*(cx*c4+sx*s4);
arg=A*A-4*B;
arg=f(arg,5);
if (arg>=0)
form=1; %assembly flag
r(1)=(-A+sigma*sqrt(arg))/2;
t3=atan2((r(1)*s1+r(4)*s4-r(2)*sx),(r(1)*c1+r(4)*c4-r(2)*cx));
theta(3)=t3;
s3=sin(t3); c3=cos(t3);
%velocity calculation
AM=[c1, r(3)*s3; s1, -r(3)*c3];
BM=[-r(2)*td(2)*sx;r(2)*td(2)*cx];
CM=AM\BM;
rd(1)=CM(1); td(3)=CM(2);
%acceleration calculation
BM=[-r(2)*tdd(2)*sx-r(2)*td(2)*td(2)*cx-r(3)*td(3)*td(3)*c3;...
r(2)*tdd(2)*cx-r(2)*td(2)*td(2)*sx-r(3)*td(3)*td(3)*s3];
CM=AM\BM;
rdd(1)=CM(1); tdd(3)=CM(2);
else
form=0;
if driver==1,
r(2:3)=[r(3) r(2)];
for j=1:4, data(j,1)=r(j).*exp(i*theta(j));end % position vectors
end
return;
end
if driver==1,
r(2:3)=[r(3) r(2)];
c2=c3;c3=cx;s2=s3;s3=sx;
td(2:3)=[td(3) td(2)];
tdd(2:3)=[tdd(3) tdd(2)];
theta(2:3)=[theta(3) theta(2)];
elseif driver==0
c2=cx;s2=sx;
end
end % end of switch
for j=1:4,
data(j,1)=r(j).*exp(i*theta(j));
data(j,2:4)=[theta(j)/d2g,td(j),tdd(j)];
end
data(:,5)=[rd(1);rdd(1);data(2,1);data(2,1)+data(3,1)];
data(1,6)=i*r(2)*td(2)*exp(i*theta(2));
data(2,6)=data(1,6)+i*r(3)*td(3)*exp(i*theta(3));
data(3,6)=r(2)*(i*tdd(2)-td(2)*td(2))*exp(i*theta(2));
data(4,6)=data(3,6)+r(3)*(i*tdd(3)-td(3)*td(3))*exp(i*theta(3));
%[Vq;Vp;Aq;Ap]
data(1:4,7)=abs(data(1:4,6));data(1:4,8)=angle(data(1:4,6))/d2g;





繪出滑塊連桿位置函數drawsldlinks程式


程式必須呼叫drawsldlinks函數繪出對應之連桿及滑塊位置。此程式未來可以作為其他多項之用途。

function [values]=drawsldlinks(r,th1,th2,sigma,driver)
%function drawlinks(r,th1,th2,sigma,driver)
%draw the positions of four-bar links
%call f4bar funcion
%r: row vector for four links
%th1: frame angle
%th2: crank angle or couple angle
%sigma: assembly mode
%driver: 0 for crank, 1 for coupler
%
% Example:
% drawsldlinks([1 2 3 4],0,60,1,0)
%clf;
[values b]=sldlink(r,th1,th2,10,0,sigma,driver);
rr=values(:,1);
rr(3)=rr(3)+rr(2);
rx=real(rr);rx(4)=0;
ry=imag(rr);ry(4)=0;
if b==1
plot([0 rx(1)],[0 0],'k-','LineWidth',4);%ground line
hold on;
plot([0 rx(1)],[0 ry(1)],'g-','LineWidth',1.5);
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)],'k-');
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');
length=max(abs(values(2:3,1)));
len=.20*length;ww=.15*length;
[coords] = sldbox(len,ww,rx(3),ry(3),th1);
plot(coords(:,1),coords(:,2),'r-','LineWidth',2);
[coords] = sldbox(len*3,0,rx(3),ry(3)-ww/2,th1);
plot(coords(:,1),coords(:,2),'r:','LineWidth',1.5);
else
fprintf('Combination of links fails at degrees %6.1f\n',th2);
end
axis equal
grid on


執行範例



此時之輸入角速度td及角加速度tdd項所代表之值為第一桿之軸向線速度及線加速度,在輸入時應特別注意。

例一:為第二桿為驅動桿



>> [val,form] = sldlink([1 2 3 4],0,60,10,0,1,0)
val =
Columns 1 through 3
2.9638 0 0
1 + 1.7321i 60 10
1.9638 + 2.2679i 49.111 -5.0922
2.4493e-016 + 4i 90 0
Columns 4 through 6
0 -5.7716 -17.321 + 10i
0 -418.87 -5.7716 +1.7764e-015i
118.15 1 + 1.7321i -100 - 173.21i
0 2.9638 + 4i -418.87 -2.8422e-014i
Columns 7 through 8
20 150
5.7716 180
200 -120
418.87 -180
form =
1


本例中,有框線者表示其為輸入值。但第一行則已經轉換為複數型式,未來複數型式要轉為x-y座表時,只要使用函數real()及imag()兩指令,即可進行轉換。值得注意的是,第一桿之數值已經變動,將val矩陣取絕對值如下:

>> abs(val)
ans =
Columns 1 through 6
2.9638 0 0 0 5.7716 20
2 60 10 0 418.87 5.7716
3 49.111 5.0922 118.15 2 200
4 90 0 0 4.9783 418.87
Columns 7 through 8
20 150
5.7716 180
200 120
418.87 180


使用繪圖程式drawsldlinks可以得到繪圖的結果與資料:

>> drawsldlinks([1 2 3 4],0,60,1,0)
ans =
Columns 1 through 3
2.9638 0 0
1 + 1.7321i 60 10
1.9638 + 2.2679i 49.111 -5.0922
2.4493e-016 + 4i 90 0
Columns 4 through 6
0 -5.7716 -17.321 + 10i
0 -418.87 -5.7716 +1.7764e-015i
118.15 1 + 1.7321i -100 - 173.21i
0 2.9638 + 4i -418.87 -2.8422e-014i
Columns 7 through 8
20 150
5.7716 180
200 -120
418.87 -180




例二:為第三桿(coupler)為驅動桿




>> drawsldlinks([1 2 3 4],0,60,1,1)
ans =
Columns 1 through 3
2.9264 0 0
1.4264 + 1.4019i 44.504 -10.516
1.5 + 2.5981i 60 10
2.4493e-016 + 4i 90 0
Columns 4 through 6
0 -11.238 14.743 - 15i
290.83 -715.46 -11.238
0 1.4264 + 1.4019i -565.46 + 259.81i
0 2.9264 + 4i -715.46 -5.6843e-014i
Columns 7 through 8
21.032 -45.496
11.238 180
622.29 155.32
715.46 -180

此時之第一桿已經與第一例略有不同。可能是計算上之誤差。其結果如下圖。



例三:以滑塊為驅動桿(slider)為驅動桿(第一桿)



>> drawsldlinks([1 2 3 4],0,60,1,2)
ans =
Columns 1 through 3
1 0 0
-0.97808 + 1.7445i 119.28 -3.4968
1.9781 + 2.2555i 48.749 -1.729
2.4493e-016 + 4i 90 0
Columns 4 through 6
0 10 6.1002 + 3.4202i
-9.0794 0 10
9.7031 -0.97808 + 1.7445i 27.799 - 12.451i
0 1 + 4i 0 +1.7764e-015i
Columns 7 through 8
6.9936 29.278
10 0
30.46 -24.127
1.7764e-015 90