用matlab分析曲柄摇杆机构.doc

完整word)用matlab分析曲柄摇杆机构
机械原理第一次作业(matlab7。0)：

求：r1旋转360°时，θ2，θ3，ω2，ω3,α2,α3和C点的加速度。
设r1=400，r2=1000,r3=700,r4=1200
1、角位移的M函数：
function y=jweiyi(x）
% Input parameters
％ x（1）=theta—1
% x（2）=theta—2 guess value
% x(3）=theta-3 guess value
% x(4)=r1
% x(5)=r2
％ x（6)=r3
％ x(7)=r4
% Output parameters
％ y(1）=theta-2
% y(2）=theta—3
theta2=x(2）;
theta3=x（3)；
％
epsilon=1.0E—6；
%
f=［x（4）*cos（x（1))+x（5)＊cos（theta2）—x(7）—x(6)*cos（theta3）；
x(4）＊sin(x(1)）+x(5)＊sin（theta2)—x（6)＊sin（theta3)];
%
while norm(f)〉epsilon
J=［—x（5)＊sin（theta2） x(6）*sin(theta3）;x(5）＊cos(theta2) —x（6）cos（theta3）］;
dth=inv（J）*（-1.0＊f）；
theta2=theta2+dth(1)；
theta3=theta3+dth(2)；
f=[x(4)＊cos（x(1）)+x(5）＊cos（theta2)-x（7）—x（6）＊cos（theta3);
x（4)＊sin（x(1）)+x（5）*sin（theta2）—x（6）＊sin（theta3）］;
norm（f)；
end;
y(1）=theta2；
y(2)=theta3；
r1旋转360°时，θ2，θ3的M文件程序：
r（1）=400;
r(2)=1000；
r(3)=700;
r(4)=1200；
dr=pi/180;
th(1）=0;
th（2)=44。0486*dr;
th（3)=96.6654*dr;
y=jweiyi（[th(1)，th(2）,th(3），r(1），r（2),r（3),r（4）])
dth=1*dr；
for i=1:360
y=jweiyi([th（1）,th(2),th(3），r(1)，r(2),r（3),r（4）])；
th23(i,:)=[th（1）/dr，th（2)/dr,th（3）/dr];
th（1）=th(1）+dth;
th(2)=y(1);
th(3)=y（2）；
end
figure（1）
plot(th23(：,1)，th23(：,2)，th23(:，1)，th23(:，3），th23(:,1),th23(:,1）)
axis(［0，360，0，360］)
grid on
title（’角位移线图’）
xlabel('曲柄转角th（1)'）
ylabel（’从动件转角th(2)，th(3)’）
text（300，50，'角th(2)'）
text（300，150，’角th（3)'）
text（200,200，'角th（1）'）

2、角速度的M函数：
function y=rrrvel(x）
% Iput parameters
％
％ x（1）=theta-1
% x（2）=theta—2
% x（3）=theta-3
％ x(4)=dtheta—1
% x（5）=r1
% x（6)=r2
％ x（7)=r3
％
％ Output parameters
％
％ y（1)=dtheta—2
％ y（2)=dtheta—3
%
A=[-x(6)*sin(x（2）） x（7）＊sin（x(3）)；
x（6）＊cos（x(2)) —x(7）＊cos（x(3））］；
B=［x(5)*sin（x（1））;—x(5）＊cos（x（1)）]*x(4)；
y=inv（A)＊B
ω2，ω3
r1旋转360°时ω2，ω3的M文件程序：
th（1)=0；
dr=pi/180；
dth(1)=10;
dth(2）=44。0486＊dr；
dth(3）=96。6654*dr;
dth1=1＊dr；
r(1)=400；
r（2）=1000;
r(3)=700;
r（4）=1200;
for j=1：360
x=rrrvel([th（1），dth（2）,dth(3），dth(1）,r(1），r（2）,r（3）］）；
y=jweiyi（［th(1)，dth(2),dth(3),r（1),r（2）,r(3）,r(4）］)；
dth（2)=y(1）;
dth(3)=y（2）；
th（1）=th（1)+dth1;
dth23(j，：)=[th(1)/dr,x（1）/dr,x（2)/dr];
end
figure（2）
plot（dth23(：,1)，dth23（：,2)，dth23(：,1)，dth23(：，3））
axis（［0,360，-400,400]）
grid on
title（'角速度线图'）
xlabel(’曲柄旋转角度'）
ylabel（'从动件角速度dth(2），dth(3)（rad/s)’)
text（250,200，'连杆2角速度')
text(50，100，'摇杆3角速度’)

3、r1旋转360°时，α2，α3程序M文件:
r（1)=400;
r（2）=1000;
r(3)=700;
r(4)=1200;
th(1）=0；
dr=pi/180；
dth（2）=44.0486＊dr;
dth(3）=96。6654*dr;
dth（1)=10;
dth2=1＊dr;
for k=1：360
x=rrrvel([th（1），dth（2)，dth（3）,dth(1),r(1）,r（2),r(3)]）;
y=jweiyi（[th(1),dth（2），dth(3）,r（1）,r（2）,r(3），r4）]);
dth(2)=y(1);
dth(3）=y(2);
w1=x（1);
w2=x(2）;
C=［-r（2)＊sin（dth（2）)，r（3）＊sin（dth(3)）;r(2）＊cos(dth（2)),—r(3)＊cos(dth(3））］；
D=-[-w1*r（2)＊cos（dth（2))，w2＊r(3)＊cos(dth(3));—w2＊r(2)＊sin（dth(2)），w2*r(3)＊sin(dth（3))]…
*[w1;w2］+10*［10＊r（1)＊cos(th(1)）;10＊r（1)*sin（th(1）)];
a=inv（C）*D；
a2(k)=a（1)；
a3（k)=a(2）;
jiao(k）=th（1)；
th（1)=th(1）+dth2；
end
figure（3)
plot（jiao/dr,a2，jiao/dr，a3）
axis（［0,360，-100，100]）
grid on
title（'角加速度线图'）
xlabel('曲柄旋转转角’)
ylabel（’从动件角加速度（rad/s^2)'）
text(200，40，'连杆2的角加速度’)
text（150，-40，'摇杆3的角加速度'）

4、C点的加速度函数：
C点的M函数文件：
function y=cdjiasudu(x）
％ x(1）=r（1);
％ x(2)=r(2);
% x（3）=th（1）；
％ x(4)=th(2）；
％ x(5)=w1;
％ x(6）=w2；
% x（7)=a（1)；
y=[—x（1）*sin（x(3)）,-x（2)＊sin（x(4））；x（1）*cos(x(3)），x（2）＊cos(x（4））］＊［0；x（7)］+［—x(1)＊x（5）＊cos(x（3))，—x(2）＊x（6)＊cos（x（4)）；-x(1)＊x（5）＊sin（x（3）)，—x(2)*x(6)＊sin(x（3）)]*［x(5)^2x(6)^2］；
r1旋转360°时, C点的加速度的M文件程序：
r（1）=400；
r(2）=1000；
r(3）=700;
r（4）=1200；
th(1)=0;
dr=pi/180;
dth(2)=44.0486*dr；
dth(3）=96.6654*dr;
dth(1)=10；
dth2=1*dr；
for k=1：360
x=rrrvel（[th（1），dth（2),dth（3),dth(1），r（1)，r（2)，r(3)］)；
y=jweiyi（［th（1），dth（2)，dth(3)，r(1)，r（2），r（3）,r(4）]）；
dth(2)=y（1）；
dth（3）=y（2）;
w1=x(1）;
w2=x(2)；
C=[—r(2)＊sin(dth(2）)，r（3)*sin(dth（3))；r(2)＊cos（dth(2)）,-r（3）＊cos（dth（3））]；
D=-[-w1*r(2)*cos(dth(2））,w2*r（3）*cos（dth(3));—w2*r(2)＊sin(dth(2)），w2*r(3)＊sin（dth（3）)]＊[w1;w2］+10＊［10＊r（1)*cos(th（1))；10＊r（1)*sin(th（1)）];
a=inv（C)D;
a2=a（1)；
ac=cdjiasudu([r（1),r（2)，th（1)，dth（2),w1,w2,a2］）;
th(1）=th（1）+dth2;
acx（k）=ac（1）；
acy（k)=ac(2）;
jiao（k）=th（1）;
end
plot(jiao/dr,acx,jiao/dr，acy)
axis（［0，360,—300000,300000］）
grid on
title(’C点加速度线图’)
xlabel（’曲柄旋转角度’）
ylabel（’C点x，y方向的加速度'）
text(300，100000，’x方向')
text(300,—100000,’y方向'）