axiolid-reference 0.3.15

//! In-plane regions of a section: faces lying in the cutting plane.
//!
//! The in-plane triangles may touch, abut and overlap. `arrangement` cuts
//! their source edges into pieces on exact points and says what covers each
//! side of a piece: `Below`, `Above`, `Both` or nothing. A region is a
//! maximal part of the plane of one class, so regions have pairwise disjoint
//! interiors, faces of one side merge whether they overlap or only share a
//! boundary segment, and `Both` is where faces of both sides coincide. The
//! pipeline, in order:
//!
//! 1. Boundary half-edges: a piece whose two sides differ is a half-edge of
//!    each non-empty class, directed with that class on its left.
//! 2. Fans: the pieces at a vertex are ordered counter-clockwise in the
//!    frame's sense by exact signs. The half-edge that follows an incoming
//!    one is the first boundary of its class met turning clockwise from the
//!    way it came, so a vertex where a region pinches is resolved per fan.
//! 3. Cycles: following half-edges closes rings. Each is an outer ring or a
//!    hole by the exact turn at its extreme vertex.
//! 4. Collinear vertices that only two boundary pieces meet at are dropped;
//!    the others, and the vertices of cut segments, stay.
//! 5. Each hole belongs to the innermost outer ring of its class that holds
//!    it, decided with exact point location.
//! 6. The rings are projected to the plane frame once per vertex, so the
//!    regions of all classes share bit-identical points, and made canonical.
//!
//! Only the returned coordinates are rounded. A ring that collapses under
//! that rounding is refused as `Degenerate`.

use std::collections::{BTreeMap, BTreeSet};

use axiolid_contracts::{ExecutionOptions, GeomError, GeomResult, Sign};
use axiolid_core::{Frame3, Point2, Point3};
use axiolid_mesh::TriMesh;
use axiolid_mesh_section_contract::{InPlaneRegion, InPlaneSide, RegionLimits};

use super::arrangement::{sort_exact, Class, InPlaneArrangement};
use super::plane_points::{Location, PlanePoints, RingBox};
use super::ring_nesting::TracedRing;
use super::{
    internal_topology_error, mesh_index, orient_sign, sequence_order, sign_times, ExactSectionPlane,
};
use crate::coplanar::{dominant_axis, project};

/// The axis dropped to reach two dimensions, and what that drop does to the
/// frame's own orientation.
#[derive(Debug, Clone, Copy)]
pub(super) struct Projection {
    axis: usize,
    frame_sign: Sign,
}

impl Projection {
    /// Pick the dominant axis of the plane normal, or the first axis that
    /// keeps the frame non-degenerate.
    ///
    /// Dropping an axis is an exact affine bijection on the plane when the
    /// projected frame points `(origin, origin + x, origin + y)` are not
    /// collinear; their exact orientation is `frame_sign`.
    pub(super) fn new(plane: &ExactSectionPlane) -> GeomResult<Self> {
        let dominant = dominant_axis(plane.normal);
        for axis in [dominant, 0, 1, 2] {
            let frame_sign = orient_sign(
                project(plane.origin, axis),
                project(plane.x_point, axis),
                project(plane.y_point, axis),
            )?;
            if frame_sign != Sign::Zero {
                return Ok(Self { axis, frame_sign });
            }
        }
        Err(GeomError::Degenerate(
            "section plane has no non-degenerate axis projection".into(),
        ))
    }

    fn flat(self, point: Point3) -> Point2 {
        project(point, self.axis)
    }

    /// The dropped axis: 0 = x, 1 = y, 2 = z.
    pub(super) fn axis(self) -> usize {
        self.axis
    }

    /// The exact orientation of the frame's own axes after the drop.
    pub(super) fn frame_sign(self) -> Sign {
        self.frame_sign
    }

    /// Exact orientation of three points of the plane in the frame's own
    /// `(x, y)` sense: positive when they turn counter-clockwise there.
    pub(super) fn orient(self, a: Point3, b: Point3, c: Point3) -> GeomResult<Sign> {
        let projected = orient_sign(self.flat(a), self.flat(b), self.flat(c))?;
        Ok(sign_times(projected, self.frame_sign))
    }
}

