You can get an approximation using Monte Carlo Method :
a = matrix(c(0 ,0 ,2 ,0 ,2 ,2 ,0 , 2, 0, 0), byrow = T, ncol = 2)
b = matrix(c(.5, 0 ,1 , 1, 1.5, 0, .5, 0), ncol = 2, byrow = T)
library(sp)
interseccao <- function(a,b, n = 100000){
xmin <- min(a[,1], b[,1])
xmax <- max(a[,1], b[,1])
ymin <- min(a[,2], b[,2])
ymax <- max(a[,2], b[,2])
x_aleatorio <- runif(n, min = xmin, max = xmax)
y_aleatorio <- runif(n, min = ymin, max = ymax)
pontos_em_a <- point.in.polygon(x_aleatorio, y_aleatorio, pol.x = a[,1], pol.y = a[,2])
pontos_em_b <- point.in.polygon(x_aleatorio, y_aleatorio, pol.x = b[,1], pol.y = b[,2])
proporcao_intersec <- mean(pontos_em_a == 1 & pontos_em_b == 1)
area_intersec <- (xmax-xmin)*(ymax - ymin)*proporcao_intersec
return(area_intersec)
}
interseccao(a,b)
[1] 0.50368
In Monte Carlo Method, random points are generated in an area of which you already know the area. Then determine the proportion of points that occur in the area you want to calculate, in this case, the intersection of the polygons. This ratio multiplied by the area in which the points were generated, returns the desired area.
Of course there should be an algorithm that computes the exact area, but if an approximation is enough for you, you are here.