Титул

Задание

Найти приближенное значение интеграла заданной функции f(x) на отрезке [a, b] по формулам прямоугольников, трапеции, Симпсона, Монте-Карло при делении отрезка на 2000 равных частей, произвести оценку погрешности методов интегрирования и сравнить точность полученных результатов: составить функцию, возвращающую значение интеграла на основе формулы метода Монте-Карло. Сравнить результаты, полученные разными методами.

№ варианта	f(x)	[a, b]
7		[1,5;2,5]

Решение

Создаем файл fun.m, в котором происходит вычисление значения функции:

function z=fun(x)

z=x.^2.*tan(x/2);

В командном окне программы МАТLAB производим вычисления подынтегральной функции на отрезке [a, b] с делением на 2000 равных частей:

>> a=1.5;

>> b=2.5;

>> N=2000;

>> i=1:N;

>> dx=(b-a)/(N-1);

>> x=a:dx:b;

>> y=fun(x);

Вычисляем значение интеграла методом правых прямоугольников:

>> m=2:N;

>> y1(m-1)=y(m);

>> Fr=sum(y1)*dx

Fr =

7.5269

Вычисляем значение интеграла методом левых прямоугольников:

>> m=1:N-1;

>> y1(m)=y(m);

>> Fl=sum(y1)*dx

Fl =

7.5185

Вычисляем значение интеграла методом трапеций:

>> s=0;

>> for i=2:N-1

s=s+y(i);

end;

>> Ft=(0.5*y(1)+s+0.5*y(N))*dx

Ft =

7.5227

Вычисляем значение интеграла методом трапеций:

>> s=0;

>> for i=2:N-1

s=s+y(i);

end;

>> Ft=(0.5*y(1)+s+0.5*y(N))*dx

Ft =

7.5227

Вычисляем значение интеграла методом Симпсона:

>> s=0;

>> for i=2:N-1

if i-2*ceil(i/2)==0

k=4;

else

k=2;

end;

s=s+k*y(i);

end;

>> Fs=(y(1)+s+y(N))*dx/3

Fs =

7.5196

Создаем файл mont_karlo.m, в котором будет реализована функция, возвращающая значение интеграла на основе формулы метода Монте-Карло:

function z=mont_karlo(f,a,b,N)

Ymin=0;

Ymax=feval(f,b);

x=a+(b-a)*rand(N,1);

y=Ymin+(Ymax-Ymin)*rand(N,1);

s=0;

for i=1:N

if y(i)<=feval(f,x(i))

s=s+1;

end;

end;

z=s*(b-a)*(Ymax-Ymin)/N;

end

Вычисляем значение интеграла с помощью метода Монте-Карло:

>> Fmk=mont_karlo('fun',a,b,N)

Fmk =

7.0631

Найдем значение интеграла с помощью функции quad, с помощью которого и будем оценивать точность полученных результатов:

>> FFt=quad('fun',a,b)

FFt =

7.5227

Определим точность каждого метода:

>> abs(Fr-FFt)/FFt

ans =

5.5583e-004

>> abs(Fl-FFt)/FFt

ans =

5.5560e-004

>> abs(Ft-FFt)/FFt

ans =

1.1442e-007

>> abs(Fs-FFt)/FFt

ans =

4.1668e-004

>> abs(Fmk-FFt)/FFt

ans =

0.0611

Как видно по полученным значением, наиболее точный результат получается по формуле трапеций.