Add naive O(n m) point location for int and float triangulation

[mikwielgus]

n is the number of triangles in triangulation, m is the number of points to locate.

Reference: #6

[NailxSharipov]

The overall direction looks good. I would like to adjust the public api.

My main request is that we keep the trait-based API discussed above. Even if the first implementation is simple and there is not much shared code today.

Something like:

pub trait PointInTriangulationLocator<P> {
    fn locate_points<T>(&self, points: &[P]) -> Vec<PointLocationInTriangulation>
    where
        P: FloatPointCompatible<T>,
        T: FloatNumber;
}

pub trait IntPointInTriangulationLocator {
    fn locate_points(&self, points: &[IntPoint]) -> Vec<PointLocationInTriangulation>;
}

[mikwielgus]

My main request is that we keep the trait-based API discussed above. Even if the first implementation is simple and there is not much shared code today.

Ok, I have added the traits and implemented them by calling the currently existing methods, so that the downstream programmer still doesn't have to import the traits to use the methods if there is no need for shared behavior. (I prefer this style myself. Not sure if you do)

[fogti]

Would it make sense to mark these forwarding methods as #[inline]?

[fogti]

Okay, in the context of this, I'm getting annoying that i_float's FloatPointCompatible is generic over T instead of having that as an associated type, even tho FloatPointCompatible in essence defines an isomorphism, because then it doesn't make sense to have FloatPointCompatible implemented for multiple different T...

[fogti]

I found out that the sorting part is fairly easy, but implementing the sliding-window / lockstep iteration is a bit more complicated. Is there some reference / other place where something like this was already implemented, such that I can use that as a basis?

[fogti]

I have something that looks like this:

diff --git i/iTriangle/src/float/locator.rs w/iTriangle/src/float/locator.rs
index b39c745..3f6ce75 100644
--- i/iTriangle/src/float/locator.rs
+++ w/iTriangle/src/float/locator.rs
@@ -1,4 +1,4 @@
-use alloc::vec::Vec;
+use alloc::{btreeset::BTreeSet, vec, vec::Vec};
 use i_overlay::i_float::fix_vec::FixVec;
 use i_overlay::i_float::float::compatible::FloatPointCompatible;
 use i_overlay::i_float::float::number::FloatNumber;
