搜索
计算方法第四次上机报告2.用欧拉方法解初值y’=10x(1-y)0=x=1Y(0)=0取步长h=0.1,保留5位有效数字，并与准确解相比较分析：该题目考察欧拉方法解初值问题程序如下：functionHeun(a,b,y0,n)h=(b-a)/n;x=a:h:b;y=y0*ones(1,n+1);forj=2:n+1yp=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);endfork=1:n+1fprintf('x[%d]=%f\ty[%d]=%f\n',k-1,x(k),k-1,y(k));endfunctionz=f(xx,yy)z=10*xx*(1-yy);运行结果：Heun(0,1,0,10)x[0]=0.000000y[0]=0.000000x[1]=0.100000y[1]=0.050000x[2]=0.200000y[2]=0.183000x[3]=0.300000y[3]=0.362740x[4]=0.400000y[4]=0.547545x[5]=0.500000y[5]=0.705905x[6]=0.600000y[6]=0.823543x[7]=0.700000y[7]=0.901184x[8]=0.800000y[8]=0.947627x[9]=0.900000y[9]=0.973290x[10]=1.000000y[10]=0.986645分析：该结果与准确结果相比比较接近，但是有一定的误差。6.用四阶龙格—库塔公式解第三题中的初值问题，取步长h=0.2，保留五位有效数字。题目目的分析：该题考查四阶龙格-库塔方法和改进欧拉方法求解精确度问题。程序：改进欧拉法：functionHeun(a,b,y0,n)h=(b-a)/n;x=a:h:b;y=y0*ones(1,n+1);forj=2:n+1yp=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);endfork=1:n+1fprintf('x[%d]=%f\ty[%d]=%f\n',k-1,x(k),k-1,y(k));end四阶龙格-库塔方法：functionodeRK4(a,b,y0,n)h=(b-a)/n;x=a:h:b;y=y0*ones(1,n+1);forj=2:n+1k1=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);endfork=1:n+1fprintf('x[%d]=%f\ty[%d]=%f\n',k-1,x(k),k-1,y(k));endfunctionz=f(xx,yy)z=-yy;结果：Heun(0,1.0,1,10)x[0]=0.000000y[0]=1.000000x[1]=0.100000y[1]=0.905000x[2]=0.200000y[2]=0.819025x[3]=0.300000y[3]=0.741218x[4]=0.400000y[4]=0.670802x[5]=0.500000y[5]=0.607076x[6]=0.600000y[6]=0.549404x[7]=0.700000y[7]=0.497210x[8]=0.800000y[8]=0.449975x[9]=0.900000y[9]=0.407228x[10]=1.000000y[10]=0.368541odeRK4(0,1.0,1,5)x[0]=0.000000y[0]=1.000000x[1]=0.200000y[1]=0.818733x[2]=0.400000y[2]=0.670324x[3]=0.600000y[3]=0.548817x[4]=0.800000y[4]=0.449335x[5]=1.000000y[5]=0.367885结果分析：由程序结果得改进欧拉法和四阶龙格-库塔法运行结果上在小数第二位后才开始有所不同，但结果相差不大，精确度都是挺高的12.将二阶方程化为一阶方程组，取h=0.1,用四阶龙格-库塔法求有y（0.2）的近似值，保留5位有效数字。解：目的分析：该题目是要将二阶方程化为一解方程组，并用四阶龙格-库塔法求近似值。程序如下：公式部分：functionz=f1(xx,yy)z=5*xx-6*yy;四阶龙格-库塔法：functionodeRK4_1(a,b,y0,n)h=(b-a)/n;x=a:h:b;y=y0*ones(1,n+1);z=-ones(1,n+1);forj=2:n+1k1=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;endfork=1:n+1fprintf('x[%d]=%f\ty[%d]=%f\n',k-1,x(k),k-1,y(k));endend结果图：结果分析：由程序结果得四阶龙格-库塔法所得y(0.2)=0.501088所给初值问题精确解为,y(0.2)=0.5009424.两结果在数值上相差不大，精确度高。