/// One source triangle lying in the plane.
#[derive(Debug, Clone, Copy)]
pub(super) struct InPlaneTriangle {
    /// Source vertices in winding order.
    pub(super) vertices: [u32; 3],
    /// Positive when the triangle is counter-clockwise in the frame, so its
    /// normal agrees with the frame's `+z`; negative otherwise.
    pub(super) facing: Sign,
}

impl InPlaneTriangle {
    pub(super) fn new(
        mesh: &TriMesh,
        projection: Projection,
        index: usize,
        vertices: [u32; 3],
    ) -> GeomResult<Self> {
        let [a, b, c] = corners(mesh, vertices)?;
        let facing = projection.orient(a, b, c)?;
        if facing == Sign::Zero {
            return Err(GeomError::Degenerate(format!(
                "in-plane source triangle {index} has zero area"
            )));
        }
        Ok(Self { vertices, facing })
    }
}

pub(super) fn corners(mesh: &TriMesh, vertices: [u32; 3]) -> GeomResult<[Point3; 3]> {
    let mut points = [Point3::ZERO; 3];
    for (slot, vertex) in points.iter_mut().zip(vertices) {
        *slot = mesh.positions[mesh_index(vertex)?];
    }
    Ok(points)
}

/// Whether `point` lies strictly between `a` and `b` on their common line.
pub(super) fn strictly_between(a: Point3, point: Point3, b: Point3) -> bool {
    let within = |a: f64, p: f64, b: f64| a.min(b) <= p && p <= a.max(b);
    point != a
        && point != b
        && within(a.x, point.x, b.x)
        && within(a.y, point.y, b.y)
        && within(a.z, point.z, b.z)
}

/// Index at which a ring starts to be canonical: at its lexicographically
/// smallest point; of several equal smallest points the rotation that orders
/// lowest wins.
pub(super) fn canonical_start(points: &[Point2]) -> usize {
    let rotated = |start: usize| -> Vec<Point2> {
        points[start..]
            .iter()
            .chain(&points[..start])
            .copied()
            .collect()
    };
    let smallest = points
        .iter()
        .min_by(|left, right| super::point_order(left, right))
        .copied();
    (0..points.len())
        .filter(|&start| Some(points[start]) == smallest)
        .min_by(|&left, &right| sequence_order(&rotated(left), &rotated(right)))
        .unwrap_or(0)
}

/// Rotate a ring to start at its lexicographically smallest point; of
/// several equal smallest points the rotation that orders lowest wins.
pub(super) fn canonical_rotation(points: Vec<Point2>) -> Vec<Point2> {
    let start = canonical_start(&points);
    points[start..]
        .iter()
        .chain(&points[..start])
        .copied()
        .collect()
}

/// The regions of a section and what the cut needs to know about them.
#[derive(Debug)]
pub(super) struct Regions {
    pub(super) regions: Vec<InPlaneRegion>,
    /// Source vertices that remain on a ring after collinear ones are dropped.
    pub(super) ring_vertices: BTreeSet<u32>,
    /// The rings of the `Both` regions on their exact points, aligned with
    /// the returned points.
    pub(super) both_rings: Vec<TracedRing>,
    /// Whether a ring of a `Both` region visits a vertex twice.
    pub(super) both_pinched: bool,
    /// Candidate pairs examined: the arrangement's and the hole containment's.
    pub(super) pairs: usize,
}

/// A direction out of a vertex along a piece.
#[derive(Debug, Clone, Copy)]
struct Spoke {
    to: usize,
    piece: usize,
    /// 0 for directions in `v > 0`, or `v = 0` and `u > 0`, in the
    /// projection's coordinates; 1 for the others.
    half: u8,
}

/// A closed boundary cycle of one class.
struct Cycle {
    class: Class,
    /// Kept vertex numbers.
    ring: Vec<usize>,
    outer: bool,
}

fn collapsed() -> GeomError {
    GeomError::Degenerate("in-plane region collapsed under projection".into())
}

/// What covers the left and the right of a piece, going out of `from`.
fn seen_from(arrangement: &InPlaneArrangement, piece: usize, from: usize) -> (Class, Class) {
    let piece = &arrangement.pieces[piece];
    if piece.ends[0] == from {
        (piece.left, piece.right)
    } else {
        (piece.right, piece.left)
    }
}

