If I were doing this:
I would try to find a way to make an N-dimensional array (where N is the dimesnion of the region of integration). Call it pop, for population. The array elements ill be (hyper)cubes of equal volume, and together they will cover the region.
The initial value of the array elements that are outside the region or which straddle the border is is -1. THe initial value of all cubes that are totally inside the region will be zero. This is a weak point of my idea, because it may be difficult to efficiently determine which cubes are totally inside the region. (One approach: check each of the 2^N corner points of each cube.. If they are all inside then the cube is totally inside. That could perhaps occasionally fail for concave regions, but probably unlikely.)
As you do the integration by choosing random points, add 1 to the appropriate element of pop, each time you choose a point. When you are done, you reshape pop into a 1D vector and you throw away all elements whose vaue is -1. This values in the cells of this vector should be uniformly distributed, if your MOnte Calro sampled the space randomly and evenly. You can use the Kolmogorov-Smirnov test or the chi-squared test to test if the values are consistent with the uniform distribution.
Good luck!
