MathGroup Archive 2004

[Date Index] [Thread Index] [Author Index]

Search the Archive

Re: Voronoi Volume calculation

  • To: mathgroup at
  • Subject: [mg50590] Re: Voronoi Volume calculation
  • From: "Steve Luttrell" <steve_usenet at>
  • Date: Sat, 11 Sep 2004 06:44:51 -0400 (EDT)
  • References: <chp91s$jln$>
  • Sender: owner-wri-mathgroup at

Here is a possible way of doing it based on using the Boole function in the 
Calculus`Integration` package. It computes the volume of a Voronoi cell very 
quickly in 2 dimensions, in a few tens of seconds in 3 dimensions, but took 
too long in 4 dimensions for me to get any answer out (I use a 2GHz Athlon 
PC running Windows XP). I presume the method will work if you are patient.

Input the package.

<< Calculus`Integration`

Define a Eucliean distance function.

distsq[x1_, x2_] := #.# &[x1 - x2];

Generate a set of points in the region [-1,1] X [-1,1]. Here I generate 
rationals rather than reals.

denom = 10;
numpoints = 10;
points = Table[Random[Integer, {-denom, denom}]/denom, {numpoints}, {2}];

Use the origin as the reference point whose Voronoi cell's volume you want 
to compute.

Compute the distance of an arbitrary point {x,y} from the origin.

d0sq = distsq[{x, y}, {0, 0}];

Compute a list of inequality constraints that ensure that {x,y} lies closer 
to the origin than to any of the above points. If you have the Delaunay 
triangulation available then (for efficiency) you can use only those points 
whose Voronoi cells share a common side with the cell at the origin; 
Mathematica can compute a 2D Delaunay triangulation but I can't find a 
higher dimensional version, so I have NOT made use of Delaunay triangulation 
information here.
boundslist = Map[Simplify[distsq[{x, y}, #] - d0sq] >= 0 &, points];

Convert this list of constraints into a single Boolean constraint.

constraint = Apply[And, boundslist];

Finally compute the 2-dimensional Voronoi cell volume. I have assumed that 
the Voronoi cell around the origin as computed using the above list of 
points is contained within the region [-1,1] X [-1,1]; and if this condition 
happens to be false (e.g. not enough points or the points are arranged in an 
unlucky configuration) then the Voronoi cell is truncated by its 
intersection with the boundary of [-1,1] X [-1,1].

Integrate[Boole[constraint], {x, -1, 1}, {y, -1, 1}]

That deals with the 2-dimensional case, and evaluates the result very 

Now gather all the code together in the 3-dimensional case to obtain the 
following code:

<< Calculus`Integration`

distsq[x1_, x2_] := #.# &[x1 - x2];

denom = 10;
numpoints = 10;
points = Table[Random[Integer, {-denom, denom}]/denom, {numpoints}, {3}];
d0sq = distsq[{x, y, z}, {0, 0, 0}];
boundslist=Map[Simplify[distsq[{x,y,z},#]-d0sq] >= 0 &,points];
constraint = Apply[And, boundslist];
Integrate[Boole[constraint], {x, -1, 1}, {y, -1, 1}, {z, -1, 1}]

The result is a rational number with an impressively large number of digits 
in its numerator and denominator.

The 4-dimensional case is analogous, but needs more time to run than my 
patience allowed.

Good kuck!

Steve Luttrell

"Ran" <efeldesh at> wrote in message 
news:chp91s$jln$1 at
> Hi everyone,
> I was trying to find the "volume" of a voronoi shape for a set of
> points in a 4-D space.
> what is the best way to do it in Mathematica?
> couldn't find any remarks on volume of these shapes at the wolfarm
> site
> thanks ahead,
> - Ran :)

  • Prev by Date: Trying to eliminate a loop
  • Next by Date: Re: Re: Log[4]==2*Log[2]
  • Previous by thread: Re: Voronoi Volume calculation
  • Next by thread: NestList