Playing around with the vd-code has highlighted some problems. Using double-precision floating point numbers for coordinates, strange things start to appear when the coordinates are around 1E-6 or so. It's not possible to accurately compute the sign of the in-circle predicate (a 4x4 determinant) and the diagram update should instead be based on topological reasoning.
The dual of the Voronoi diagram is the Delaunay triangulation. Here I've modified my earlier VD code to also output the dual (shown in red). If this can be developed further to do constrained Delaunay triangulations (arbitrary pockets with islands) then it will be useful for cutter-location surfaces in opencamlib.
The algorithm works by incrementally constructing the diagram while adding the point-generators one by one. This initial configuration is used at the beginning:
The diagram will be correct if new generator points are placed inside the orange circle (this way we avoid edges that extend to infinity). Once about 500 randomly chosen points are added the diagram looks like this:
Because of floating-point errors it gets harder to construct the correct diagram when more and more points are added. My code now crashes at slightly more than 500 generators (due to an assert() that fails, not a segfault - yay!). It boils down to calculating the sign of a 4-by-4 determinant, which due to floating-point error goes wrong at some point. That's why the Sugihara-Iri algorithm is based on strictly preserving the correct topology of the diagram, and in fact they show in their papers that their algorithm doesn't crash when a random error is added to the determinant-calculations, and it even works (in some sense) when the determinant calculation is completely replaced by a random number generator. Their 1992 paper constructs the diagram for one million generators using 32-bit float
numbers, while my naive attempt using double here crashes after 535 points... some work to do still then!
Note how within one voronoi-region the offset is determined by the associated generator (a line, point, or arc). It's very easy to figure out the offset within one region: points and arcs have circular offsets, while lines have line offsets. The complete offset path is a result of connecting the offset from one region with the next region.