The area of a polygon (or any plane shape) can be evaluated by Monte Carlo integration. The process involves 'random' sampling of (x,y) points in the xy plane, and testing whether they fall inside or outside the shape. That is illustrated schematically below for a circle.
Inputs to your function will be N, the number of points to sample, and polygonX & polygonY, which together constitute an ordered list of polygon vertices. N will always be a positive integer. polygonX & polygonY will always describe a single, continuous outline; however, the polygon may be self-intersecting.
For polygons that are self-intersecting, you must find the area of the enclosed point set. In the case of the cross-quadrilateral, it is treated as two simple triangles (cf. three triangles in the first figure below), and the clockwise/counterclockwise ordering of points is irrelevant. A self-intersecting pentagram (left & middle of the second figure below) will therefore have the same area as a concave decagon (right).
% Input N = 1000 polygonX = [1 0 -1 0] polygonY = [0 1 0 -1] % Output A = 2.036
Explanation: the above polygon is a square of side-length √2 centred on the origin, but rotated by 45°. The exact area is hence 2. As the value of N is moderate, the estimate by Monte Carlo integration is moderately accurate (an error of 0.036 in this example, corresponding to 509 of the 1000 sampled points falling within the polygon). Of course, if the code were run again then a slightly different value of A would probably be returned, such as A=1.992 (corresponding to 498 of the 1000 sampled points falling within the polygon), due to the random qualities of the technique.