数值分析常微分方程边值问题.pptVIP

数值分析常微分方程边值问题第一页，共十六页，2022年，8月28日 初值问题欧拉公式: yn+1 = yn+ h f(xn , yn)修改的欧拉公式(Heun方法):k1 = f(xn , yn) , k2 = f( xn+1 , yn+ h k1)(n = 0, 1, 2, ·······, N )(n = 0, 1, 2, ·······, N )2/16第二页，共十六页，2022年，8月28日 二阶Range-Kutta公式一般形式yn+1= yn+ h[c1k1+ c2k2]k1=f (xn, yn)k2=f(xn+λ2h, yn+μ21hk1)( n=0,1,····,N )c1=0.5, c2 = 0.5, λ2=1, μ21= 1 c1= 0, c2 = 1, λ2=0.5, μ21=0.5 yn+1= yn+ 0.5h[ k1+ k2] k1=f (xn, yn) k2=f(xn+h, yn+hk1) yn+1= yn+ h k2k1=f (xn, yn)k2=f(xn+.5h, yn+.5hk1)3/16第三页，共十六页，2022年，8月28日 三阶Range-Kutta公式一般形式yn+1= yn+ h[k1+4k2+k3]/6k1=f(xn,yn), k2=f(xn+0.5h, yn+0.5hk1)k3=f(xn+h, yn – hk1+2hk2)四阶Range-Kutta公式一般形式yn+1= yn+ h[k1+2k2+2k3+k4]/6k1=f(xn,yn), k2=f(xn+0.5h, yn+0.5hk1)k3=f(xn+0.5h, yn+0.5hk2), k4=f(xn+h, yn+hk3)4/16第四页，共十六页，2022年，8月28日 例1数值实验:几种不同求数值解公式的误差比较 n 10 20 30 40 h 0.2 0.1 0.0667 0.05RK4 6.862e-005 3.747e-006 7.071e-007 2.186e-007RK3 0.0012 1.529e-004 4.517e-005 1.906e-005RK2 0.0123 0.0026 0.0011 5.9612e-004Euler 0.1059 0.0521 0.0342 0.02565/16第五页，共十六页，2022年，8月28日 n=input(input n:= );f=inline(y-x.*y.^2);h=2/n;x=h:h:2;y0=1;k1=f(0,y0);k2=f(.5*h,y0+0.5*h*k1);k3=f(.5*h,y0+.5*h*k2);k4=f(h,y0+h*k3);y(1)=y0+h*(k1+2*k2+2*k3+k4)/6;for k=1:n-1 xk=x(k);yk=y(k);k1=f(xk,yk); k2=f(xk+h/2,yk+0.5*h*k1); k3=f(xk+h/2,yk+0.5*h*k2); k4=f(xk+h,yk+h*k3); y(k+1)=yk+h*(k1+2*k2+2*k3+k4)/6;endt=-1:h/2:5;yy=1./(t-1+2*exp(-t));plot(0,1,ro,x,y,ro,t,yy)y1=1./(x-1+2*exp(-x));error=max(abs(y-y1))input n:= 20error = 3.7475e-0066/16第六页，共十六页，2022年，8月28日 例2. 捕食者与被捕食者模型0≤ x ≤ 15 y1(0)=20, y2(0)=20记7/16第七页，共十六页，2022年，8月28日 x(1)=20;y(1)=20;for n=1:length(t)-1 tn=t(n);xn=x(n);yn=y(n); k11=f1(xn,yn);k21=f2(xn,yn); k12=f1(xn+.5*h*k11,yn+.5*h*k21); k22=f2(xn+.5*h*k11,yn+.5*h*k21); k13=f1(xn+.5*h*k12,yn+.5*h*k22); k23=f2(xn+.5*h*k12,yn+.5*h*k22);