牛头刨床运动分析机械原理剖析Word格式.docx
《牛头刨床运动分析机械原理剖析Word格式.docx》由会员分享，可在线阅读，更多相关《牛头刨床运动分析机械原理剖析Word格式.docx（13页珍藏版）》请在冰豆网上搜索。
（5）
其中
为刨刀的水平速度，
为滑块2相对于杆3的速度。
由于每个
对应的
已求出，方程组式（5）的系数矩阵均为常数，采用按列选主元的高斯消去法可求解（式5）可解得角速度ω3、ω4、
加速度方程
把（5）对时间求导得矩阵式：
（6）
同样采用按列选主元的高斯消去法可求解（6）可得角加速度
%主程序开始
clear;
clc;
l1=180;
%L1=lab
l3=960;
%l3=lCD
l4=160;
%l4=lED
h=900;
h1=460;
h2=110;
du=180/pi;
omega1=1;
alpha1=0;
theta1=linspace（0,35*pi/18,36）;
%定义常量和已知参数，l1代表杆1的长度，l3代表杆3的长度，l4代表杆4的长度,h表示EG的长度，h1表示AE的竖直距离,h2表示AE的水平距离，theta1表示角θ1的不同值。
theta3=zeros（1,36）;
theta4=zeros（1,36）;
s3=zeros（1,36）;
s5=zeros（1,36）;
test=zeros（1,36）;
vBe=zeros（1,36）;
vc=zeros（1,36）;
omega1=ones（1,36）;
omega3=zeros（1,36）;
omega4=zeros（1,36）;
aBe=zeros（1,36）;
ac=zeros（1,36）;
alpha1=zeros（1,36）;
alpha3=zeros（1,36）;
alpha4=zeros（1,36）;
A=zeros（4,4）;
dA=zeros（4,1）;
%定义最终的结果数据，当θ1取不同值时,theta3表示θ3的值，theta4表示θ4的值，s3表示BD的长度，s5表示GC的长度，vBe表示B点在杆3上运动的速度，vc表示杆5的运动速度，即牛头刨刀的速度，omega3表示杆3的转动角速度，omega4表示杆4的转动角速度，aBe表示B点在杆3上运动的角加速度，ac表示杆5的加速度，即牛头刨刀的加速度，alpha3表示杆3的角加速度，alpha4表示杆4的角加速度，矩阵A,dA表示线性方程组的系数矩阵
i=0;
%i为循环变量，在循环结构中使用
symsTHETA1THETA4%定义符号变量，为以下计算做准备
fun1=（（h1+l1*sin（THETA1）-l4*sin（THETA4））^2+（h2+l1*cos（THETA1）-l4*cos（THETA4））^2）*（l4^2*sin（THETA4）^2+h^2-2*h*l4*sin（THETA4））-l3^2*（h1+l1*sin（THETA1）-l4*sin（THETA4））^2;
%定义迭代法中要求解的关于THETA4的方程。
x0=0;
%定义在牛顿迭代法中的变量THEA4的初值。
fori=1:
36%用循环结构求当theta1取不同值时，theta3值。
fun2=subs（fun1,THETA1,theta1（i））;
%把不同的THETA1的值代入要求解的方程
[theta4（i）,EA,it]=NEWTON（fun2,'
THETA4'
x0,0.0001,1000）;
%用牛顿迭代法求得THEATA4，并赋值到theta4的数组中
x0=theta4（i）;
%把这次计算的解作为下一次计算的初值。
end
36
%用循环结构求当theta1的值取不同值时，theta3、s3、s5的取值。
因为theta3的值可能的取值范围为[0,π]，对theta3求解时应分以下两种情况讨论
ifsign（h2+l1*cos（theta1（i））-l4*cos（theta4（i）））>
0%theta3<
π/2
theta3（i）=asin（（h-l4*sin（theta4（i）））/l3）;
elsetheta3（i）=pi-asin（（h-l4*sin（theta4（i）））/l3）;
%theta3>
end
test（i）=h1+l1*sin（theta1（i））-l4*sin（theta4（i））;
s5（i）=l4*cos（theta4（i））+l3*cos（theta3（i））;
s3（i）=（h1+l1*sin（theta1（i））-l4*sin（theta4（i）））/sin（theta3（i））;
end
%用循环结构求当theta1的值取不同值时vBe、omega3、omega4、vc的值。
A（1,1）=cos（theta3（i））;
A（1,2）=-s3（i）*sin（theta3（i））;
A（1,3）=-l4*sin（theta4（i））;
A（2,1）=sin（theta3（i））;
A（2,2）=s3（i）*cos（theta3（i））;
A（2,3）=l4*cos（theta4（i））;
A（3,2）=-l3*sin（theta3（i））;
A（3,3）=-l4*sin（theta4（i））;
A（3,4）=-1;
A（4,2）=l3*cos（theta3（i））;
A（4,3）=l4*cos（theta4（i））;
dA（1,1）=-omega1（1,1）*l1*sin（theta1（i））;
dA（2,1）=omega1（1,2）*l1*cos（theta1（i））;
x=gauss（A,dA）;
%用按列选主元的高斯消去法求解
vBe（i）=x
（1）;
omega3（i）=x
（2）;
omega4（i）=x（3）;
vc（i）=x（4）;
%把求得的结构赋值给各物理量
36%用循环结构求当theta1的值取不同值时aBe、alpha3、alpha4、vc的值。
A（1,3）=-l4*sin（theta4（i））;
A（2,1）=sin（theta3（i））;
A（2,2）=s3（i）*cos（theta3（i））;
dA（1,1）=-omega3（i）*sin（theta3（i））*vBe（i）*2-s3（i）*omega3（i）^2*cos（theta3（i））-l4*omega4（i）^2*cos（theta4（i））-l1*cos（theta1（i））;
dA（2,1）=omega3（i）*cos（theta3（i））*vBe（i）*2-s3（i）*omega3（i）^2*sin（theta3（i））-l4*omega4（i）^2*sin（theta4（i））-l1*sin（theta1（i））;
dA（3,1）=-l3*omega3（i）^2*cos（theta3（i））-l4*omega4（i）^2*cos（theta4（i））;
dA（4,1）=-l3*omega3（i）^2*sin（theta3（i））-l4*omega4（i）^2*sin（theta4（i））;
%构造速度方程的系数矩阵
aBe（i）=x
（1）;
alpha3（i）=x
（2）;
alpha4（i）=x（3）;
ac（i）=x（4）;
%主程序结束
%出图程序
figure
（1）;
i=1:
10:
360;
%l3角位移图
subplot（2,2,1）;
plot（i,theta3*du,'
r'
）;
title（'
角位移图'
xlabel（'
曲柄转角\theta_1/\circ'
ylabel（'
角位移/\circ'
gridon;
holdon;
text（200,110,'
\theta3'
%l4角位移图
subplot（2,2,2）;
plot（i,theta4*du,'
text（150,10,'
\theta4'
%滑块2位移
subplot（2,2,3）;
plot（i,s3,'
位置'
滑块位置\s3/\circ'
毫米/\circ'
text（150,500,'
s3'
%c点位移
subplot（2,2,4）;
plot（i,s5,'
text（150,0,'
s5'
figure
（2）;
%l3角速度
plot（i,omega3,'
角速度图'
）
角速度/rad\cdots^{-1}'
\omega_3'
%l4角速度
plot（i,omega4,'
text（150,0.2,'
\omega_4'
%c点速度
plot（i,vc,'
速度图'
速度mm/s'
text（200,0,'
Vc'
%l4的速度
plot（i,vBe,'
text（200,-100,'
Vbe'
figure（3）;
%l3角加速度图
plot（i,alpha3,'
角加速度图'
角加速度/rad\cdots^{-2}'
text（200,-0.1,'
\alpha_3'
%l4角加速度图
plot（i,alpha4,'
text（200,-0.5,'
\alpha_4'
%c点加速度图
plot（i,ac,'
加速度图'
加速度/m\cdots^{-2}'
text（200,150,'
ac'
plot（i,aBe,'
text（200,100,'
Abe'
%牛顿迭代法的函数定义
function[r,ea,iter]=NEWTON（fun,x,x0,es,maxit）
%定义函数名和输入输出的参数。
输出参数为r,ea,iter。
其中，r代表方程的解，ea代表最终解r代入方程的误差值，iter代表在运算过程中迭代的次数。
fun,x,x0,es,maxit为输入参数。
其中，fun代表要求解的方程，x代表要求解的未知数名称，x0代表求解过程取的初值，es表示求解要求的精度，maxit表示最大迭代步数。
ifnargin<
3,error（'
请输入包括函数名，函数变量，变量初始值在内的至少三个参数'
）,end
4||isempty（es）,es=0.005;
5||isempty（maxit）,maxit=300;
%对调用函数时输入的参数进行检查，如果调用时输入的参数不足3个，则报错，如果用户没输入es和maxit的值，则设置默认的求解精度为0.005，默认的迭代步数为300。
iter=0;
%初始迭代次数为0。
t=sym（'
t'
'
real'
%定义符号变量，来取代输入函数的变量。
y=sym（'
y'
%定义符号变量，来取代输入函数。
y=subs（fun,x,t）;
y1=inline（y）;
%定义内联函数，表达式与y的表达式一致，方便计算。
y2=inline（t-y/diff（y,t））;
%定义内联函数，用于牛顿迭代。
r=x0;
fori=0:
maxit
%定义一个循环结构来产生迭代过程，当循环次数大于最大迭代次数时，循环结束。
r=y2（r）;
%牛顿迭代
iter=iter+1;
%每循环一次，迭代次数iter加一。
ifabs（y1（r））<
=es||iter>
=maxit,break;
%如果误差小于允许误差或循环次数大于最大迭代次数，迭代停止。
ea=y1（r）;
%把最终迭代误差赋值给ea。
%NEWTON函数定义完成。
%按列选主元的高斯消去法的函数定义
functionx=gauss（A,b）
输出参数为列向量x，即线性方程的解向量。
A,b为输入参数。
其中，A为与解向量x维数相同的方阵，b为与解向量维数相同的列向量，这个函数的作用是解线性方程组‘AX=b’。
[m,n]=size（A）;
%获得A的行数和列数，其中m代表矩阵的行数，n代表矩阵的列数。
ifm~=n,error（'
A必须是方阵'
end%检查A是否是方阵。
B=[A,b];
%把方阵A和向量b组成增广矩阵B。
fork=1:
n-1%用嵌套的循环结构进行消元。
[big,i]=max（abs（B（k:
n,k）））;
%找出B（k，k）、B（k+1，k）、B（k+2，k）…B（n,k）中的最大值。
u=i+k-1;
%u为第k列中的最大元素所在的列。
ifu~=k
a=B（k,:
B（k,:
）=B（u,:
B（u,:
）=a;
%把第u行的元素与第k行的元素位置互换。
end
fori=k+1:
n
factor=B（i,k）/B（k,k）;
B（i,k:
n+1）=B（i,k:
n+1）-factor*B（k,k:
n+1）;
%用初等行变换对第k+1行至第n行进行消元。
x=zeros（n,1）;
x（n）=B（n,n+1）;
fori=n-1:
-1:
1
x（i）=（B（i,n+1）-B（i,i+1:
n）*x（i+1:
n））/B（i,i）;
%的回代过程
%gauss函数定义完成。
位置-时间曲线
速度-时间曲线
加速度-时间曲线
升级会员 