/// The pieces at each vertex that has a boundary, in counter-clockwise order
/// in the frame's sense. Both `orient` and the half-plane split are in the
/// projection's coordinates up to the frame sign, which reverses the order
/// within a half but keeps the cyclic order.
fn fans(
    arrangement: &InPlaneArrangement,
    exact: &PlanePoints,
    on_boundary: &[usize],
    options: &ExecutionOptions,
) -> GeomResult<BTreeMap<usize, Vec<Spoke>>> {
    let mut fans = BTreeMap::<usize, Vec<Spoke>>::new();
    for (index, piece) in arrangement.pieces.iter().enumerate() {
        for (from, to) in [
            (piece.ends[0], piece.ends[1]),
            (piece.ends[1], piece.ends[0]),
        ] {
            if on_boundary[from] == 0 {
                continue;
            }
            let (here, there) = (
                arrangement.vertices[from].exact,
                arrangement.vertices[to].exact,
            );
            let upper = exact.compare(there, here, 1)? == std::cmp::Ordering::Greater
                || (exact.compare(there, here, 1)? == std::cmp::Ordering::Equal
                    && exact.compare(there, here, 0)? == std::cmp::Ordering::Greater);
            fans.entry(from).or_default().push(Spoke {
                to,
                piece: index,
                half: u8::from(!upper),
            });
        }
    }
    for (&vertex, fan) in &mut fans {
        options.check_cancelled()?;
        let here = arrangement.vertices[vertex].exact;
        sort_exact(fan, |left, right| {
            if left.half != right.half {
                return Ok(left.half.cmp(&right.half));
            }
            let (a, b) = (
                arrangement.vertices[left.to].exact,
                arrangement.vertices[right.to].exact,
            );
            match exact.orient(here, a, b)? {
                Sign::Positive => Ok(std::cmp::Ordering::Less),
                Sign::Negative => Ok(std::cmp::Ordering::Greater),
                // Two pieces leave a vertex the same way only if they
                // overlap, which the arrangement has cut apart.
                _ => Err(internal_topology_error()),
            }
        })?;
    }
    Ok(fans)
}

/// The vertex the boundary of `class` leaves `at` for, after it came from
/// `from`: the first boundary met turning clockwise from `from`.
fn next_vertex(
    arrangement: &InPlaneArrangement,
    fans: &BTreeMap<usize, Vec<Spoke>>,
    from: usize,
    at: usize,
    class: Class,
) -> GeomResult<usize> {
    let fan = fans.get(&at).ok_or_else(internal_topology_error)?;
    let index = fan
        .iter()
        .position(|spoke| spoke.to == from)
        .ok_or_else(internal_topology_error)?;
    for step in 1..=fan.len() {
        let spoke = fan[(index + fan.len() - step) % fan.len()];
        let (left, right) = seen_from(arrangement, spoke.piece, at);
        // The sector between two neighbouring pieces is the left of one and
        // the right of the other.
        if left != class {
            return Err(internal_topology_error());
        }
        if right != class {
            return Ok(spoke.to);
        }
    }
    Err(internal_topology_error())
}

/// The eight lexicographic orders tried for the extreme vertex of a ring:
/// which axis is major, and whether each runs ascending.
const EXTREME_ORDERS: [(usize, bool, bool); 8] = [
    (0, true, true),
    (0, true, false),
    (0, false, true),
    (0, false, false),
    (1, true, true),
    (1, true, false),
    (1, false, true),
    (1, false, false),
];

