欧拉法,梯形法,龙格库塔法

数值分析实验报告（6）

实验目的：用欧拉法、梯形法和龙格—库塔法计算例题，分析三种方法的优点和缺点，掌握规律。

实验要求：用不三种方法来解同一个例子，比较误差。

1.欧拉法：

MATLAB程序：

function e=aa(f,a,b,n,x,y0,g)

f=inline(f);

g=inline(g);

h=(b-a)/n;

y(1)=y0;x0=0;

disp(sprintf('i x(i) y(i) yx(i) wucha(i)'));

for i=1:n+1

x(i)=x0+(i-1)*h;

y(i+1)=y(i)+h*f(x(i),y(i));

yx(i)=g(x(i));

wucha(i)=abs(yx(i)-y(i));

disp(sprintf('%d %f %f %f %f',i-1,x(i),y(i),yx(i),wucha(i))); end

运行结果：

>> aa('y-2*x/y',0,1,10,0,1,'(1+2*x)^(1/2)')

i x(i) y(i) yx(i) wucha(i)

0 0.000000 1.000000 1.000000 0.000000

1 0.100000 1.100000 1.095445 0.004555

2 0.200000 1.191818 1.183216 0.008602

3 0.300000 1.277438 1.264911 0.012527

4 0.400000 1.358213 1.341641 0.016572

5 0.500000 1.435133 1.414214 0.020919

6 0.600000 1.508966 1.483240 0.025727

7 0.700000 1.580338 1.549193 0.031145

8 0.800000 1.649783 1.612452 0.037332

9 0.900000 1.717779 1.673320 0.044459

10 1.000000 1.784771 1.732051 0.052720

2.梯形方法

MATLAB程序：

function e=tixingaa(f,a,b,n,x,y0,g)

f=inline(f);

g=inline(g);

h=(b-a)/n;

y(1)=y0;x0=0;

disp(sprintf('i x(i) y(i) yx(i) wucha(i)'));

for i=1:n+1

x(i)=x0+(i-1)*h;

x(i+1)=h+x(i);

s(i+1)=y(i)+h*f(x(i),y(i));

y(i+1)=y(i)+h/2*(f(x(i),y(i))+f(x(i+1),s(i+1)));

yx(i)=g(x(i));

wucha(i)=abs(yx(i)-y(i));

disp(sprintf('%d %f %f %f %f',i-1,x(i),y(i),yx(i),wucha(i))); end

运行结果：

>> tixing('y-2*x/y',0,1,10,0,1,'(1+2*x)^(1/2)')

i x(i) y(i) yx(i) wucha(i)

0 0.000000 1.000000 1.000000 0.000000

1 0.100000 1.095909 1.095445 0.000464

2 0.200000 1.184097 1.183216 0.000881

3 0.300000 1.266201 1.264911 0.001290

4 0.400000 1.343360 1.341641 0.001719

5 0.500000 1.416402 1.414214 0.002188

6 0.600000 1.485956 1.483240 0.002716

7 0.700000 1.552514 1.549193 0.003321

8 0.800000 1.616475 1.612452 0.004023

9 0.900000 1.678166 1.673320 0.004846

10 1.000000 1.737867 1.732051 0.005817

3.三阶R-K方法：

MATLAB程序：

function e=R3(f,a,b,n,x,y0,g)

f=inline(f);

g=inline(g);

h=(b-a)/n;

y(1)=y0;x0=0;

disp(sprintf('i x(i) y(i) yx(i) wucha(i)'));

for i=1:n+1

x(i)=x0+(i-1)*h;

x(i+1)=h+x(i);

k(1)=f(x(i),y(i));

k(2)=f(x(i)+(h/2),y(i)+(h/2)*k(1));

k(3)=f(x(i)+h,y(i)-h*k(1)+2*h*k(2));

y(i+1)=y(i)+(h/6)*(k(1)+4*k(2)+k(3));

yx(i)=g(x(i));

wucha(i)=abs(yx(i)-y(i));

disp(sprintf('%d %f %f %f %f',i-1,x(i),y(i),yx(i),wucha(i))); end

运行结果：

>> R3('y-2*x/y',0,1,5,0,1,'(1+2*x)^(1/2)')

i x(i) y(i) yx(i) wucha(i)

0 0.000000 1.000000 1.000000 0.000000

1 0.200000 1.183244 1.183216 0.000028

2 0.400000 1.341729 1.341641 0.000088

3 0.600000 1.483408 1.483240 0.000169

4 0.800000 1.612727 1.612452 0.000276

5 1.000000 1.732472 1.732051 0.000421

实验结果：经过用三种方法解题可以看出，如果我们用的方法的精度越高，误差越小，我们算出的近似值越接近精确解。可是精度越高计算量也越多。