R - Ellipse Area with Montecarlo Method

Question

I need to calculate the area of the eclipse (a=6 b=3) with the Montecarlo Method. Also I have to make a plot (a diagram) of the result with the inside points red and the out ones black. At the end I have to compare the "Montecarlo result" with the "Regular Result"

The equation is (x^2)/36+(y^2)/9=1

The method must have 100000 replies.

This is what I do. Obviously it doesn't work.

set.seed(157619)

n <- 100000

xmin <- (-6)
xmax <- (+6)

ymin <- (-3)
ymax <- (+3)

rx <- (xmax-xmin)/2
ry <- (ymax-ymin)/2


outa <- runif(n,min=xmin,max=xmax)
outb <- runif(n,min=ymin,max=ymax)

dx <- outa*2
dy <- outb*2

ly <- dy<=(ry^2); my <- dy>(ry^2)
lx <- dx<=(ry^2); mx <- dx>(rx^2)

This is an example code that work for the circle:

n <- 200
xmin <- -1; xmax <- 1
r <- (xmax-xmin)/2
out <- runif(n,min=xmin,max=xmax)
x <- matrix(out,ncol=2)
d <- x[,1]^2 + x[,2]^2
l <- d<=(r^2); m <- d>(r^2)
win.graph(7,7.8) # così è quadrato
plot(c(xmin,xmax),c(xmin,xmax),type="n")
plot(x[l,1],x[l,2])
points(x[m,1],x[m,2],col="red",pch=19)
(p <- sum(l)/length(l))
p*4

Answer 1

I suspect this is homework, but here we go:

set.seed(42)
n <- 1e5
xmax <- 6
ymax <- 3

x <- runif(n, 0, xmax)
y <- runif(n, 0, ymax)

inside <- (x^2)/36+(y^2)/9 <= 1

plot(x, y, pch=16, cex=0.5, col=inside+1)

mean(inside) * (xmax*ymax) *4
#[1] 56.54376
pi*6*3
#[1] 56.54867

Answer 2

set.seed(1)
n = 1000
a = 6
b = 3
x.samp = runif(n, -a, a)
y.samp = runif(n, -b, b)

p.in = (x.samp/a)^2 + (y.samp/b)^2 <= 1

S = 4*a*b*sum(p.in)/n
print(S)

plot(x.samp, y.samp, col = p.in + 1)