欧拉法matlab程序.doc

1.Euler法function[x,y]=naeuler(dyfun,xspan,y0,h)x=xspan(1):h:xspan(2);y(1)=y0;forn=1:length(x)-1y(n+1)=y(n)+h*f(dyfun,x(n),y(n));endx=x ;y=y ;x1=0:0.2:1;y1=(1+2*x1).^0.5;plot(x,y,x1,y1)>>dyfun=inline( y-2*x/y );[x,y]=naeuler(dyfun,[0,1],1,0.2);[x,y]ans=01.00000.20001.20000.40001.37330.60001.53150.80001.68111.00001.82692.隐式Euler法function[x,y]=naeulerb(dyfun,xspan,y0,h)x=xspan(1):h:xspan(2);y(1)=y0;forn=1:length(x)-1y(n+1)=iter(dyfun,x(n+1),y(n),h);endx=x ;y=y ;x1=0:0.2:1;y1=(1+2*x1).^0.5;plot(x,y,x1,y1)functiony=iter(dyfun,x,y,h)y0=y;e=1e-4;K=1e+4;y=y+h*f(dyfun,x,y);y1=y+2*e;k=1;whileabs(y-y1)>ey1=y;y=y0+h*f(dyfun,x,y);k=k+1;ifk>Kerror( 迭代发散 );endend>>dyfun=inline( y-2*x/y );[x,y]=naeulerb(dyfun,[0,1],1,0.2);[x,y]ans=01.00000.20001.16410.40001.30140.60001.41460.80001.50191.00001.55613.改进Euler法function[x,y]=naeuler2(dyfun,xspan,y0,h)x=xspan(1):h:xspan(2);y(1)=y0;forn=1:length(x)-1k1=f(dyfun,x(n),y(n));y(n+1)=y(n)+h*k1;k2=f(dyfun,x(n+1),y(n+1));y(n+1)=y(n)+h*(k1+k2)/2;endx=x ;y=y ;x1=0:0.2:1;y1=(1+2*x1).^0.5;plot(x,y,x1,y1)>>dyfun=inline( y-2*x/y );[x,y]=naeuler2(dyfun,[0,1],1,0.2);[x,y]ans=01.00000.20001.18670.40001.34830.60001.49370.80001.62791.00001.7542