单自由度机械系统动力学——牛头刨床运动例题

PAGE\*MERGEFORMAT#单自由度机械系统动力学作业 题 快递公司问题件快递公司问题件货款处理关于圆的周长面积重点题型关于解方程组的题及答案关于南海问题 目：图1所示为一牛头刨床。各构件长度为：L1110mm,L3540mm,L4135mm；尺寸H580mm,H〔380mm。导杆3重量G3200N，质心S3位于导杆中心，导杆绕2&的转动惯量J31.1kgm。滑枕5的重量G5700N。其余构件重量均可不计。电动机型号为Y100L2-4，电动机轴至曲柄1的传动比i23.833,电动机转子及传动齿轮等折算到曲柄上的转动惯量Ji133.3kgm2。刨床的平均传动效率0.85。空行程时作用在滑枕上的摩擦阻力Ff50N，切削某工件时的切削力和摩擦阻力如图2所示。求空载启动后曲柄的稳态运动规律；求开始刨削工件的加载过程，直至稳态。图1牛头刨床图2牛头刨床加工某工件时的负载图解：运动 分析 定性数据统计分析pdf销售业绩分析模板建筑结构震害分析销售进度分析表京东商城竞争战略分析 可以用解析法列出各杆角速度、各杆质心速度的 表 关于同志近三年现实表现材料材料类招标技术评分表图表与交易pdf视力表打印pdf用图表说话 pdf 达式。但为简便起见，现调用改自课本附录I中的Matlab子程序来进行计算。图1中给出了构件和运动副的编号。先调用子程序crank分析点②的运动学参数，再调用子程序vosc进行滑块2一导杆3这一杆组的运动学分析，然后再调用子程序vguide进行小连杆4一滑枕5这一杆组的运动学分析。这一段的Matlab程序如下：crank(1,2,L(1),TH(1),W(1));vosc(2,3,4,L(3));vguide(4,5,L(4));其中：L(i)、TH(i)、W(i)分别表示第i个杆的长度、位置角、角速度。等效转动惯量和等效力矩取曲柄1为等效构件，等效转动惯量为Je32G3Vs3、2G5/V5\2(a)J1J3()~g()式中：g为重力加速度，Vs3为导杆3质心的速度，V5为滑枕的速度。等效驱动力矩可由电动机机械特性导出,设Mm、Mde分别为电动机输出力矩和等效驱(b)(c)(d)(e)(f)动力矩，两者有如下关系:MdeiMm式中i为电动机轴和曲轴间的传动比。电动机轴转速m和曲柄转速1间有如下关系：m1一i将式（b）和式（c）代入电动机机械特性2Mmabmcm可得Mdeaibi21ci3:a、b、c的值代入式（e）,得到等效驱动将传动比i的值和课本例题3.2.2中求出的系数力矩TOC\o"1-5"\h\z一-2Mde148466076.81580.261等效阻力矩M怡中只计入滑枕上的摩擦阻力Ff和切削阻力F「，以及导杆的重力G3：Mre（G3Vs3y|（Ff日沮|）/（i）（g）式中Vs3y为导杆3重心S3的y向速度。等效力矩Me为MeMdeMre（h）等效力矩和等效转动惯量均随机构位置而变化。需将曲柄运动周期分成k个等份（k可取为60），对每一机构位置计算等效力矩Me和等效转动惯量Je。（3）运动方程的求解本题属于等效力矩同时为等效构件转角和角速度的函数，而等效力矩Me的表达式中与可以分离，即可以表达为两个函数的和，其中一个等效驱动力矩Mde为角速度的函数，另一个等效阻力矩Mre为转角的函数。这样采用能量形式的运动方程求解更为简便、快速。已知等效驱动力矩Mde的表达式为式(f)，等效阻力矩Mre的表达式为式(g)，设MreMr(),为已知量。由能量形式的运动方程,对从1到2的区间，可以写出n21J2-Je22-Je11222(MdeMre)d1(i)式中,Jei为与角i相对应的位置的角速度和等效转动惯量。用梯形公式求积分，式-J22e22(i)可写为2e11~2~[Mde(1)Mre(1)Mde(2)Mre(2)](j)用式(e)和MreMr()代入，得到(Je2b2[Mr(1)Mr(2)2ab112]Je112(k)这是一个以2为未知数的元二次方程，如果1已知，2便可很容易地求出。同理，对第i个区间，即i1之间的区间，可以有如下递推公式:21BiCi0(l)式中:ABi[(Je)ibc]Ci{[Mr(i)Mr(i1)2a22-i](Je)ii}用此递推公式，当已知初始条件1时便可逐步求出各位置的角速度。计算空载启动后的稳态响应不必取初值0,为在计算中迅速收敛，可任意取一接近电动机额定角速度的初值，如取16.5rad/s，则不到两周便求出稳态解，如图3所示。可以看出，在空载下速度波动很小。开始刨削后的加载过程的初值可取空载稳态1360时的1值。加载过程如图4所示，最后得到切削时的稳态响应如图5中曲线(2)所示，可以看出，负载的波动导致了较大的速度波动。将图5、图4与图2对比，可以很清楚地看出，工作循环中的两段有切削力的部分基本上与1变化中两次降速的位置相对应。6.636.626.616.66.596.586.576.566.556.54501001502002503003504001图3空载启动后曲柄的稳态运动规律6.86.66.56.46.36.26.165.901002003004005006007008001图4开始刨削工件的加载过程6.86.76.66.56.46.36.26.165.90501001502002503003504001图5空载与切削时的稳态响应Matlab程序：[main.m]globalPVP%&点位置与速度为全局变量P=zeros(5,2);VP=zeros(5,2);P(3,2)=-0.38;P(5,2)=0.2;Je=zeros(1,61);Mre=zeros(1,61);Mre0=zeros(1,61);DeltaPhi=pi/30;%准备工作，先计算各个位置时的等效转动惯量Je，等效阻力矩Mre%因为本题等效转动惯量与等效阻力矩均只与机构位置有关，与角速度无关，设曲柄角速度为1进行计算fork=1:60crank(1,2,0.11,2*pi-(k-1)*DeltaPhi,1);W3=vosc(2,3,4,0.54);vguide(4,5,0.135);Vs3=sqrt(VP(4,1)A2+VP(4,2)A2)/2;Je(k)=133.3+1.1*W52+200/10*Vs3A2+700/10*VP(5,1)人2;if((k>=33&&k<=43)||(k>=50&&k<=59))F=9500;elseF=50;endMre(k)=(-200*VP(4,2)/2-abs(F*VP(5,1)))/0.85;%刨削工件时的阻力Mre0(k)=(-200*VP(4,2)/2-abs(50*VP(5,1)))/0.85;盼载阻力矩endJe(61)=Je(1);%第61点值与第1点值相同，只是为了方便后面的迭代计算Mre(61)=Mre(1);Mre0(61)=Mre0(1);n=0;%己录迭代次数，其实没什么用w=zeros(1,61);w(1)=6.5;w(61)=1;whileabs(w(61)-w(1))/w(61)>=1e-4ifn==0n=1;elsew(1)=w(61);%更新原点比较值n=n+1;endfork=1:60A=Je(k+1)-DeltaPhi*(-580.26);B=-DeltaPhi*6076.8;58C=-(DeltaPhi*(Mre0(k)+Mre0(k+1)+2*(-14846)+6076.8*w(k)+(-0.26)*w(k)A2)+Je(k)*w(k)A2);w(k+1)=(-B+sqrt(BA2-4*A*C))/(2*A);endendPhi=0:6:360;plot(Phi,w);%绘制空载稳态TOC\o"1-5"\h\zxlabel('\phi_1');ylabel('\omega1/(rad/s)');figure(3);%用来绘制空载与工作状态稳态对比plot(Phi,w);%绘制空载稳态xlabel('\phi_1');ylabel('\omega_1/(rad/s)');w0=w(61);w=zeros(1,121);趣加载后两个周期的数据w(1)=w0;fork=1:120sk=mod(k-1,60)+1;%sk取值范围为1~60A=Je(sk+1)-DeltaPhi*(-580.26);B=-DeltaPhi*6076.8C=-(DeltaPhi*(Mre(sk)+Mre(sk+1)+2*(-14846)+6076.8*w(k)+(-580.26)*w(k)A2)+Je(sk)*w(k)A2);w(k+1)=(-B+sqrt(BA2-4*A*C))/(2*A);endPhi=0:6:120*6;figure;plot(Phi,w);%绘制加载过程xlabel('\phi_1');ylabel('\omega_1/(rad/s)');while里的语%计算加载后达到的稳态响应，其实之前的计算值已经满足要求精度了，所以句不会执行w(1:61)=w(61:121);whileabs(w(61)-w(1))/w(61)>=1e-4|w(1)=w(61);fork=1:60A=Je(k+1)-DeltaPhi*(-580.26);B=-DeltaPhi*6076.8;C=-(DeltaPhi*(Mre(k)+Mre(k+1)+2*(-14846)+6076.8*w(k)+(-580.26)*w(k)A2)+Je(k)*w(k)A2);w(k+1)=(-B+sqrt(BA2-4*A*C))/(2*A);endendPhi=0:6:360;figure(3);holdon;plot(Phi,w(1:61));%绘制工作状态稳态，与空载对比[crank.m]functioncrank(N1,N2,R,TH,W)globalPVPVP(N1,1)=0;VP(N1,2)=0;RX=R*cos(TH);RY=R*sin(TH);P(N2,1)=P(N1,1)+RX;P(N2,2)=P(N1,2)+RY;VP(N2,1)=-RY*WVP(N2,2)=RX*W[vosc.m]function[W]=vosc(N1,N2,N3,R)globalPVPTH=posc(N1,N2,N3,R);R2=sqrt((P(N2,1)-P(N1,1))A2+(P(N2,2)-P(N1,2))A2);W=((VP(N1,2)-VP(N2,2))*cos(TH)-(VP(N1,1)-VP(N2,1))*sin(TH))/R2;VP(N3,1)=VP(N2,1)-W*R*sin(TH);VP(N3,2)=VP(N2,2)+WR*cos(TH);[posc.m]function[TH]=posc(N1,N2,N3,R)globalPTH=atan2(P(N1,2)-P(N2,2),P(N1,1)-P(N2,1));P(N3,1)=P(N2,1)+R*cos(TH);P(N3,2)=P(N2,2)+R*sin(TH);[vguide.m]functionvguide(N1,N2,R)globalPVPTH=pi-asin((P(N2,2)-P(N1,2))/R);W=-VP(N1,2)/(R*cos(TH));VP(N2,1)=VP(N1,1)-R*W*sin(TH);