clear all

F = @(x) x.^cos(x);

%CAŁKI
a=0;
b=10;

X = [a:0.01:b];
Y = F(X);

F_integral = integral(F, a, b);
sprintf('Wbudowana funkcja integral: %d', F_integral)

%Metoda prostokątów
n = 10;
x = linspace(a,b,n+1);
h = (b-a)/n;

F_rect = 0;
for i=1:length(x)-1
   F_rect = F_rect + F(x(i));
end
F_rect = F_rect * h;
sprintf('Metoda prostokątów: %d', F_rect)

%Metoda trapezów
F_trap = 0;
for i=2:length(x)-1
   F_trap = F_trap + F(x(i));
end
F_trap = (F_trap + (F(x(1)) + F(x(n+1)))/2) * h;
sprintf('Metoda trapezów: %d', F_trap)

%Metoda Simpsona
even = 0;
for i=2:2:n
   even = even + F(x(i));
end
odd = 0;
for i=3:2:n-1
   odd = odd + F(x(i));
end
F_simp = (F(a) + 4*even + 2*odd + F(b)) * (h/3);
sprintf('Metoda Simpsona: %d', F_simp)

%Monte Carlo
N = 10000;
x = b * rand(N, 1);
y = max(Y) * rand(N, 1);

z = y(y<=F(x));
figure(1)
fplot(F, [0 10])
hold on
for i=1:N
   if y(i) <= F(x(i))
       plot(x(i),y(i),'og')
   else
       plot(x(i),y(i),'or')
   end
end
F_monte = (length(z)/N) * (b-a) * (max(Y)-min(Y));

sprintf('Metoda Monte Carlo: %d', F_monte)