计算方法上机作业

1、精选优质文档-倾情为你奉上计算方法第四次上机报告2.用欧拉方法解初值y=10x(1-y) 0<=x<=1Y(0)=0取步长h=0.1,保留5位有效数字，并与准确解相比较分析：该题目考察欧拉方法解初值问题程序如下：function Heun(a,b,y0,n)h=(b-a)/n;x=a:h:b;y=y0*ones(1,n+1);for j=2:n+1 yp=y(j-1)+h*f(x(j-1),y(j-1); yc=y(j-1)+h*f(x(j),yp); y(j)=1/2*(yp+yc);endfor k=1:n+1 fprintf('x%d=%fty%d=%fn',

2、k-1,x(k),k-1,y(k);endfunction z=f(xx,yy)z=10*xx*(1-yy);运行结果：>> Heun(0,1,0,10)x0=0.y0=0.x1=0.y1=0.x2=0.y2=0.x3=0.y3=0.x4=0.y4=0.x5=0.y5=0.x6=0.y6=0.x7=0.y7=0.x8=0.y8=0.x9=0.y9=0.x10=1.y10=0.>>分析： 该结果与准确结果相比比较接近，但是有一定的误差。6.用四阶龙格库塔公式解第三题中的初值问题，取步长h=0.2，保留五位有效数字。题目目的分析：该题考查四阶龙格-库塔方法和改进欧拉方法求解

3、精确度问题。程序：改进欧拉法：function Heun(a,b,y0,n)h=(b-a)/n;x=a:h:b;y=y0*ones(1,n+1);for j=2:n+1 yp=y(j-1)+h*f(x(j-1),y(j-1); yc=y(j-1)+h*f(x(j),yp); y(j)=1/2*(yp+yc);endfor k=1:n+1 fprintf('x%d=%fty%d=%fn',k-1,x(k),k-1,y(k);end四阶龙格-库塔方法：function odeRK4(a,b,y0,n)h=(b-a)/n;x=a:h:b;y=y0*ones(1,n+1);for j=

4、2:n+1 k1=f(x(j-1),y(j-1); k2=f(x(j-1)+h/2,y(j-1)+h/2*k1); k3=f(x(j-1)+h/2,y(j-1)+h/2*k2); k4=f(x(j-1)+h,y(j-1)+h*k3); y(j)=y(j-1)+h/6*(k1+k4)+h/3*(k2+k3);endfor k=1:n+1 fprintf('x%d=%fty%d=%fn',k-1,x(k),k-1,y(k);endfunction z=f(xx,yy)z=-yy;结果：>> Heun(0,1.0,1,10)x0=0.y0=1.x1=0.y1=0.x2=0

5、.y2=0.x3=0.y3=0.x4=0.y4=0.x5=0.y5=0.x6=0.y6=0.x7=0.y7=0.x8=0.y8=0.x9=0.y9=0.x10=1.y10=0.>>>> odeRK4(0,1.0,1,5)x0=0.y0=1.x1=0.y1=0.x2=0.y2=0.x3=0.y3=0.x4=0.y4=0.x5=1.y5=0.>>结果分析：由程序结果得改进欧拉法和四阶龙格-库塔法运行结果上在小数第二位后才开始有所不同，但结果相差不大，精确度都是挺高的12.将二阶方程y''-5y+6y=0 (y0=1,y'0=-1)化为一阶

6、方程组，取h=0.1,用四阶龙格-库塔法求有y（0.2）的近似值，保留5位有效数字。解：目的分析： 该题目是要将二阶方程化为一解方程组，并用四阶龙格-库塔法求近似值。程序如下：公式部分：function z=f1(xx,yy)z=5*xx-6*yy;四阶龙格-库塔法：function odeRK4_1(a,b,y0,n)h=(b-a)/n;x=a:h:b;y=y0*ones(1,n+1);z=-ones(1,n+1);for j=2:n+1 k1=z(j-1); l1=f1(z(j-1),y(j-1); k2=z(j-1)+h/2*l1; l2=f1(z(j-1)+h/2*l1,y(j-1)+h/2*k1); k3=z(j-1)+h/2*l2; l3=f1(z(j-1)+h*l2/2,y(j-1)+h/2*k2); k4=z(j-1)+h*l3; l4=f1(z(j-1)+h*l3,y(j-1)+h*k3); y(j)=y(j-1)+h*(k1+2*k2+2*k3+k4)/6; z(j)=z(j-1)+h*(l1+2*l2+2*l3+l4)/6;endfor k=1:n+1 fprintf('x%d=%fty%d=%fn',k-1,x(k)