/// Whether the cycle is an outer boundary (counter-clockwise in the frame),
/// by the exact turn at its extreme vertex.
///
/// The minimum of a lexicographic order is a convex corner. The first order
/// whose minimum the ring visits once decides; a ring that touches itself at
/// every extreme (a weakly simple ring) is refused.
fn is_outer(
    arrangement: &InPlaneArrangement,
    exact: &PlanePoints,
    ring: &[usize],
) -> GeomResult<bool> {
    for (major, major_up, minor_up) in EXTREME_ORDERS {
        let order = |left: usize, right: usize| -> GeomResult<std::cmp::Ordering> {
            let (a, b) = (
                arrangement.vertices[ring[left]].exact,
                arrangement.vertices[ring[right]].exact,
            );
            let first = exact.compare(a, b, major)?;
            let first = if major_up { first } else { first.reverse() };
            if first.is_ne() {
                return Ok(first);
            }
            let second = exact.compare(a, b, 1 - major)?;
            Ok(if minor_up { second } else { second.reverse() })
        };
        let mut smallest = 0;
        for index in 1..ring.len() {
            if order(index, smallest)?.is_lt() {
                smallest = index;
            }
        }
        let mut visits = 0;
        for index in 0..ring.len() {
            if order(index, smallest)?.is_eq() {
                visits += 1;
            }
        }
        if visits != 1 {
            continue;
        }
        let count = ring.len();
        let at = |index: usize| arrangement.vertices[ring[index % count]].exact;
        return match exact.orient(at(smallest + count - 1), at(smallest), at(smallest + 1))? {
            Sign::Positive => Ok(true),
            Sign::Negative => Ok(false),
            _ => Err(GeomError::Degenerate(
                "in-plane region boundary has a spike at its extreme vertex".into(),
            )),
        };
    }
    Err(GeomError::Degenerate(
        "in-plane region boundary visits every extreme vertex twice".into(),
    ))
}

/// What simplifying a ring may keep: vertices where other boundary meets, and
/// the protected ones.
struct Keep<'a> {
    /// Boundary pieces at each vertex.
    on_boundary: &'a [usize],
    protected: &'a BTreeSet<u32>,
}

/// Drop vertices that lie strictly between their neighbours on a line and
/// where only this ring's two boundary pieces meet, except protected ones; a
/// spike is refused.
fn simplify(
    arrangement: &InPlaneArrangement,
    exact: &PlanePoints,
    ring: &mut Vec<usize>,
    keep: &Keep<'_>,
    options: &ExecutionOptions,
) -> GeomResult<()> {
    let collapsed =
        || GeomError::Degenerate("in-plane region collapsed under simplification".into());
    let mut index = 0usize;
    let mut clean = 0usize;
    while clean < ring.len() {
        options.check_cancelled()?;
        if ring.len() < 3 {
            return Err(collapsed());
        }
        let at = index % ring.len();
        let (previous, vertex, next) = (
            ring[(at + ring.len() - 1) % ring.len()],
            ring[at],
            ring[(at + 1) % ring.len()],
        );
        let point = |vertex: usize| arrangement.vertices[vertex].exact;
        if exact.orient(point(previous), point(vertex), point(next))? != Sign::Zero {
            clean += 1;
            index = at + 1;
            continue;
        }
        if !exact.strictly_between(point(previous), point(vertex), point(next))? {
            return Err(GeomError::Degenerate(
                "in-plane region boundary has a zero-width spike".into(),
            ));
        }
        let protected = arrangement.vertices[vertex]
            .sources
            .iter()
            .any(|source| keep.protected.contains(source));
        if protected || keep.on_boundary[vertex] != 2 {
            clean += 1;
            index = at + 1;
            continue;
        }
        ring.remove(at);
        clean = 0;
        index = (at + ring.len().saturating_sub(1)) % ring.len().max(1);
    }
    if ring.len() < 3 {
        return Err(collapsed());
    }
    Ok(())
}

/// Whether `inner` lies in `outer`: its vertices not on the boundary agree.
/// `Degenerate` when they do not, or when none decides.
fn within(exact: &PlanePoints, inner: &[usize], outer: &[usize]) -> GeomResult<bool> {
    let (mut inside, mut outside) = (false, false);
    for &node in inner {
        match exact.locate(node, outer)? {
            Location::Inside => inside = true,
            Location::Outside => outside = true,
            Location::Boundary => {}
        }
    }
    match (inside, outside) {
        (true, false) => Ok(true),
        (false, true) => Ok(false),
        _ => Err(GeomError::Degenerate(
            "in-plane region rings cross or touch without a deciding vertex".into(),
        )),
    }
}

