%==========================================================================
%SIMULACIÖN DE MONTE CARLO: Estimación del área de una circunferencia de 
%radio unitario (A=pi*R^2=pi*1=pi)
%Para estimar el área, se coloca circunferencia dentro de un recinto de
%área conocida, en este caso un cuadrado de lado 2 (Área=4). El centro de 
%la circunferencia coincide con el cero de los ejes cartesianos. Se genera 
%un conjunto de números aleatorios tal que -1<x<1 y -1<y<1. Los mismos 
%conforman N pares ordenados. Luego se cuantifica la cantidad de puntos que 
%caen dentro de la circunferencia (ND), los cuáles deben cumplir la 
%condición x^2+y^2<1. Finalmente se estima el área buscada de la siguiente 
%forma:
%                                 Aest=4*ND/N
%
%donde N es la cantidad total de puntos generados
%==========================================================================
close all; clear; clc
%Variables aleatorias X e Y entre -1 y 1
N=500;
x=2*rand(1,N)-1;
y=2*rand(1,N)-1;
%Bucle de simulación con k = Cantidad de puntos considerados 
AEst=zeros(1,N);
Error=zeros(1,N);
ATeo=pi*(ones(1,N));
n=1:N;
%Se estima E(x) considerando cada ves más muestras
for k=1:N;
    %Cantidad de puntos dentro
    ND=sum((x(1:k).^2+y(1:k).^2)<1);
    %Cantidad de puntos fuera
    NF=k-ND;
    %Estimación de la esperanza
    AEst(k)=4*ND/k;
    Error(k)=abs(pi-AEst(k));
    %Visualización
    subplot(311),bar([0,1],[NF,ND]);
    legend(['Valor esperado con ' num2str(k) ' intentos = ' num2str(AEst(k))]);
    axis([0 1 0 N]);
    subplot(312),plot(n,AEst,'.-',n,ATeo,'r--');
    axis([1 N 0 4]);
    subplot(313),plot(n,Error,'.-');
    pause (0.01);
end