常微分方程作业欧拉法与改进欧拉法.docx

1、常微分方程作业欧拉法与改进欧拉法常微分方程作业欧拉法与改进欧拉法P77 31.利用改进欧拉方法计算下列初值问题，并画出近似解的草图：代码：%改进欧拉法function Euler(t0,y0,inv,h)n=round(inv(2)-inv(1)/h;t(1)=t0;y(1)=y0;for i=1:n y1(i+1)=y(i)+h*fun(t(i),y(i); t(i+1)=t(i)+h; y(i+1)=y(i)+1/2*h*(fun(t(i),y(i)+ fun(t(i+1),y1(i+1)endplot(t,y,*r)function y=fun(t,y);y=y+1;调用：Euler(0

2、,3,0,2,0.5)得到解析解：hold on;y=dsolve(Dy=y+1,(y(0)=3),t);ezplot(y,0,2)图像：代码：function Euler1(t0,y0,inv,h)n=round(inv(2)-inv(1)/h;t(1)=t0;y(1)=y0;for i=1:n y1(i+1)=y(i)+h*fun(t(i),y(i); t(i+1)=t(i)+h; y(i+1)=y(i)+1/2*h*(fun(t(i),y(i)+ fun(t(i+1),y1(i+1)endplot(t,y,*r)function y=fun(t,y);y=y2-4*t;调用：Euler1

3、(0,0.5,0,2,0.2)图像：代码：function Euler2(t0,y0,inv,h)n=round(inv(2)-inv(1)/h;t(1)=t0;y(1)=y0;for i=1:n y1(i+1)=y(i)+h*fun(t(i),y(i); t(i+1)=t(i)+h; y(i+1)=y(i)+1/2*h*(fun(t(i),y(i)+ fun(t(i+1),y1(i+1)endplot(t,y,*r)function y=fun(t,y);y=(3-y)*(y+1);调用：Euler2(0,4,0,5,1)得到解析解：hold on;y=dsolve(Dy=(3-y)*(y+

4、1),y(0)=4,t);ezplot(y)图像：代码：function Euler2(t0,y0,inv,h)n=round(inv(2)-inv(1)/h;t(1)=t0;y(1)=y0;for i=1:n y1(i+1)=y(i)+h*fun(t(i),y(i); t(i+1)=t(i)+h; y(i+1)=y(i)+1/2*h*(fun(t(i),y(i)+ fun(t(i+1),y1(i+1)endplot(t,y,*r)function y=fun(t,y);y=(3-y)*(y+1);调用：Euler2(0,4,0,5,0.5)得到解析解：hold on;y=dsolve(Dy=

5、(3-y)*(y+1),y(0)=4,t);ezplot(y)图像：14.考虑满足初始条件（x(0)，y(0)）=(1,1)的下列方程组：选定时间步长t=0.25，n=5.用改进欧拉方法求两个方程组的近似解；（1）代码：function Euler4(t0,int,n,h)t=t0;x(1)=int(1);y(1)=int(2);for i=1:n x1(i+1)=x(i)+h*xfun(t(i),x(i),y(i); y1(i+1)=y(i)+h*yfun(t(i),x(i),y(i); t(i+1)=t(i)+h; x(i+1)=x(i)+1/2*h*(xfun(t(i),x(i),y(i

6、)+xfun(t(i+1),x1(i+1),y1(i+1); y(i+1)=y(i)+1/2*h*(yfun(t(i),x(i),y(i)+yfun(t(i+1),x1(i+1),y1(i+1);endplot(t,x,o-r)hold onplot(t,y,*-g)hold onplot(x,y)function x=xfun(t,x,y);x=y;function y=yfun(t,x,y);y=-2*x-3*y;调用函数：Euler4(0,1,1,5,0.25)图像：（2）代码：function Euler5(t0,int,n,h)t=t0;x(1)=int(1);y(1)=int(2)

7、;for i=1:n x1(i+1)=x(i)+h*xfun(t(i),x(i),y(i); y1(i+1)=y(i)+h*yfun(t(i),x(i),y(i); t(i+1)=t(i)+h; x(i+1)=x(i)+1/2*h*(xfun(t(i),x(i),y(i)+xfun(t(i+1),x1(i+1),y1(i+1); y(i+1)=y(i)+1/2*h*(yfun(t(i),x(i),y(i)+yfun(t(i+1),x1(i+1),y1(i+1);endplot(t,x,o-r)hold onplot(t,y,*-g)hold onplot(x,y)function x=xfun(t,x,y);x=y+y2;function y=yfun(t,x,y);y=-x+0.2*y-x*y+1.2*y2;调用函数：Euler5(0,1,1,5,0.25)图像：