Help!
I've got a wacky idea I'd like to explore (pseudo-preliminary-pre-research stage) which will require some code, namely for Voronoi diagrams. What I'd like seems like it should be available. I think I want the following:
1) I input a list of points -- 2-D is fine, but hey, if the code can handle 3-D, more time-wasting fun.
1a) Update: I suppose what I really want is for code that works on the "2-D torus" -- i.e say on the unit square [0,1]^2 where the boundaries wrap around, so that there's symmetry. That would be ideal, but from what I understand "torus" versions of the algorithms are harder to find, so maybe I'll just have to deal with regular 2-D and truncate somehow.
2) For the output, all I really want is the area of each cell corresponding to each point. Sure, if more information is available -- things like number of sides of the cell, etc. -- again, more time-wasting fun for me. But for now cell area seems most interesting. (I don't need pretty pictures.)
3) Now, the painful part -- I expect to be doing a lot of incremental inserting and deleting of points. So I'd like to have code that's very fast if say I delete a point from my list and insert a new point elsewhere. The code I've seen available seems to (if I'm understanding right) just implement the algorithms where you have a list of points. Some of them use incremental insertion and can therefore probably be modified easily to handle point insertions, but not point deletions. Somehow I think dealing with point sets of thousands or more and re-computing from scratch when I delete a point seems like it will be way too slow for my explorations. I realize handling deletions is harder, but that there exist algorithms for it, so I was surprised I can't easily find relevant code. But maybe I'm just looking at the wrong place.
I suppose whatever language the code is I can figure out how to use/wrap/rewrite it to my needs, so I shouldn't quibble about those sorts of issues.
Strangely, I can't seem to find anything that meets these needs. Feel free to set me straight via comments or mail. Or, if it seems to be something I'll need to implement myself, feel free to point me to the most helpful relevant papers.
The first libraries that come to mind are CGAL (cgal.org) and Qhull (qhull.org).
ReplyDeleteAFAIK neither supports dynamic addition or deletion of points. But Voronoi diagrams are O(n log n), so if you have only a few thousand points and a multi-gigahertz machine, recomputing everything from scratch may be plenty fast.
http://cran.r-project.org/web/packages/deldir/index.html
ReplyDeleteI used it once three years ago. He won't get the volumes
of each cell though (but for each tile he can bound it by the
volume of the sphere passing by the d+1 edges further
from the center)
tile.list
tile.centroids
deldir
try C++ numerical recipes 3rd edition page 1145
ReplyDeleteI'll check, but I don't think any of these things support dynamic deletion of points. (I don't think the C++ numerical recipes does either, though I'm finding it a bit hard to tell -- may need a more careful read.) And Kevin, I intend (as usual) to be thrashing the code, so even at few thousand points, I think the difference will be something like 3 orders of magnitude in time, so that's the difference between finishing experiments in a day or in a year...
ReplyDeleteI forgot one other desiderata -- I'd really like for the code to do this for points "on the unit torus", just because I like symmetry, and to make the problem even harder. Will add to the original post.
How about http://math.lbl.gov/voro++/about.html?
ReplyDeleteThey seem to support tori, though I am not sure about its dynamic behavior. Beyond: It maybe problematic if you're aiming at exact constructions.
Anon:
ReplyDeleteI had a student look at this and he had been working with Voro++. It doesn't handle deletions on its own but he suggested the code could be modified to handle it. Though, truthfully, it seemed a bit slow -- I might be able to handle experiments with 1000 points, but not 10000.
But that might be the way to go.
Isn't dynamic insertion and deletion necessarily Omega(n) at least? Consider n-1 points on a circle, and repeatedly inserting and deleting the center of the circle. Each insertion changes all n-1 Voronoi cells.
ReplyDeleteEric --
ReplyDeleteMy cases won't be worst cases, so I'd like something that has good "average case" performance (like randomized incremental algorithms should on insertions, but also for deletions).
Similarly, vor++ seems maybe a little loose on how it handles degenerate cases -- but I don't really care about degenerate cases.
What are your accuracy requirements? I suspect one could come with some approximate version based on kd trees based nearest neighbor data structures.
ReplyDeleteIt would be easier if you worked with Delaunay triangulations instead of Voronoi diagrams, particularly if all you want is the area of the Voronoi region, which can be extracted quite readily from the DT. Deleting a point from a 2-d triangulation is quite simple if you have access to primitives for 2-2 flips and 3-1 flips. Try Shewchuk's code.
ReplyDeleteYep. Try any reasonable implementation of Delaunay triangulation. Any such implementation would have insertions and deletions. CGAL is a good starting point. I would just replicate the points 9 times if you want tori topology. Sariel.
ReplyDeleteThanks -- I was wondering if Delaunay triangulations might prove easier to deal with, but wasn't clear on how hard getting cell areas and things might be. I'll take a look.
ReplyDeleteTo add to that: computing the area of a voronoi region is a primitive in "Natural Neighbor" interpolation. The code for that is given here:
ReplyDeletehttp://code.google.com/p/nn-c/
Someone created a dynamic voroni in d3.js. You can add, remove, and drag nodes around the screen: http://bl.ocks.org/1734663
ReplyDeleteTo compute Voronoi diagrams on the torus, build a 3x3 grid of squares, dump a copy of your point set in each square, compute the Voronoi diagram of the 9n points, and then clip to the center square.
ReplyDeleteThis comment has been removed by the author.
ReplyDeleteI was wondering if Delaunay triangulations might prove easier to deal with
ReplyDeleteNope. Any data structure that represents a Voronoi diagram can be interpreted as a data structure that represents a Delaunay triangulation and vice versa. Thus, any algorithm that computes Voronoi diagrams also computes Delaunay triangulations and vice versa. Duality is your friend.
Jeff --
ReplyDeleteBoth your comments are good/useful in theory. In practice, however, making 9 copies of the points is going to be a 10x slowdown; always nice to avoid that if possible (if there's an easy way to do so). Similarly, just because the two are duals doesn't mean that one isn't easier to work with than another. (Or, in particular, the point is I don't want to "work with" either; I just want something that gives me areas of Voronoi cells out of the box.)
In my experience, MATLAB is the path of least resistance for quickly (inefficiently) testing out ideas. Apparently there are ways to add/delete points without rebuilding everything but when I did these kinds of tests, I just rebuilt everytime since I was going to implement things correctly later and I am not sure how much functionality octave supports.
ReplyDeletehttp://www.mathworks.com/products/matlab/examples.html?file=/products/demos/shipping/matlab/demoDelaunayTri.html#13
CGAL has a structure for performing dynamic insertion/deletion that is fast most of the time. But just getting started with CGAL takes some effort...
http://www.cgal.org/Manual/latest/doc_html/cgal_manual/Triangulation_2/Chapter_main.html#Section_37.10
About reflecting points: can't you compute the Voronoi diagram and then only reflect points with Voronoi cells that touch the boundary and add them to the triangulation. If you are already writing/finding the code to clip Voronoi cells to the boundary (and incrementally adding points), that should be not too bad. And "usually" you only have to reflect O(n^{1/2}) points (in 2D).
This might be useful!! nice way of visualizing an insertion.. a deletion is highly probably as easy!
ReplyDeletehttp://bl.ocks.org/4060366