计算方法上机作业

计算方法第四次上机报告

2.用欧拉方法解初值

y’=10x(1-y) 0<=x<=1

Y(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);

end

for k=1:n+1

fprintf('x[%d]=%f\ty[%d]=%f\n',k-1,x(k),k-1,y(k)); end

function z=f(xx,yy)

z=10*xx*(1-yy);

运行结果：

>> Heun(0,1,0,10)

x[0]=0.000000 y[0]=0.000000

x[1]=0.100000 y[1]=0.050000

x[2]=0.200000 y[2]=0.183000

x[3]=0.300000 y[3]=0.362740

x[4]=0.400000 y[4]=0.547545

x[5]=0.500000 y[5]=0.705905

x[6]=0.600000 y[6]=0.823543

x[7]=0.700000 y[7]=0.901184

x[8]=0.800000 y[8]=0.947627

x[9]=0.900000 y[9]=0.973290

x[10]=1.000000 y[10]=0.986645

>>

分析：

该结果与准确结果相比比较接近，但是有一定的误差。

6.用四阶龙格—库塔公式解第三题中的初值问题，取步长h=0.2，保留五位有效数字。

题目目的分析：

该题考查四阶龙格-库塔方法和改进欧拉方法求解精确度问题。

程序：

改进欧拉法：

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);

end

for k=1:n+1

fprintf('x[%d]=%f\ty[%d]=%f\n',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=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);

end

for k=1:n+1

fprintf('x[%d]=%f\ty[%d]=%f\n',k-1,x(k),k-1,y(k)); end

function z=f(xx,yy)

z=-yy;

结果：

>> Heun(0,1.0,1,10)

x[0]=0.000000 y[0]=1.000000

x[1]=0.100000 y[1]=0.905000

x[2]=0.200000 y[2]=0.819025

x[3]=0.300000 y[3]=0.741218

x[4]=0.400000 y[4]=0.670802

x[5]=0.500000 y[5]=0.607076

x[6]=0.600000 y[6]=0.549404

x[7]=0.700000 y[7]=0.497210

x[8]=0.800000 y[8]=0.449975

x[9]=0.900000 y[9]=0.407228

x[10]=1.000000 y[10]=0.368541

>>

>> odeRK4(0,1.0,1,5)

x[0]=0.000000 y[0]=1.000000

x[1]=0.200000 y[1]=0.818733

x[2]=0.400000 y[2]=0.670324

x[3]=0.600000 y[3]=0.548817

x[4]=0.800000 y[4]=0.449335

x[5]=1.000000 y[5]=0.367885

>>

结果分析：

由程序结果得改进欧拉法和四阶龙格-库塔法运行结果上在小数第二位后才开始有所不同，但结果相差不大，精确度都是挺高的

12.将二阶方程

化为一阶方程组，取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);