/// Count one candidate pair against the budget.
fn count_pair(pairs: &mut usize, max_pairs: usize) -> GeomResult<()> {
    *pairs += 1;
    if *pairs > max_pairs {
        return Err(GeomError::BudgetExceeded {
            resource: "section in-plane pairs",
        });
    }
    Ok(())
}

/// Regions of the arranged in-plane triangles: one per maximal part of a
/// class, with its holes.
///
/// `protected` vertices (nodes of cut segments) are never simplified away,
/// so a chain end stays a vertex of the ring it meets.
pub(super) fn build_regions(
    mesh: &TriMesh,
    frame: Frame3,
    arrangement: &InPlaneArrangement,
    exact: &PlanePoints,
    protected: &BTreeSet<u32>,
    limits: RegionLimits,
    options: &ExecutionOptions,
) -> GeomResult<Regions> {
    let mut pairs = arrangement.pairs;

    // Boundary half-edges, with the class on their left.
    let mut half_edges = BTreeSet::<(Class, usize, usize)>::new();
    let mut on_boundary = vec![0usize; arrangement.vertices.len()];
    for piece in &arrangement.pieces {
        if piece.left == piece.right {
            continue;
        }
        let [from, to] = piece.ends;
        on_boundary[from] += 1;
        on_boundary[to] += 1;
        if piece.left != Class::Empty {
            half_edges.insert((piece.left, from, to));
        }
        if piece.right != Class::Empty {
            half_edges.insert((piece.right, to, from));
        }
    }
    let fans = fans(arrangement, exact, &on_boundary, options)?;

    // Cycles, as vertex numbers.
    let mut seen = BTreeSet::<(Class, usize, usize)>::new();
    let mut cycles = Vec::<Cycle>::new();
    let keep = Keep {
        on_boundary: &on_boundary,
        protected,
    };
    let mut ring_vertices = BTreeSet::new();
    for &start in &half_edges {
        if seen.contains(&start) {
            continue;
        }
        let (class, _, _) = start;
        let mut ring = Vec::new();
        let mut edge = start;
        loop {
            options.check_cancelled()?;
            seen.insert(edge);
            ring.push(edge.1);
            let after = next_vertex(arrangement, &fans, edge.1, edge.2, class)?;
            edge = (class, edge.2, after);
            if edge == start {
                break;
            }
            if !half_edges.contains(&edge) || seen.contains(&edge) || ring.len() > half_edges.len()
            {
                return Err(internal_topology_error());
            }
        }
        let outer = is_outer(arrangement, exact, &ring)?;
        simplify(arrangement, exact, &mut ring, &keep, options)?;
        for &vertex in &ring {
            ring_vertices.extend(arrangement.vertices[vertex].sources.iter().copied());
        }
        cycles.push(Cycle { class, ring, outer });
    }

    // Exact rings, and the owner of each hole.
    let exact_ring = |cycle: &Cycle| -> Vec<usize> {
        cycle
            .ring
            .iter()
            .map(|&vertex| arrangement.vertices[vertex].exact)
            .collect()
    };
    let rings: Vec<Vec<usize>> = cycles.iter().map(exact_ring).collect();
    let boxes: Vec<RingBox> = rings.iter().map(|ring| exact.ring_box(ring)).collect();
    let mut owner = vec![None::<usize>; cycles.len()];
    for (hole, cycle) in cycles.iter().enumerate() {
        if cycle.outer {
            continue;
        }
        let mut containers = Vec::new();
        for (candidate, other) in cycles.iter().enumerate() {
            if !other.outer
                || other.class != cycle.class
                || !boxes[hole].possibly_within(&boxes[candidate])
            {
                continue;
            }
            options.check_cancelled()?;
            count_pair(&mut pairs, limits.max_in_plane_pairs)?;
            if within(exact, &rings[hole], &rings[candidate])? {
                containers.push(candidate);
            }
        }
        let mut innermost = Vec::new();
        for &candidate in &containers {
            let mut holds_all = true;
            for &other in &containers {
                if other == candidate {
                    continue;
                }
                count_pair(&mut pairs, limits.max_in_plane_pairs)?;
                if !within(exact, &rings[candidate], &rings[other])? {
                    holds_all = false;
                    break;
                }
            }
            if holds_all {
                innermost.push(candidate);
            }
        }
        match innermost.as_slice() {
            [only] => owner[hole] = Some(*only),
            _ => {
                return Err(GeomError::Degenerate(
                    "in-plane hole has no single innermost outer boundary".into(),
                ))
            }
        }
    }

    // Output points, rounded once per vertex.
    let mut points = BTreeMap::<usize, Point2>::new();
    let mut point_of = |vertex: usize| -> GeomResult<Point2> {
        if let Some(&point) = points.get(&vertex) {
            return Ok(point);
        }
        let point = arrangement.output_point(mesh, frame, exact, vertex)?;
        points.insert(vertex, point);
        Ok(point)
    };
    // A ring as plane-local points and the exact points they come from, in
    // canonical rotation.
    let mut traced = |cycle: &Cycle| -> GeomResult<TracedRing> {
        let local = cycle
            .ring
            .iter()
            .map(|&vertex| point_of(vertex))
            .collect::<GeomResult<Vec<_>>>()?;
        if local.len() < 3
            || (0..local.len()).any(|index| local[index] == local[(index + 1) % local.len()])
        {
            return Err(collapsed());
        }
        let start = canonical_start(&local);
        let rotate = |items: &[usize]| -> Vec<usize> {
            items[start..]
                .iter()
                .chain(&items[..start])
                .copied()
                .collect()
        };
        Ok(TracedRing {
            nodes: rotate(
                &cycle
                    .ring
                    .iter()
                    .map(|&vertex| arrangement.vertices[vertex].exact)
                    .collect::<Vec<_>>(),
            ),
            points: canonical_rotation(local),
        })
    };

    let mut regions = Vec::new();
    let mut both_rings = Vec::new();
    let mut both_pinched = false;
    let mut vertices = 0usize;
    for (index, cycle) in cycles.iter().enumerate() {
        if !cycle.outer {
            continue;
        }
        let Some(side) = cycle.class.side() else {
            return Err(internal_topology_error());
        };
        if regions.len() >= limits.max_regions {
            return Err(GeomError::BudgetExceeded {
                resource: "section in-plane regions",
            });
        }
        let outer = traced(cycle)?;
        let mut holes = Vec::new();
        for (hole, owned) in owner.iter().enumerate() {
            if *owned == Some(index) {
                holes.push(traced(&cycles[hole])?);
            }
        }
        holes.sort_by(|left, right| sequence_order(&left.points, &right.points));
        vertices += outer.points.len() + holes.iter().map(|ring| ring.points.len()).sum::<usize>();
        if vertices > limits.max_region_vertices {
            return Err(GeomError::BudgetExceeded {
                resource: "section in-plane region vertices",
            });
        }
        if side == InPlaneSide::Both {
            for ring in std::iter::once(&outer).chain(&holes) {
                let distinct: BTreeSet<usize> = ring.nodes.iter().copied().collect();
                both_pinched |= distinct.len() != ring.nodes.len();
            }
            both_rings.push(TracedRing {
                nodes: outer.nodes.clone(),
                points: outer.points.clone(),
            });
            for ring in &holes {
                both_rings.push(TracedRing {
                    nodes: ring.nodes.clone(),
                    points: ring.points.clone(),
                });
            }
        }
        regions.push(InPlaneRegion::new(
            side,
            outer.points,
            holes.into_iter().map(|ring| ring.points).collect(),
        ));
    }
    if owner
        .iter()
        .zip(&cycles)
        .any(|(owned, cycle)| !cycle.outer && owned.is_none())
    {
        return Err(GeomError::Degenerate(
            "in-plane component has a hole but no outer boundary".into(),
        ));
    }
    regions.sort_by(|left, right| {
        left.side
            .cmp(&right.side)
            .then_with(|| sequence_order(&left.outer, &right.outer))
    });
    Ok(Regions {
        regions,
        ring_vertices,
        both_rings,
        both_pinched,
        pairs,
    })
}

