機(jī)械機(jī)構(gòu)計算機(jī)輔助分析

1、機(jī)械機(jī)構(gòu)計算機(jī)輔助分析專業(yè)：農(nóng)業(yè)生物環(huán)境與能源工程姓名：遲明路v二分法原理 通過每次把f(x)的零點(diǎn)所在小區(qū)間收縮一半的方法，使區(qū)間的兩個端點(diǎn)逐步迫近函數(shù)的零點(diǎn)，以求得零點(diǎn)的近似值，這種方法叫做二分法。 二分法也稱對分法或區(qū)間套法，是逐次逼近一個實(shí)根的簡單方法。 02xx21.用二分法求解方程二分法程序框圖 Matlab函數(shù)M文件如下：vfunction x=erfen(fx,xa,xb,n,eps)vx=xa;fa=eval(fx);vx=xb;fb=eval(fx);v disp( n xa xb xc fc );vfor i=1:n xc=(xa+xb)/2;x=xc;fc=eval(f

2、x);v X=i,xa,xb,xc,fc;vdisp(X),vif fc*fa0vxb=xc;v else v xa=xc;v endv if (xb-xa)/2eps,break,endvendvfx=x2-x-2;verfen(fx,-1.5,0.5,20,10(-5)vans = -1.0000vfx=x2-x-2;verfen(fx,1.5,2.5,20,10(-5)vans = 2.0000v牛頓法原理 牛頓法是解非線性方程（組）的一種重要迭代方法。牛頓法的一般公式：)( / )f(n1nnnxfxxx 牛頓法也常稱作切線法，其求單實(shí)根的收斂速度較快，但求重實(shí)根的收斂速度較快。2、用

3、牛頓法求解方程0e-ex牛頓法程序框圖 Matlab函數(shù)M文件如下：vfunction x=newton(fx,dfx,x0,eps,N)vx=x0; f0=eval(fx); df0=eval(dfx);vn=0;v disp( n xn xn+1 fn+1 );vwhile n=Nv x1=x0-f0/df0;v x=x1;f1=eval(fx);v X=n,x0,x1,f1;v disp(X);v if abs(x0-x1)epsv fprintf(The procedure was successful.)v returnv else v n=n+1;v x0=x1;f0=f1;v e

4、ndvendvif n=N+1v fprintf(the method failed after N iterations. )v endvfx=exp(1)-exp(x);vdfx=-exp(x);vnewton(fx,dfx,0.5,10(-5),100) vThe procedure was successful.vans = 1.0000 梯形法原理 給定步長，確定梯形的高之后，求出所有梯形的面積之和，再根據(jù)定積分的幾何意義可知曲頂高梯形的面積就是定積分的值。復(fù)化梯形公式為：niiixfxfhfT112)()()(nabh/ )( 其中,n是區(qū)間的等分?jǐn)?shù)。ba,3、用梯形法求解10si

