Implementation of Simulation of Simplicity by Edelsbrunner and Mücke

Simulation of simplicity is a technique for ignoring degeneracies when calculating geometric predicates, such as the orientation of one point with respect to a list of points. Each point p_ i is perturbed by some polynomial in ε, a sufficiently small positive number. Specifically, coordinate p_(i,j) is perturbed by ε^(3^(d*i - j)), where d is more than the number of dimensions.

Predicates

Orientation

The orientation of 2 points p_0, p_1 in 1-dimensional space is positive if p_0 is to the right of p_1 and negative otherwise. We don't consider the case where p_0 = p_1 because of the perturbations.

The orientation of n points p_0, ..., p_(n - 1) in (n - 1)-dimensional space is the same as the orientation of p_1, ..., p_(n - 1) when looked at from p_0. In particular, the orientation of 3 points in 2-dimensional space is positive iff they form a left turn.

Orientation predicates for 1, 2, and 3 dimensions are implemented. They return whether the orientation is positive.

In Hypersphere

The in-circle of 4 points measures whether the last point is inside the circle that goes through the first 3 points. Those 3 points are not collinear because of the perturbations.

The in-sphere of 5 points measures whether the last point is inside the sphere that goes through the first 4 points. Those 4 points are not coplanar because of the perturbations.

Usage

use simplicity::{nalgebra, orient_2d};
use nalgebra::Vector2;

let points = vec![
    Vector2::new(0.0, 0.0),
    Vector2::new(1.0, 0.0),
    Vector2::new(1.0, 1.0),
    Vector2::new(0.0, 1.0),
    Vector2::new(2.0, 0.0),
];

// Positive orientation
let result = orient_2d(&points, |l, i| l[i], 0, 1, 2);
assert!(result);

// Negative orientation
let result = orient_2d(&points, |l, i| l[i], 0, 3, 2);
assert!(!result);

// Degenerate orientation, tie broken by perturbance
let result = orient_2d(&points, |l, i| l[i], 0, 1, 4);
assert!(result);
let result = orient_2d(&points, |l, i| l[i], 4, 1, 0);
assert!(!result);

Because the predicates take an indexing function, this can be used for arbitrary lists without having to implement Index for them:

let points = vec![
    (Vector2::new(0.0, 0.0), 0.8),
    (Vector2::new(1.0, 0.0), 0.4),
    (Vector2::new(2.0, 0.0), 0.6),
];

let result = orient_2d(&points, |l, i| l[i].0, 0, 1, 2);