#[cfg(test)]
mod tests {
    use axiolid_core::{Tolerance, Vec3};

    use super::super::classify_vertices;
    use super::*;

    /// Eight triangles fan around `c = (1,1,0)` over `[0,2]^2`; the edge
    /// midpoints `(1,0)` and `(1,2)` are lifted to `z = 1`, so the four
    /// triangles touching them leave the plane. The remaining four form two
    /// groups that meet only at `c`: a pinch. Only the in-plane triangles
    /// are given; their missing neighbours make the boundary.
    fn built(protected: &BTreeSet<u32>) -> Regions {
        let positions = vec![
            Point3::new(1.0, 1.0, 0.0), // 0: c
            Point3::new(0.0, 0.0, 0.0), // 1: r0
            Point3::new(2.0, 0.0, 0.0), // 2: r2
            Point3::new(2.0, 1.0, 0.0), // 3: r3, collinear between r2 and r4
            Point3::new(2.0, 2.0, 0.0), // 4: r4
            Point3::new(0.0, 2.0, 0.0), // 5: r6
            Point3::new(0.0, 1.0, 0.0), // 6: r7, collinear between r6 and r0
        ];
        let mesh = TriMesh::new(positions, vec![]);
        let frame = Frame3 {
            origin: Point3::ZERO,
            x: Vec3::X,
            y: Vec3::Y,
            z: Vec3::Z,
        };
        let plane = ExactSectionPlane::new(frame).expect("plane");
        let projection = Projection::new(&plane).expect("projection");
        let classifications = classify_vertices(&mesh, plane).expect("classes");
        let triangles: Vec<InPlaneTriangle> = [[0, 2, 3], [0, 3, 4], [0, 5, 6], [0, 6, 1]]
            .into_iter()
            .enumerate()
            .map(|(index, vertices)| {
                InPlaneTriangle::new(&mesh, projection, index, vertices).expect("triangle")
            })
            .collect();
        let options = ExecutionOptions::new(Tolerance::METRE);
        let mut exact = PlanePoints::new(&plane, projection).expect("points");
        let arrangement = InPlaneArrangement::build(
            &mesh,
            &classifications,
            &triangles,
            &mut exact,
            1000,
            &options,
        )
        .expect("arrangement");
        build_regions(
            &mesh,
            frame,
            &arrangement,
            &exact,
            protected,
            RegionLimits::new(10, 100, 100),
            &options,
        )
        .expect("regions")
    }

