欧拉近似方法求常微分方程.docVIP

欧拉近似方法求常微分方程 朱翼 1、编程实现以下科学计算算法，并举一例应用之。 “欧拉近似方法求常微分方程” 算法说明： 欧拉法是简单有效的常微分方程数值解法，欧拉法有多种形式的算法，其中简单欧拉法是一种单步递推算法。其基本原理为 对简单的一阶方程的初值问题： y’=f(x,y) 其中 y(x0 )=y0 欧拉法等同于将函数微分转换为数值微分，由欧拉公式可得 yn+1 =y n+hf(x n ,y n) 程序代码： function [tout,yout]=myeuler(ypfun,t0,tfinal,y0,tol,trace) %初始化 pow=1/3; if nargin5,tol=1.e-3;end if nargin6,trace=0;end t=t0; hmax=(tfinal-t)/16; h=hmax/8; y=y0(:); chunk=128; tout=zeros(chunk,1); yout=zeros(chunk,length(y)); k=1; tout(k)=t; yout(k,:)=y.; if trace %绘图 clc,t,h,y end while (ttfinal)(t+ht) %主循环 if t+htfinal,h=tfinal-t;end % Compute the slopes f=feval(ypfun,t,y);f=f(:); %估计误差并设定可接受误差 delta=norm(h*f,inf); tau=tol*max(norm(y,inf),1.0); %当误差可接受时重写解 if delta=tau t=t+h; y=y+h*f; k=k+1; if klength(tout) tout=[tout;zeros(chunk,1)]; yout=[yout;zeros(chunk,length(y))]; end tout(k)=t; yout(k,:)=y.; end if trace home,t,h,y end % Update the step size if delta~=0.0 h=min(hmax,0.9*h*(tau/delta)^pow); end end if (ttfinal) dish(Singularity likely.) t end tout=tout(1:k); yout=yout(1:k,:); 流程图： 开始 开始 [tout,yout]=my [tout,yout]=myeuler(ypfun,t0,tfinal,y0,tol,trace) Pow=1/3 Pow=1/3 Ynargin5 Y nargin5 Ntol=1.e-3 N tol=1.e-3 Y Y nargin nargin6 trace=0 trace=0 t=t0;hmax=(tfinal-t)/16;h=hm t=t0;hmax=(tfinal-t)/16;h=hmax/8;y=y0(:);chunk=128; tout=zeros(chunk,1);yout=zeros(chunk,length(y)); k=1;tout(k)=t;yout(k,:)=y.; N N Y Y trace= 1 trace= 1 N输出t,h,y N 输出t,h,y N N ttfinal ttfinal t+ht Y Y Y Y t+htfinal t+htfinal h=tfinal-t h=tfinal-t N N f=feval(ypfun,t,y); f=feval(ypfun,t,y); f=f(:); delta=norm(h*f,inf); tau=tol*max(norm(y,inf)1.0); tau=tol*max(norm(y,inf)1.0); Y Y delta =tau delta =tau Nt=t+h; N t=t+h; y=y+h*f; k=k+l; k klength(tout) tout=[tout;zeros(chunk,1)]; tout=[tout;zeros(chunk,1)]; yout=[yout;zeros(chunk,lengt