5、nxdx梯形法程序框圖 Matlab函數(shù)M文件如下：vfunction y=T_quad(a,b,n)v h=(b-a)/n;v x1=a;v S=0;v for i=1:nv T=1/2*h*(sin(x1)+sin(x1+h);v S=S+T;v x1=x1+h;vendv SvT_quad(0,1,100)vS = 0.4597v歐拉法原理 歐拉（Euler）法（亦稱切線法）是基于查分法的最簡單的顯式單步方法。盡管其精度較低，但算法簡單，起源最早，也是其它較復(fù)雜數(shù)值解法的一個很好的向?qū)А?歐拉法公式：),(1jjjjtyhfyy) 1,.,2 , 1 , 0(nj4、歐拉法和改進(jìn)求歐拉法

6、求：微分方程的解y(0)=1,h=0.2,區(qū)間0,1。yx2y y歐拉法程序框圖 Matlab函數(shù)M文件如下：vfunction x,y=euler(fun,x0,y0,xfine,h) v N=(xfine-x0)/h;v x(1)=x0;y(1)=y0;v for n=1:Nv x(n+1)=x(n)+h;v y(n+1)=y(n)+h*feval(fun,x(n),y(n);vendvfun=inline(y-2*x/y,x,y);vx,y=euler(fun,0,1,1,0.2)vx = 0 0.2000 0.4000 0.6000 0.8000 1.0000vy =1.0000 1.

7、2000 1.3733 1.5315 1.6811 1.8269v改進(jìn)歐拉法：改進(jìn)歐拉法程序框圖 Matlab函數(shù)M文件如下：vfunction x,y=euler_g(fun,x0,y0,xfine,h) v N=(xfine-x0)/h;v x(1)=x0;y(1)=y0;v for n=1:Nv x(n+1)=x(n)+h;v ym=y(n)+h*feval(fun,x(n),y(n);v y(n+1)=y(n)+h/2*(feval(fun,x(n),y(n)v+feval(fun,x(n+1),ym);vendvfun=inline(y-2*x/y,x,y);vx,y=euler_g

8、(fun,0,1,1,0.2)vx = 0 0.2000 0.4000 0.6000 0.8000 1.0000vy = 1.0000 1.1867 1.3483 1.4937 1.6279 1.7542v黃金分割法程序框圖輸入?yún)^(qū)間a,b，精度a1=b-0.618*(b-a); y1=f(a1);a2=a+0.618*(b-a); y2=f(a2);a2-a1y1停a=0.5(a1+a2)b=a2;a2=a1;y2=y1;a1=b-0.618*(b-a)y1=f(a1);a=a1;a1=a2;y1=y2;a2=a+0.618*(b-a)y2=f(a2);yynv 5、用黃金分割法求解方程f(x

9、)=3x4-16x3+30 x2-24x+8在區(qū)間-10，10，=0.01 v Matlab函數(shù)M文件如下：vfunction y=gold(f,a,b,eps)v while abs(b-a)epsv x1=a+0.618*(b-a);v x2=a+0.382*(b-a);v f1=3*x14-16*x13+30*x12-24*x1+8 ;v f2=3*x24-16*x23+30*x22-24*x2+8 ;v if f1f2v b=x1;x1=x2;v x2=a+0.382*(b-a);v elsev a=x2;x2=x1;v x1=a+0.618*(b-a); v endv endvx=(

10、x1+x2)/2運(yùn)行結(jié)果如下：vf=3*x4-16*x3+30*x2-24*x+8;vgold(f,-10,10,0.01)vx = 1.9984 6、斐波那契法求：f(x)=x+20/x在0.2，1 上的極值，=0.001 vfunction y=fibonacci(a,b,eps,n)v F(1)=1;F(2)=1;v for i=3:n v F(i)=F(i-1)+F(i-2);v x1=a+F(i-1)/F(i)*(b-a);v x2=a+F(i-2)/F(i)*(b-a);v while abs(x2-x1)epsv f1=x1+20/x1;v f2=x2+20/x2; v if f

11、1f2v b=x1;x1=x2;v x2=a+F(i-2)/F(i)*(b-a);v elsev a=x2;x2=x1;v x1=a+F(i-1)/F(i)*(b-a);v endv endv endv x=(x1+x2)/2vfibonacci(0.2,1,0.001,50)vx =v 0.99897、用進(jìn)退法求f(x)= 3x4-16x3+30 x2-24x+8 ，a0=0.5，h=0.5，=0.001 Matlab函數(shù)M文件如下：vfunction y=jintui(a0,h,eps)v x1=a0;v x2=x1+h; v i=1;j=1; vwhile(abs(x2-x1)eps)v

12、 f1=3*x14-16*x13+30*x12-24*x1+8;v f2=3*x24-16*x23+30*x22-24*x2+8;v if(f1f2)v x1=x2;v x2=x1+2i*h;v i=i+1;v elsev x2=x1;v x1=x2-h/4j;v j=j+1;v endv endvx=(x1+x2)/2vf= 3*x4-16*x3+30*x2-24*x+8v運(yùn)行結(jié)果如下：vjintui(0.5,0.5,0.001)vx =v 1.9998vf =v 3.5751e-007v8、質(zhì)點(diǎn)在振動分離篩上的分析v(1)推導(dǎo)出k1,k2,k3。v(2)分析=-25，-10，0，10，25

13、，推導(dǎo)出=030的變化范圍。左移右移躍起的k的變化范圍。繪制=030時，k的變化規(guī)律。 v（1）篩面上任一點(diǎn)位移： R為OA的長度。v速度方程：v加速度方程：v物料在篩面上受到的牽連慣性力為： 當(dāng)曲柄在左半軸周（2,3象限)時，其慣性力的方向指向x的負(fù)方向，摩擦力的方向指向D。當(dāng)物料向下移動且沒有躍起時，則： 為物料相對篩面的加速度。 整理得：v由于-1 coswt 0,可得：0tansin)cos(cos0cos)sin(cosF22NNrFamgawtmRwmaamgawtmRwwtRwxcos2 wtRwxsinwtRw cosS2wtRcosxra0)sin(g)cos(cos)sin

14、(coscos22aawtRwawtmRwamgFN12)cos()sin(kaagRwv當(dāng)曲柄在右半軸周（1,4象限)時，慣性力的方向指向x的正方向，摩擦力的方向指向C。當(dāng)物料向上移動且沒有躍起時，則：v由于-1 coswt 0,可得：v物料躍起條件：v當(dāng) 時，由 式可以看出，只有coswt 時，只有coswt0時，才有可能 ，即物料在曲柄轉(zhuǎn)至右半周時躍起。v（2 2）c clear; lear; vclcclc; ;va=input(a=input(輸入角度：輸入角度：); );vb b=a=a* *pi/180pi/180; ; vc c=linspace(0.0001,30=linsp

15、ace(0.0001,30* *pi/180,100); pi/180,100); vd d=30=30* *(pi/180); (pi/180); vk1=sin(k1=sin(d d- -b b)./cos()./cos(c c- -b b+ +d d); ); vk2=sin(k2=sin(d d+ +b b)./cos()./cos(c c- -b b- -d d); );vk3=cos(k3=cos(b b)./abs(sin()./abs(sin(c c- -b b););0FNa322)sin(cos0cos1cos)sin(cosRkaagRwwtagawtwc1=c*180/

16、pi;plot(c1,k1,-b,c1,k2,-g,c1,k3,-r);xlabel(c);ylabel(k);title(k1,k2,k3曲線圖曲線圖);legend(k1 ,k2 ,k3 ); axis(0 30 0 8); v9（1）對于空間Vx,Vy,Vz的表達(dá)式。v（2）聯(lián)合收割谷物清選（小石頭，飽滿谷物，不成熟谷物，短秸桿）.其中V0=5,=10m/s,=10清涼室高度是0.3m,u=0-30m/s，速度分別是25 m/s， 16 m/s ，14 m/s的位移曲線。v（1）氣流的速度記為u，VL記為質(zhì)點(diǎn)的飄浮速度，則質(zhì)點(diǎn)相對于氣流的相對速度v的表達(dá)式為：v轉(zhuǎn)換成三維微分方程：uVv

17、FgmdtVdmgVuVVVVgdtdVVuVVVVgdtdVuVuVVVVgdtdVzxzyLzyxzyLyxxzyLx222222222222)()()()(v質(zhì)點(diǎn)初始速度v0在x、y、z方向上的分量：v將時間分為若干小等分，每等分為t，求出t時間后的速度，再依據(jù)已知速度求出下一個t后的速度，以此類推，可以求所有點(diǎn)的速度。v 令：sinsincoscoscos000000vvvvvvzyxgVuVVVVgVVVfVuVVVVgVVVfuVuVVVVgVVVfzxzyLzyxyxzyLzyxxxzyLzyx222232222222221)(),()(),()()(),(),(),(2),()

18、,(2),(),(2t111111111111332211iiiiiiiiiiiiiiiiiiiiiiiizyxzyxzzzyxzyxyyzyxzyxxxVVVfVVVftVVVVVfVVVftVVVVVfVVVfVVv在微小的t時間內(nèi),其速度可看作初始速度和終結(jié)速度的平均值,物料的位移量為：v質(zhì)點(diǎn)的位移：v（2）氣流分別為25m/s，16m/s，14m/s的位移曲線iiiiiiiiizzzyyyxxx1112/)(2/)(2/)(111tVVztVVytVVxiiiiiizziyyixxifor j=1:1:4 u=0; if j=1 xL=24; elseif j=2 xL=10; els

19、eif j= 3 xL=6.5; else xL=3; end clear; clc; w=30;H=0.3;L=1;v0=5;a=10*pi/180; b=a;g=9.8; t=0.001;h=w/2000; for i=1:1:2000 vx=v0*cos(a); vy=-v0*sin(a); u=u+h; x1=0; y1=0; while x1=0.3 & y1=1 f1=-g/xL2*sqrt(vx+u*sin(b)2+(vy-u*cos(b)2)*(vx+u*sin(b)+g; f2=-g/xL2*sqrt(vx+u*sin(b)2+(vy-u*cos(b)2)*(vy-u

20、*cos(b); vx1=vx+t*f1; vy1=vy+t*f2; f3=-g/xL2*sqrt(vx1+u*sin(b)2+(vy1-u*cos(b)2)*(vx1+u*sin(b)+g; f4=-g/xL2*sqrt(vx1+u*sin(b)2+(vy1-u*cos(b)2)*(vy1-u*cos(b); vx2=vx+t/2*(f1+f3); vy2=vy+t/2*(f2+f4); dx=(vx2+vx)*t/2; dy=(vy2+vy)*t/2; vx=vx1; vy=vy1; x2=x1+dx; y2=y1+dy; x1=x2; y1=y2; end if j=1 v1(i)=y2

21、; elseif j=2 v2(i)=y2; elseif j=3 v3(i)=y2; else v4(i)=y2; end endendwhile 1 fengsu=input(請輸入風(fēng)速請輸入風(fēng)速:); if fengsu0 break; else num=input(請輸入輸出參數(shù)個數(shù)請輸入輸出參數(shù)個數(shù):); n=round(fengsu/h); if n=0n=1; end out(1)=v1(n); out(2)=v2(n); out(3)=v3(n); out(4)=v4(n); if num=1 location=out(1) end if num=2 location=out(

22、1:2) end if num=3 location=out(1:3) end if num=4 location=out(1:4) end end end運(yùn)行結(jié)果請輸入風(fēng)速:25請輸入輸出參數(shù)個數(shù):4location = -0.0311 0.0786 0.3931 1.0029請輸入風(fēng)速:16請輸入輸出參數(shù)個數(shù):4location = -0.0423 0.0008 0.0884 1.0061請輸入輸出參數(shù)個數(shù):4location = -0.0442 -0.0111 0.0520 1.0023 10、曲柄滑塊機(jī)構(gòu)運(yùn)動學(xué)分析clear;clc;w=input(角速度(轉(zhuǎn)數(shù)):r=);R=input(OA桿長度:R=);L=in