Algorithm fails to generate correct Delaunay triangulation

[JackStade]

I was using this library for a project of mine, and I came across a bug caused by the algorithm generating a triangulation that badly fails to meet the Delaunay condition.

Here is code that demonstrates the bug:

use delaunator::{triangulate, Point};

fn main {
    let points_raw = [
        [-1.0849334460551594, -1.0849334460551594], // 0
        [-1.5265847189170323, -1.5265847189170323], // 1
        [-2.095815826893781, -1.274850699584578], // 2
        [-2.770865690566727, -0.9763199041819535], // 3
        [-3.6843898126113546, -0.5723274031100705], // 4
        [-4.403658177176129, -0.25424163117441645], // 5
        [-5.02567003534954, 0.020833892521026076], // 6
        [-5.701084620214019, 0.3195259804640507], // 7
        [-6.561463270870451, 0.7000156846325452], // 8
        [-7.31105511135829, 1.0315115642859167], // 9
        [-8.0, 1.336187228503853], // 10
        [-7.339371017948894, 1.1048861305058357], // 11
        [-6.616986689211032, 0.8519630935920505], // 12
        [-5.816767071935726, 0.5717881676590966], // 13
        [-5.121128245254447, 0.3282293338915143], // 14
        [-4.512948676142796, 0.11529195763725286], // 15
        [-3.850960067633096, -0.11648517623155441], // 16
        [-3.1594534122386113, -0.3585972436874558], // 17
        [-2.288934450277422, -0.663385554827794], // 18
        [-1.6751076145244035, -0.8783001664294665], // 19
    ];
    let mut points = Vec::with_capacity(20);
    for [x, y] in &points_raw {
        points.push(Point { x: *x, y: *y });
    }
    let tris = triangulate(&points).unwrap();
    println!("{:?}", tris.triangles);
}

The output is:

[5, 6, 15, 6, 14, 15, 14, 16, 15, 15, 16, 5, 6, 7, 14, 16, 4, 5, 5, 4, 6, 7, 13, 14, 16, 17, 4, 14, 17, 16, 9, 8, 7, 7, 8, 13, 17, 3, 4, 4, 3, 6, 8, 12, 13, 19, 18, 17, 17, 18, 3, 6, 9, 7, 8, 9, 12, 18, 2, 3, 9, 11, 12, 12, 11, 13, 13, 11, 14, 14, 19, 17, 18, 19, 2, 19, 1, 2, 2, 10, 3, 10, 0, 3]

Which contains the triangle [10, 0, 3]. By plotting these points we can see that this triangle definitely shouldn't be part of the Delaunay triangulation, as the circumcircle of these points contains most of the other points:

It looks like there are other triangles in this triangulation that fail the condition. I have idea why this specific set of points causes a problem.

[timothee-haudebourg]

This seems to be quite an old issue... Has this been fixed yet? I believe this is one of the most downloaded delaunay triangulation crate (just after spade), so if nothing is planned to fix this bug, maybe users should be warned that this crate is bugged and no longer developed.

[mourner]

@timothee-haudebourg it's not a trivial issue to fix. You're welcome to submit a PR if you have time to spend on this — it would involve switching to Shewchuk's robust version of the orient2d check. Equivalent issue in JS (solved just recently) mapbox/delaunator#61

Also this kind of bug (caused by floating point errors) is present in a lot of computational geometry software, and usually has workarounds (such as scaling/rounding the input coordinates), so I wouldn't call it unusable because of it.

[timothee-haudebourg]

I see, so it is just a precision problem. Although it can be deduced by the precision of the floats values used in the example, that wasn't quite clear. Sorry if I sounded a bit harsh, I know it can be difficult to find time working on such projects, but it would not be the first time I come across a buggy dead crate without any such indication for the users. Your response gives just the information I needed, it would have been nice to give it earlier to help users navigate the sea of Rust crates. Apart from that, thanks for your work!

[mourner]

@timothee-haudebourg this crate definitely needs more love — I wrote it as a small Rust-learning project (porting my own JS library), but can't devote much free time to maintain it. New contributors are very welcome though.

[andreesteve]

@JackStade , @timothee-haudebourg - it would be super helpful if you could try this scenario against the branch from #19 - I was not able to repro the issue on current master on my box; but if this is related to the missing robustness checks, it should help. I am wondering if there is some other factor in play here (compiler version maybe?).

Here is the triangulation I am getting right now (does not include [10, 0, 3]):

EDIT: so I noticed that the triangulation being generated is putting some of the points that visually we would expect to be on the hull inside of the hull. For example, 19 is inside of the hull (black dot) and 0 is on the hull (red dot). Then degenerated triangle [19, 0, 14] is being formed.

[andreesteve]

So I was learning WASM and decided to build an inspector on the browser to compare the rust library running in WASM with the JS library.

This is this data set running against the master branch

And this is the data set for the robust branch

I still have no clue why when I compile rust to native, it cannot repro the issue on master. I can only imagine that the instructions that WASM compiles to may be slightly different than the ones for the native target. I wonder if the compiler, when it targets native, it is able to identify some float error accumulation patterns and replace with ones that accumulates less errors. Frankly I do not know enough about floats!

@mourner - this makes me feel better about #19

[mourner]

Closed by #19