# This script computes Pi using Monte Carlo method. For this exercise we assume to draw random points inside the square on the [-1,1] unit (We consider a circle of radius R equal to 1 inscribed in a square with side length 2).

R <- 1
pi.estimate<-function(N){
x <- runif(N, min= -R, max= R)             # Generate two vectors for random points in unit circle
y <- runif(N, min= -R, max= R)
is.inside <- (x^2 + y^2) <= R^2             # Test if draws are inside the unit circle
pi.estimate <- 4 * sum(is.inside) / N      # Estimate Pi
pi.estimate
}
pi.estimate(20)              #examples
pi.estimate(200)
pi.estimate(2000)
pi.estimate(2000000)

# We can also graphically represent our simulation

N <- 20000
R <- 1
x <- runif(N, min= -R, max= R)
y <- runif(N, min= -R, max= R)
is.inside <- (x^2 + y^2) <= R^2
pi.estimate <- 4 * sum(is.inside) / N
print(pi.estimate)

plot.new()
plot.window(xlim = 1.1 * R * c(-1, 1), ylim = 1.1 * R * c(-1, 1))
points(x[ is.inside], y[ is.inside], pch = '.', col = "green")
points(x[!is.inside], y[!is.inside], pch = '.', col = "blue")