@@ -7,7 +7,7 @@ use i_overlay::i_float::triangle::Triangle;
 use crate::{
     float::triangulation::Triangulation,
     int::triangulation::IndexType,
-    location::{PointLocationInTriangulation, TriangleIndex},
+    location::{LocationInfo, PointLocationInTriangulation, TriangleIndex},
 };
 
 pub trait PointInTriangulationLocator<P> {
@@ -22,7 +22,50 @@ impl<P, I: IndexType> Triangulation<P, I> {
     where
         P: FloatPointCompatible<T>,
     {
-        let mut result = Vec::with_capacity(points.len());
+        let mut tmp: Vec<LocationInfo> = points
+            .iter()
+            .enumerate()
+            .map(|(ptidx, _)| LocationInfo::Point(ptidx))
+            .chain(self.indices.iter().enumerate().map(|(idxidx, ptidx)| {
+                LocationInfo::TrianglePoint {
+                    point: ptidx,
+                    index: idxidx,
+                }
+            }))
+            .collect();
+        tmp.sort_by_key(|i| {
+            (
+                to_fix_vec::<P, T>(match i {
+                    LocationInfo::Point(ptidx) => &points[ptidx],
+                    LocationInfo::TrianglePoint { point, .. } => &self.points[point.into_usize()],
+                }),
+                *i,
+            )
+        });
+
+        let mut result = vec![PointLocationInTriangulation::Outside; points.len()];
+        let mut active_triangles = vec![0u8; self.indices.len() / 3];
+        let mut active_triangles2 = BTreeSet::new();
+
+        for locinf in tmp {
+            match locinf {
+                LocationInfo::Point(ptidx) => {
+                    for &active in active_triangles2.iter() {
+                        
+                    }
+                },
+                LocationInfo::TrianglePoint { index, .. } => {
+                    let index = index / 3;
+                    let cur = &mut active_triangles[index];
+                    *cur = (*cur + 1) % 3;
+                    match *cur {
+                        0 => active_triangles2.remove(&index);
+                        1 => active_triangles2.insert(index);
+                        _ => {}
+                    }
+                }
+            }
+        }
 
         for point in points.iter() {
             let mut vertex_hits = Vec::new();
diff --git i/iTriangle/src/location.rs w/iTriangle/src/location.rs
index 4f5e2e2..19e4c87 100644
--- i/iTriangle/src/location.rs
+++ w/iTriangle/src/location.rs
@@ -24,3 +24,9 @@ pub enum PointLocationInTriangulation {
     OnEdge(TriangleIndex, TriangleIndex),
     OnVertex(Vec<TriangleIndex>),
 }
+
+#[derive(Clone, Copy, PartialEq, Eq, PartialOrd, Ord)]
+pub(crate) enum LocationInfo<I> {
+    Point(usize),
+    TrianglePoint { point: I, index: usize },
+}

But I think I need to search forwards and backwards to find all involved triangles (because that code only would work reliably for being in the triangle, but not for the case where a point is in the same place as a triangle end-point).

[NailxSharipov]

Okay, in the context of this, I'm getting annoying that i_float's FloatPointCompatible is generic over T instead of having that as an associated type, even tho FloatPointCompatible in essence defines an isomorphism, because then it doesn't make sense to have FloatPointCompatible implemented for multiple different T...

Thanks again for your PR in i_float. I've already migrated all depends(i_shape, i_overlay and i_triangle).

[NailxSharipov]

I found out that the sorting part is fairly easy, but implementing the sliding-window / lockstep iteration is a bit more complicated. Is there some reference / other place where something like this was already implemented, such that I can use that as a basis?

The closest reference is ShapeBinder in iOverlay:

https://github.com/iShape-Rust/iOverlay/blob/main/iOverlay/src/bind/solver.rs

More specifically, look at:

fn private_solve<S: KeyExpCollection<VSegment, i32, ContourIndex>>(
    mut scan_list: S,
    shape_count: usize,
    anchors: Vec<IdSegment>,
    segments: Vec<IdSegment>,
) -> BindSolution

The problem there is slightly different, but the sweep logic is very similar.

In iOverlay we have inner contours / holes and outer contours. For each hole, we need to find the parent outer contour that contains it.

The mapping to this problem is roughly:

iOverlay anchors -> points we want to locate.
iOverlay outer contour bottom edges -> triangle bottom edges.
scan_list -> active edges intersecting the current sweep position.
first_less(...) -> find the closest active edge below the query point.

The algorithm there works like this:

For every hole, take its left-bottom anchor segment.
It is not just a point, but a tiny/real edge starting at the leftmost point.
This helps with degenerate cases when several contours start at the same vertex.
Collect candidate bottom edges from outer contours.

Sort both lists:

• anchors by left endpoint,
• edges by left endpoint, with angle as an additional tie-breaker.

Sweep from left to right.
Add all edges that start before the current anchor.
Keep them in scan_list.
For each anchor, query the nearest active edge below it.

For iTriangle, the same structure can be adapted:

auxiliary point -> query anchor,
triangle bottom edge -> active scan edge,
triangle index -> stored payload in the active set.

[mikwielgus]

Would it make sense to mark these forwarding methods as #[inline]?

I don't know how much of a difference this makes, but sure, done.

It would be nice if this was merged and we continued the non-naive implementation in another PR.

[NailxSharipov]

Would it make sense to mark these forwarding methods as #[inline]?

I don't know how much of a difference this makes, but sure, done.

It would be nice if this was merged and we continued the non-naive implementation in another PR.

inline is useless here, it's not a hot loop fn.

Some important things must be done before merging:

• float must be a wrapper over int.
• carefully process border cases:
  • point on shared/outer edge
  • point on shared/outer vertex

we can extend PointLocationInTriangulation and add this cases

pub enum PointLocationInTriangulation {
    Outside,
    InsideTriangle(TriangleIndex),
    InsideTriangleOnEdge(TriangleIndex),
    InsideTriangleOnVertex(TriangleIndex),
    OnSharedEdge(TriangleIndex, TriangleIndex),
    OnSharedVertex(Vec<TriangleIndex>),
}

If you found my optimization tips complicated, we can skip them for this PR.

[mikwielgus]

I would like to use slightly different naming for clarity and add one more special case for vertices that are shared by exactly two triangles:

pub enum PointLocationInTriangulation {
    Outside,
    InsideTriangle(TriangleIndex),
    OnExteriorEdge(TriangleIndex),
    OnInteriorEdge(TriangleIndex, TriangleIndex),
    OnUnsharedExteriorVertex(TriangleIndex),
    OnSharedExteriorVertex(TriangleIndex, TriangleIndex),
    OnInteriorVertex(Vec<TriangleIndex>),
}

I think the terms exterior and interior are good for the distinction between what lies on the on the outermost boundary and what is inside the triangulation mesh.

[NailxSharipov]

I think the terms exterior and interior are good for the distinction between what lies on the on the outermost boundary and what is inside the triangulation mesh.

Ok. I like it

[mikwielgus]

float must be a wrapper over int.

Can you give me a hint how you want this to be done?

[NailxSharipov]

float must be a wrapper over int.

Can you give me a hint how you want this to be done?

Sorry, looks like I don't submit after review. So You haven't seen my comments I made 3 days ago

[NailxSharipov]

Do we actually need an interior/exterior vertex distinction in this API?

Right now this part does not look correct. See the picture.

A valid way to decide whether a vertex is exterior is to inspect the edges incident to that vertex. If any incident edge is exterior/boundary, then the vertex also exterior.

However, IntTriangulation does not contain this topology information.
We do have enough information when working with RawIntTriangulation and IntDelaunay, because they keep triangle adjacency/topology.

So maybe the locator API should be built over RawIntTriangulation and IntDelaunay.

[mikwielgus]

A valid way to decide whether a vertex is exterior is to inspect the edges incident to that vertex. If any incident edge is exterior/boundary, then the vertex also exterior.

I myself don't need immediate information whether vertex is interior or exterior from this function -- this can be derived by inspecting the incident triangles afterwards.

I did mistakenly add these additional variants, OnUnsharedExteriorVertex, OnSharedExteriorVertex, OnInteriorVertex, to expand on these two enum variants which you suggested:

InsideTriangleOnVertex(TriangleIndex),
OnSharedVertex(Vec<TriangleIndex>)

But I think we don't need to distinguish these two either, do we? I have opted to remove this distinction altogether now. Now it's just

#[derive(Clone, Debug, Eq, PartialEq)]
pub enum PointLocationInTriangulation {
    Outside,
    InsideTriangle(TriangleIndex),
    OnExteriorEdge(TriangleIndex),
    OnInteriorEdge(TriangleIndex, TriangleIndex),
    OnVertex(Vec<TriangleIndex>),
}

Is this alright for you? Or do you want these two separate *On*Vertex enum variants back?

[NailxSharipov]

I think we should write a helper trait for the IntTriangle

fn locate_point(p, a, b, c) -> TriangleLocation (combine code form is_contain_point_exclude_borders and is_contain_point)

enum TriangleLocation {
    Outside,
    Inside,
    OnVertex,
    OnEdge
}