    fn regions(protected: &BTreeSet<u32>) -> Vec<InPlaneRegion> {
        built(protected).regions
    }

    fn ring(points: &[(f64, f64)]) -> Vec<Point2> {
        points.iter().map(|&(x, y)| Point2::new(x, y)).collect()
    }

    #[test]
    fn a_pinch_gives_two_regions_and_drops_collinear_vertices() {
        // The boundary of the left group is c, r6, r7, r0 = (1,1), (0,2),
        // (0,1), (0,0); (0,1) lies between its neighbours and goes. The ring
        // is counter-clockwise from its smallest point (0,0). The right group
        // is (1,1), (2,0), (2,1), (2,2) without (2,1). Sorted by outer ring.
        let found = regions(&BTreeSet::new());
        assert_eq!(found.len(), 2);
        assert_eq!(found[0].side, InPlaneSide::Below);
        assert_eq!(found[0].outer, ring(&[(0.0, 0.0), (1.0, 1.0), (0.0, 2.0)]));
        assert_eq!(found[1].outer, ring(&[(1.0, 1.0), (2.0, 0.0), (2.0, 2.0)]));
        assert!(found.iter().all(|region| region.holes.is_empty()));
    }

    #[test]
    fn the_kept_ring_vertices_exclude_the_dropped_collinear_one() {
        // Vertex 6 is (0,1), the collinear one on the left group.
        assert!(!built(&BTreeSet::new()).ring_vertices.contains(&6));
        assert!(built(&BTreeSet::from([6])).ring_vertices.contains(&6));
    }

    #[test]
    fn a_protected_collinear_vertex_is_kept() {
        let found = regions(&BTreeSet::from([6]));
        assert_eq!(
            found[0].outer,
            ring(&[(0.0, 0.0), (1.0, 1.0), (0.0, 2.0), (0.0, 1.0)])
        );
        assert_eq!(found[1].outer, ring(&[(1.0, 1.0), (2.0, 0.0), (2.0, 2.0)]));
    }
}