From 8d9de0e1028daed981246174182a39dd917b72bc Mon Sep 17 00:00:00 2001
From: Tavian Barnes <tavianator@tavianator.com>
Date: Sun, 19 Apr 2020 16:38:24 -0400
Subject: metric/vp: Implement vantage-point trees

---
 src/metric/vp.rs | 168 +++++++++++++++++++++++++++++++++++++++++++++++++++++++
 1 file changed, 168 insertions(+)
 create mode 100644 src/metric/vp.rs

(limited to 'src/metric')

diff --git a/src/metric/vp.rs b/src/metric/vp.rs
new file mode 100644
index 0000000..8d5b091
--- /dev/null
+++ b/src/metric/vp.rs
@@ -0,0 +1,168 @@
+//! [Vantage-point trees](https://en.wikipedia.org/wiki/Vantage-point_tree).
+
+use super::{Metric, NearestNeighbors, Neighborhood};
+
+use std::iter::FromIterator;
+
+/// A node in a VP tree.
+#[derive(Debug)]
+struct VpNode<T> {
+    /// The vantage point itself.
+    item: T,
+    /// The radius of this node.
+    radius: f64,
+    /// The subtree inside the radius, if any.
+    inside: Option<Box<Self>>,
+    /// The subtree outside the radius, if any.
+    outside: Option<Box<Self>>,
+}
+
+impl<T: Metric> VpNode<T> {
+    /// Create a new VpNode.
+    fn new(mut items: Vec<T>) -> Option<Box<Self>> {
+        if items.is_empty() {
+            return None;
+        }
+
+        let item = items.pop().unwrap();
+
+        items.sort_by_cached_key(|a| item.distance(a));
+
+        let mid = items.len() / 2;
+        let outside: Vec<T> = items.drain(mid..).collect();
+
+        let radius = items.last().map(|l| item.distance(l).into()).unwrap_or(0.0);
+
+        Some(Box::new(Self {
+            item,
+            radius,
+            inside: Self::new(items),
+            outside: Self::new(outside),
+        }))
+    }
+}
+
+trait VpSearch<'a, T, U, N> {
+    /// Recursively search for nearest neighbors.
+    fn search(&'a self, neighborhood: &mut N);
+
+    /// Search the inside subtree.
+    fn search_inside(&'a self, distance: f64, neighborhood: &mut N);
+
+    /// Search the outside subtree.
+    fn search_outside(&'a self, distance: f64, neighborhood: &mut N);
+}
+
+impl<'a, T, U, N> VpSearch<'a, T, U, N> for VpNode<T>
+where
+    T: 'a,
+    U: Metric<&'a T>,
+    N: Neighborhood<&'a T, U>,
+{
+    fn search(&'a self, neighborhood: &mut N) {
+        let distance = neighborhood.consider(&self.item).into();
+
+        if distance <= self.radius {
+            self.search_inside(distance, neighborhood);
+            self.search_outside(distance, neighborhood);
+        } else {
+            self.search_outside(distance, neighborhood);
+            self.search_inside(distance, neighborhood);
+        }
+    }
+
+    fn search_inside(&'a self, distance: f64, neighborhood: &mut N) {
+        if let Some(inside) = &self.inside {
+            if neighborhood.contains(distance - self.radius) {
+                inside.search(neighborhood);
+            }
+        }
+    }
+
+    fn search_outside(&'a self, distance: f64, neighborhood: &mut N) {
+        if let Some(outside) = &self.outside {
+            if neighborhood.contains(self.radius - distance) {
+                outside.search(neighborhood);
+            }
+        }
+    }
+}
+
+/// A [vantage-point tree](https://en.wikipedia.org/wiki/Vantage-point_tree).
+#[derive(Debug)]
+pub struct VpTree<T> {
+    root: Option<Box<VpNode<T>>>,
+}
+
+impl<T: Metric> FromIterator<T> for VpTree<T> {
+    fn from_iter<I: IntoIterator<Item = T>>(items: I) -> Self {
+        Self {
+            root: VpNode::new(items.into_iter().collect::<Vec<_>>()),
+        }
+    }
+}
+
+impl<T, U> NearestNeighbors<T, U> for VpTree<T>
+where
+    T: Metric,
+    U: Metric<T>,
+{
+    fn search<'a, 'b, N>(&'a self, mut neighborhood: N) -> N
+    where
+        T: 'a,
+        U: 'b,
+        N: Neighborhood<&'a T, &'b U>,
+    {
+        if let Some(root) = &self.root {
+            root.search(&mut neighborhood);
+        }
+        neighborhood
+    }
+}
+
+/// An iterator that moves values out of a VP tree.
+#[derive(Debug)]
+pub struct IntoIter<T> {
+    stack: Vec<Box<VpNode<T>>>,
+}
+
+impl<T> IntoIter<T> {
+    fn new(node: Option<Box<VpNode<T>>>) -> Self {
+        Self {
+            stack: node.into_iter().collect(),
+        }
+    }
+}
+
+impl<T> Iterator for IntoIter<T> {
+    type Item = T;
+
+    fn next(&mut self) -> Option<T> {
+        self.stack.pop().map(|node| {
+            self.stack.extend(node.inside);
+            self.stack.extend(node.outside);
+            node.item
+        })
+    }
+}
+
+impl<T> IntoIterator for VpTree<T> {
+    type Item = T;
+    type IntoIter = IntoIter<T>;
+
+    fn into_iter(self) -> Self::IntoIter {
+        IntoIter::new(self.root)
+    }
+}
+
+#[cfg(test)]
+mod tests {
+    use super::*;
+
+    use crate::metric::tests::test_nearest_neighbors;
+
+    #[test]
+    fn test_vp_tree() {
+        test_nearest_neighbors(VpTree::from_iter);
+    }
+}
-- 
cgit v1.2.3


From 1c560791902a4ef72efa671106d8f6d97fea50c1 Mon Sep 17 00:00:00 2001
From: Tavian Barnes <tavianator@tavianator.com>
Date: Sun, 19 Apr 2020 16:40:29 -0400
Subject: metric/kd: Implement k-d trees

---
 src/metric.rs    |   1 +
 src/metric/kd.rs | 217 +++++++++++++++++++++++++++++++++++++++++++++++++++++++
 2 files changed, 218 insertions(+)
 create mode 100644 src/metric/kd.rs

(limited to 'src/metric')

diff --git a/src/metric.rs b/src/metric.rs
index 549db67..95191da 100644
--- a/src/metric.rs
+++ b/src/metric.rs
@@ -1,5 +1,6 @@
 //! [Metric spaces](https://en.wikipedia.org/wiki/Metric_space).
 
+pub mod kd;
 pub mod vp;
 
 use ordered_float::OrderedFloat;
diff --git a/src/metric/kd.rs b/src/metric/kd.rs
new file mode 100644
index 0000000..ab0cd2e
--- /dev/null
+++ b/src/metric/kd.rs
@@ -0,0 +1,217 @@
+//! [k-d trees](https://en.wikipedia.org/wiki/K-d_tree).
+
+use super::{Metric, NearestNeighbors, Neighborhood};
+
+use ordered_float::OrderedFloat;
+
+use std::iter::FromIterator;
+
+/// A point in Cartesian space.
+pub trait Cartesian {
+    /// Returns the number of dimensions necessary to describe this point.
+    fn dimensions(&self) -> usize;
+
+    /// Returns the value of the `i`th coordinate of this point (`i < self.dimensions()`).
+    fn coordinate(&self, i: usize) -> f64;
+}
+
+/// Blanket [Cartesian] implementation for references.
+impl<'a, T: Cartesian> Cartesian for &'a T {
+    fn dimensions(&self) -> usize {
+        (*self).dimensions()
+    }
+
+    fn coordinate(&self, i: usize) -> f64 {
+        (*self).coordinate(i)
+    }
+}
+
+/// Standard cartesian space.
+impl Cartesian for [f64] {
+    fn dimensions(&self) -> usize {
+        self.len()
+    }
+
+    fn coordinate(&self, i: usize) -> f64 {
+        self[i]
+    }
+}
+
+/// A node in a k-d tree.
+#[derive(Debug)]
+struct KdNode<T> {
+    /// The value stored in this node.
+    item: T,
+    /// The left subtree, if any.
+    left: Option<Box<Self>>,
+    /// The right subtree, if any.
+    right: Option<Box<Self>>,
+}
+
+trait KdSearch<'a, T, U, N> {
+    /// Recursively search for nearest neighbors.
+    fn search(&'a self, i: usize, neighborhood: &mut N);
+
+    /// Search the left subtree.
+    fn search_left(&'a self, i: usize, distance: f64, neighborhood: &mut N);
+
+    /// Search the right subtree.
+    fn search_right(&'a self, i: usize, distance: f64, neighborhood: &mut N);
+}
+
+impl<'a, T, U, N> KdSearch<'a, T, U, N> for KdNode<T>
+where
+    T: 'a + Cartesian,
+    U: Cartesian + Metric<&'a T>,
+    N: Neighborhood<&'a T, U>,
+{
+    fn search(&'a self, i: usize, neighborhood: &mut N) {
+        neighborhood.consider(&self.item);
+
+        let distance = neighborhood.target().coordinate(i) - self.item.coordinate(i);
+        let j = (i + 1) % self.item.dimensions();
+        if distance <= 0.0 {
+            self.search_left(j, distance, neighborhood);
+            self.search_right(j, -distance, neighborhood);
+        } else {
+            self.search_right(j, -distance, neighborhood);
+            self.search_left(j, distance, neighborhood);
+        }
+    }
+
+    fn search_left(&'a self, i: usize, distance: f64, neighborhood: &mut N) {
+        if let Some(left) = &self.left {
+            if neighborhood.contains(distance) {
+                left.search(i, neighborhood);
+            }
+        }
+    }
+
+    fn search_right(&'a self, i: usize, distance: f64, neighborhood: &mut N) {
+        if let Some(right) = &self.right {
+            if neighborhood.contains(distance) {
+                right.search(i, neighborhood);
+            }
+        }
+    }
+}
+
+impl<T: Cartesian> KdNode<T> {
+    /// Create a new KdNode.
+    fn new(i: usize, mut items: Vec<T>) -> Option<Box<Self>> {
+        if items.is_empty() {
+            return None;
+        }
+
+        items.sort_unstable_by_key(|x| OrderedFloat::from(x.coordinate(i)));
+
+        let mid = items.len() / 2;
+        let right: Vec<T> = items.drain((mid + 1)..).collect();
+        let item = items.pop().unwrap();
+        let j = (i + 1) % item.dimensions();
+        Some(Box::new(Self {
+            item,
+            left: Self::new(j, items),
+            right: Self::new(j, right),
+        }))
+    }
+}
+
+/// A [k-d tree](https://en.wikipedia.org/wiki/K-d_tree).
+#[derive(Debug)]
+pub struct KdTree<T> {
+    root: Option<Box<KdNode<T>>>,
+}
+
+impl<T: Cartesian> FromIterator<T> for KdTree<T> {
+    /// Create a new k-d tree from a set of points.
+    fn from_iter<I: IntoIterator<Item = T>>(items: I) -> Self {
+        Self {
+            root: KdNode::new(0, items.into_iter().collect()),
+        }
+    }
+}
+
+impl<T, U> NearestNeighbors<T, U> for KdTree<T>
+where
+    T: Cartesian,
+    U: Cartesian + Metric<T>,
+{
+    fn search<'a, 'b, N>(&'a self, mut neighborhood: N) -> N
+    where
+        T: 'a,
+        U: 'b,
+        N: Neighborhood<&'a T, &'b U>,
+    {
+        if let Some(root) = &self.root {
+            root.search(0, &mut neighborhood);
+        }
+        neighborhood
+    }
+}
+
+/// An iterator that the moves values out of a k-d tree.
+#[derive(Debug)]
+pub struct IntoIter<T> {
+    stack: Vec<Box<KdNode<T>>>,
+}
+
+impl<T> IntoIter<T> {
+    fn new(node: Option<Box<KdNode<T>>>) -> Self {
+        Self {
+            stack: node.into_iter().collect(),
+        }
+    }
+}
+
+impl<T> Iterator for IntoIter<T> {
+    type Item = T;
+
+    fn next(&mut self) -> Option<T> {
+        self.stack.pop().map(|node| {
+            self.stack.extend(node.left);
+            self.stack.extend(node.right);
+            node.item
+        })
+    }
+}
+
+impl<T> IntoIterator for KdTree<T> {
+    type Item = T;
+    type IntoIter = IntoIter<T>;
+
+    fn into_iter(self) -> Self::IntoIter {
+        IntoIter::new(self.root)
+    }
+}
+
+#[cfg(test)]
+mod tests {
+    use super::*;
+
+    use crate::metric::tests::{test_nearest_neighbors, Point};
+    use crate::metric::SquaredDistance;
+
+    impl Metric<[f64]> for Point {
+        type Distance = SquaredDistance;
+
+        fn distance(&self, other: &[f64]) -> Self::Distance {
+            self.0.distance(other)
+        }
+    }
+
+    impl Cartesian for Point {
+        fn dimensions(&self) -> usize {
+            self.0.dimensions()
+        }
+
+        fn coordinate(&self, i: usize) -> f64 {
+            self.0.coordinate(i)
+        }
+    }
+
+    #[test]
+    fn test_kd_tree() {
+        test_nearest_neighbors(KdTree::from_iter);
+    }
+}
-- 
cgit v1.2.3


From 5de377b2b00a927a4f6463c1c5a5fd18606ad006 Mon Sep 17 00:00:00 2001
From: Tavian Barnes <tavianator@tavianator.com>
Date: Thu, 30 Apr 2020 22:51:06 -0400
Subject: metric/kd: Prune k-d tree searches more aggressively

---
 src/metric/kd.rs | 116 +++++++++++++++++++++++++++++++------------------------
 1 file changed, 65 insertions(+), 51 deletions(-)

(limited to 'src/metric')

diff --git a/src/metric/kd.rs b/src/metric/kd.rs
index ab0cd2e..db1b2bd 100644
--- a/src/metric/kd.rs
+++ b/src/metric/kd.rs
@@ -7,7 +7,7 @@ use ordered_float::OrderedFloat;
 use std::iter::FromIterator;
 
 /// A point in Cartesian space.
-pub trait Cartesian {
+pub trait Cartesian: Metric<[f64]> {
     /// Returns the number of dimensions necessary to describe this point.
     fn dimensions(&self) -> usize;
 
@@ -26,6 +26,15 @@ impl<'a, T: Cartesian> Cartesian for &'a T {
     }
 }
 
+/// Blanket [Metric<[f64]>](Metric) implementation for [Cartesian] references.
+impl<'a, T: Cartesian> Metric<[f64]> for &'a T {
+    type Distance = T::Distance;
+
+    fn distance(&self, other: &[f64]) -> Self::Distance {
+        (*self).distance(other)
+    }
+}
+
 /// Standard cartesian space.
 impl Cartesian for [f64] {
     fn dimensions(&self) -> usize {
@@ -37,6 +46,21 @@ impl Cartesian for [f64] {
     }
 }
 
+/// Marker trait for cartesian metric spaces.
+pub trait CartesianMetric<T: ?Sized = Self>:
+    Cartesian + Metric<T, Distance = <Self as Metric<[f64]>>::Distance>
+{
+}
+
+/// Blanket [CartesianMetric] implementation for cartesian spaces with compatible metric distance
+/// types.
+impl<T, U> CartesianMetric<T> for U
+where
+    T: ?Sized,
+    U: ?Sized + Cartesian + Metric<T, Distance = <U as Metric<[f64]>>::Distance>,
+{
+}
+
 /// A node in a k-d tree.
 #[derive(Debug)]
 struct KdNode<T> {
@@ -48,54 +72,6 @@ struct KdNode<T> {
     right: Option<Box<Self>>,
 }
 
-trait KdSearch<'a, T, U, N> {
-    /// Recursively search for nearest neighbors.
-    fn search(&'a self, i: usize, neighborhood: &mut N);
-
-    /// Search the left subtree.
-    fn search_left(&'a self, i: usize, distance: f64, neighborhood: &mut N);
-
-    /// Search the right subtree.
-    fn search_right(&'a self, i: usize, distance: f64, neighborhood: &mut N);
-}
-
-impl<'a, T, U, N> KdSearch<'a, T, U, N> for KdNode<T>
-where
-    T: 'a + Cartesian,
-    U: Cartesian + Metric<&'a T>,
-    N: Neighborhood<&'a T, U>,
-{
-    fn search(&'a self, i: usize, neighborhood: &mut N) {
-        neighborhood.consider(&self.item);
-
-        let distance = neighborhood.target().coordinate(i) - self.item.coordinate(i);
-        let j = (i + 1) % self.item.dimensions();
-        if distance <= 0.0 {
-            self.search_left(j, distance, neighborhood);
-            self.search_right(j, -distance, neighborhood);
-        } else {
-            self.search_right(j, -distance, neighborhood);
-            self.search_left(j, distance, neighborhood);
-        }
-    }
-
-    fn search_left(&'a self, i: usize, distance: f64, neighborhood: &mut N) {
-        if let Some(left) = &self.left {
-            if neighborhood.contains(distance) {
-                left.search(i, neighborhood);
-            }
-        }
-    }
-
-    fn search_right(&'a self, i: usize, distance: f64, neighborhood: &mut N) {
-        if let Some(right) = &self.right {
-            if neighborhood.contains(distance) {
-                right.search(i, neighborhood);
-            }
-        }
-    }
-}
-
 impl<T: Cartesian> KdNode<T> {
     /// Create a new KdNode.
     fn new(i: usize, mut items: Vec<T>) -> Option<Box<Self>> {
@@ -115,6 +91,40 @@ impl<T: Cartesian> KdNode<T> {
             right: Self::new(j, right),
         }))
     }
+
+    /// Recursively search for nearest neighbors.
+    fn search<'a, U, N>(&'a self, i: usize, closest: &mut [f64], neighborhood: &mut N)
+    where
+        T: 'a,
+        U: CartesianMetric<&'a T>,
+        N: Neighborhood<&'a T, U>,
+    {
+        neighborhood.consider(&self.item);
+
+        let target = neighborhood.target();
+        let ti = target.coordinate(i);
+        let si = self.item.coordinate(i);
+        let j = (i + 1) % self.item.dimensions();
+
+        let (near, far) = if ti <= si {
+            (&self.left, &self.right)
+        } else {
+            (&self.right, &self.left)
+        };
+
+        if let Some(near) = near {
+            near.search(j, closest, neighborhood);
+        }
+
+        if let Some(far) = far {
+            let saved = closest[i];
+            closest[i] = si;
+            if neighborhood.contains_distance(target.distance(closest)) {
+                far.search(j, closest, neighborhood);
+            }
+            closest[i] = saved;
+        }
+    }
 }
 
 /// A [k-d tree](https://en.wikipedia.org/wiki/K-d_tree).
@@ -135,7 +145,7 @@ impl<T: Cartesian> FromIterator<T> for KdTree<T> {
 impl<T, U> NearestNeighbors<T, U> for KdTree<T>
 where
     T: Cartesian,
-    U: Cartesian + Metric<T>,
+    U: CartesianMetric<T>,
 {
     fn search<'a, 'b, N>(&'a self, mut neighborhood: N) -> N
     where
@@ -143,8 +153,12 @@ where
         U: 'b,
         N: Neighborhood<&'a T, &'b U>,
     {
+        let target = neighborhood.target();
+        let dims = target.dimensions();
+        let mut closest: Vec<_> = (0..dims).map(|i| target.coordinate(i)).collect();
+
         if let Some(root) = &self.root {
-            root.search(0, &mut neighborhood);
+            root.search(0, &mut closest, &mut neighborhood);
         }
         neighborhood
     }
-- 
cgit v1.2.3


From 9699f4657ecaaf4361448f249e4f2e210a854af4 Mon Sep 17 00:00:00 2001
From: Tavian Barnes <tavianator@tavianator.com>
Date: Thu, 23 Apr 2020 09:55:26 -0400
Subject: metric/vp: Flatten the tree representation

---
 src/metric/vp.rs | 133 +++++++++++++++++++++----------------------------------
 1 file changed, 51 insertions(+), 82 deletions(-)

(limited to 'src/metric')

diff --git a/src/metric/vp.rs b/src/metric/vp.rs
index 8d5b091..fae62e5 100644
--- a/src/metric/vp.rs
+++ b/src/metric/vp.rs
@@ -11,78 +11,62 @@ struct VpNode<T> {
     item: T,
     /// The radius of this node.
     radius: f64,
-    /// The subtree inside the radius, if any.
-    inside: Option<Box<Self>>,
-    /// The subtree outside the radius, if any.
-    outside: Option<Box<Self>>,
+    /// The size of the subtree inside the radius.
+    inside_len: usize,
 }
 
 impl<T: Metric> VpNode<T> {
     /// Create a new VpNode.
-    fn new(mut items: Vec<T>) -> Option<Box<Self>> {
-        if items.is_empty() {
-            return None;
+    fn new(item: T) -> Self {
+        Self {
+            item,
+            radius: 0.0,
+            inside_len: 0,
         }
+    }
 
-        let item = items.pop().unwrap();
-
-        items.sort_by_cached_key(|a| item.distance(a));
-
-        let mid = items.len() / 2;
-        let outside: Vec<T> = items.drain(mid..).collect();
+    /// Build a VP tree recursively.
+    fn build(slice: &mut [VpNode<T>]) {
+        if let Some((node, children)) = slice.split_first_mut() {
+            let item = &node.item;
+            children.sort_by_cached_key(|n| item.distance(&n.item));
 
-        let radius = items.last().map(|l| item.distance(l).into()).unwrap_or(0.0);
+            let (inside, outside) = children.split_at_mut(children.len() / 2);
+            if let Some(last) = inside.last() {
+                node.radius = item.distance(&last.item).into();
+            }
+            node.inside_len = inside.len();
 
-        Some(Box::new(Self {
-            item,
-            radius,
-            inside: Self::new(items),
-            outside: Self::new(outside),
-        }))
+            Self::build(inside);
+            Self::build(outside);
+        }
     }
-}
 
-trait VpSearch<'a, T, U, N> {
     /// Recursively search for nearest neighbors.
-    fn search(&'a self, neighborhood: &mut N);
-
-    /// Search the inside subtree.
-    fn search_inside(&'a self, distance: f64, neighborhood: &mut N);
-
-    /// Search the outside subtree.
-    fn search_outside(&'a self, distance: f64, neighborhood: &mut N);
-}
+    fn recurse<'a, U, N>(slice: &'a [VpNode<T>], neighborhood: &mut N)
+    where
+        T: 'a,
+        U: Metric<&'a T>,
+        N: Neighborhood<&'a T, U>,
+    {
+        let (node, children) = slice.split_first().unwrap();
+        let (inside, outside) = children.split_at(node.inside_len);
 
-impl<'a, T, U, N> VpSearch<'a, T, U, N> for VpNode<T>
-where
-    T: 'a,
-    U: Metric<&'a T>,
-    N: Neighborhood<&'a T, U>,
-{
-    fn search(&'a self, neighborhood: &mut N) {
-        let distance = neighborhood.consider(&self.item).into();
+        let distance = neighborhood.consider(&node.item).into();
 
-        if distance <= self.radius {
-            self.search_inside(distance, neighborhood);
-            self.search_outside(distance, neighborhood);
+        if distance <= node.radius {
+            if !inside.is_empty() && neighborhood.contains(distance - node.radius) {
+                Self::recurse(inside, neighborhood);
+            }
+            if !outside.is_empty() && neighborhood.contains(node.radius - distance) {
+                Self::recurse(outside, neighborhood);
+            }
         } else {
-            self.search_outside(distance, neighborhood);
-            self.search_inside(distance, neighborhood);
-        }
-    }
-
-    fn search_inside(&'a self, distance: f64, neighborhood: &mut N) {
-        if let Some(inside) = &self.inside {
-            if neighborhood.contains(distance - self.radius) {
-                inside.search(neighborhood);
+            if !outside.is_empty() && neighborhood.contains(node.radius - distance) {
+                Self::recurse(outside, neighborhood);
             }
-        }
-    }
-
-    fn search_outside(&'a self, distance: f64, neighborhood: &mut N) {
-        if let Some(outside) = &self.outside {
-            if neighborhood.contains(self.radius - distance) {
-                outside.search(neighborhood);
+            if !inside.is_empty() && neighborhood.contains(distance - node.radius) {
+                Self::recurse(inside, neighborhood);
             }
         }
     }
@@ -90,15 +74,13 @@ where
 
 /// A [vantage-point tree](https://en.wikipedia.org/wiki/Vantage-point_tree).
 #[derive(Debug)]
-pub struct VpTree<T> {
-    root: Option<Box<VpNode<T>>>,
-}
+pub struct VpTree<T>(Vec<VpNode<T>>);
 
 impl<T: Metric> FromIterator<T> for VpTree<T> {
     fn from_iter<I: IntoIterator<Item = T>>(items: I) -> Self {
-        Self {
-            root: VpNode::new(items.into_iter().collect::<Vec<_>>()),
-        }
+        let mut nodes: Vec<_> = items.into_iter().map(VpNode::new).collect();
+        VpNode::build(nodes.as_mut_slice());
+        Self(nodes)
     }
 }
 
@@ -113,36 +95,23 @@ where
         U: 'b,
         N: Neighborhood<&'a T, &'b U>,
     {
-        if let Some(root) = &self.root {
-            root.search(&mut neighborhood);
+        if !self.0.is_empty() {
+            VpNode::recurse(&self.0, &mut neighborhood);
         }
+
         neighborhood
     }
 }
 
 /// An iterator that moves values out of a VP tree.
 #[derive(Debug)]
-pub struct IntoIter<T> {
-    stack: Vec<Box<VpNode<T>>>,
-}
-
-impl<T> IntoIter<T> {
-    fn new(node: Option<Box<VpNode<T>>>) -> Self {
-        Self {
-            stack: node.into_iter().collect(),
-        }
-    }
-}
+pub struct IntoIter<T>(std::vec::IntoIter<VpNode<T>>);
 
 impl<T> Iterator for IntoIter<T> {
     type Item = T;
 
     fn next(&mut self) -> Option<T> {
-        self.stack.pop().map(|node| {
-            self.stack.extend(node.inside);
-            self.stack.extend(node.outside);
-            node.item
-        })
+        self.0.next().map(|n| n.item)
     }
 }
 
@@ -151,7 +120,7 @@ impl<T> IntoIterator for VpTree<T> {
     type IntoIter = IntoIter<T>;
 
     fn into_iter(self) -> Self::IntoIter {
-        IntoIter::new(self.root)
+        IntoIter(self.0.into_iter())
     }
 }
 
-- 
cgit v1.2.3


From e9a81a6d0df149252164003975addf175d5c6f4b Mon Sep 17 00:00:00 2001
From: Tavian Barnes <tavianator@tavianator.com>
Date: Thu, 23 Apr 2020 09:55:13 -0400
Subject: metric/kd: Flatten the tree representation

---
 src/metric/kd.rs | 113 ++++++++++++++++++++++++++-----------------------------
 1 file changed, 54 insertions(+), 59 deletions(-)

(limited to 'src/metric')

diff --git a/src/metric/kd.rs b/src/metric/kd.rs
index db1b2bd..2caf4a3 100644
--- a/src/metric/kd.rs
+++ b/src/metric/kd.rs
@@ -66,61 +66,71 @@ where
 struct KdNode<T> {
     /// The value stored in this node.
     item: T,
-    /// The left subtree, if any.
-    left: Option<Box<Self>>,
-    /// The right subtree, if any.
-    right: Option<Box<Self>>,
+    /// The size of the left subtree.
+    left_len: usize,
 }
 
 impl<T: Cartesian> KdNode<T> {
     /// Create a new KdNode.
-    fn new(i: usize, mut items: Vec<T>) -> Option<Box<Self>> {
-        if items.is_empty() {
-            return None;
+    fn new(item: T) -> Self {
+        Self { item, left_len: 0 }
+    }
+
+    /// Build a k-d tree recursively.
+    fn build(slice: &mut [KdNode<T>], i: usize) {
+        if slice.is_empty() {
+            return;
         }
 
-        items.sort_unstable_by_key(|x| OrderedFloat::from(x.coordinate(i)));
-
-        let mid = items.len() / 2;
-        let right: Vec<T> = items.drain((mid + 1)..).collect();
-        let item = items.pop().unwrap();
-        let j = (i + 1) % item.dimensions();
-        Some(Box::new(Self {
-            item,
-            left: Self::new(j, items),
-            right: Self::new(j, right),
-        }))
+        slice.sort_unstable_by_key(|n| OrderedFloat::from(n.item.coordinate(i)));
+
+        let mid = slice.len() / 2;
+        slice.swap(0, mid);
+
+        let (node, children) = slice.split_first_mut().unwrap();
+        let (left, right) = children.split_at_mut(mid);
+        node.left_len = left.len();
+
+        let j = (i + 1) % node.item.dimensions();
+        Self::build(left, j);
+        Self::build(right, j);
     }
 
     /// Recursively search for nearest neighbors.
-    fn search<'a, U, N>(&'a self, i: usize, closest: &mut [f64], neighborhood: &mut N)
-    where
+    fn recurse<'a, U, N>(
+        slice: &'a [KdNode<T>],
+        i: usize,
+        closest: &mut [f64],
+        neighborhood: &mut N,
+    ) where
         T: 'a,
         U: CartesianMetric<&'a T>,
         N: Neighborhood<&'a T, U>,
     {
-        neighborhood.consider(&self.item);
+        let (node, children) = slice.split_first().unwrap();
+        neighborhood.consider(&node.item);
 
         let target = neighborhood.target();
         let ti = target.coordinate(i);
-        let si = self.item.coordinate(i);
-        let j = (i + 1) % self.item.dimensions();
+        let ni = node.item.coordinate(i);
+        let j = (i + 1) % node.item.dimensions();
 
-        let (near, far) = if ti <= si {
-            (&self.left, &self.right)
+        let (left, right) = children.split_at(node.left_len);
+        let (near, far) = if ti <= ni {
+            (left, right)
         } else {
-            (&self.right, &self.left)
+            (right, left)
         };
 
-        if let Some(near) = near {
-            near.search(j, closest, neighborhood);
+        if !near.is_empty() {
+            Self::recurse(near, j, closest, neighborhood);
         }
 
-        if let Some(far) = far {
+        if !far.is_empty() {
             let saved = closest[i];
-            closest[i] = si;
+            closest[i] = ni;
             if neighborhood.contains_distance(target.distance(closest)) {
-                far.search(j, closest, neighborhood);
+                Self::recurse(far, j, closest, neighborhood);
             }
             closest[i] = saved;
         }
@@ -129,16 +139,14 @@ impl<T: Cartesian> KdNode<T> {
 
 /// A [k-d tree](https://en.wikipedia.org/wiki/K-d_tree).
 #[derive(Debug)]
-pub struct KdTree<T> {
-    root: Option<Box<KdNode<T>>>,
-}
+pub struct KdTree<T>(Vec<KdNode<T>>);
 
 impl<T: Cartesian> FromIterator<T> for KdTree<T> {
     /// Create a new k-d tree from a set of points.
     fn from_iter<I: IntoIterator<Item = T>>(items: I) -> Self {
-        Self {
-            root: KdNode::new(0, items.into_iter().collect()),
-        }
+        let mut nodes: Vec<_> = items.into_iter().map(KdNode::new).collect();
+        KdNode::build(nodes.as_mut_slice(), 0);
+        Self(nodes)
     }
 }
 
@@ -153,40 +161,27 @@ where
         U: 'b,
         N: Neighborhood<&'a T, &'b U>,
     {
-        let target = neighborhood.target();
-        let dims = target.dimensions();
-        let mut closest: Vec<_> = (0..dims).map(|i| target.coordinate(i)).collect();
+        if !self.0.is_empty() {
+            let target = neighborhood.target();
+            let dims = target.dimensions();
+            let mut closest: Vec<_> = (0..dims).map(|i| target.coordinate(i)).collect();
 
-        if let Some(root) = &self.root {
-            root.search(0, &mut closest, &mut neighborhood);
+            KdNode::recurse(&self.0, 0, &mut closest, &mut neighborhood);
         }
+
         neighborhood
     }
 }
 
 /// An iterator that the moves values out of a k-d tree.
 #[derive(Debug)]
-pub struct IntoIter<T> {
-    stack: Vec<Box<KdNode<T>>>,
-}
-
-impl<T> IntoIter<T> {
-    fn new(node: Option<Box<KdNode<T>>>) -> Self {
-        Self {
-            stack: node.into_iter().collect(),
-        }
-    }
-}
+pub struct IntoIter<T>(std::vec::IntoIter<KdNode<T>>);
 
 impl<T> Iterator for IntoIter<T> {
     type Item = T;
 
     fn next(&mut self) -> Option<T> {
-        self.stack.pop().map(|node| {
-            self.stack.extend(node.left);
-            self.stack.extend(node.right);
-            node.item
-        })
+        self.0.next().map(|n| n.item)
     }
 }
 
@@ -195,7 +190,7 @@ impl<T> IntoIterator for KdTree<T> {
     type IntoIter = IntoIter<T>;
 
     fn into_iter(self) -> Self::IntoIter {
-        IntoIter::new(self.root)
+        IntoIter(self.0.into_iter())
     }
 }
 
-- 
cgit v1.2.3


From a4a75059f302de2a00971f1f485fcf4389710628 Mon Sep 17 00:00:00 2001
From: Tavian Barnes <tavianator@tavianator.com>
Date: Sun, 19 Apr 2020 16:55:17 -0400
Subject: metric/forest: Implement dynamized forests

---
 src/metric.rs        |   1 +
 src/metric/forest.rs | 152 +++++++++++++++++++++++++++++++++++++++++++++++++++
 2 files changed, 153 insertions(+)
 create mode 100644 src/metric/forest.rs

(limited to 'src/metric')

diff --git a/src/metric.rs b/src/metric.rs
index 95191da..7a5f5f7 100644
--- a/src/metric.rs
+++ b/src/metric.rs
@@ -1,5 +1,6 @@
 //! [Metric spaces](https://en.wikipedia.org/wiki/Metric_space).
 
+pub mod forest;
 pub mod kd;
 pub mod vp;
 
diff --git a/src/metric/forest.rs b/src/metric/forest.rs
new file mode 100644
index 0000000..f23c451
--- /dev/null
+++ b/src/metric/forest.rs
@@ -0,0 +1,152 @@
+//! [Dynamization](https://en.wikipedia.org/wiki/Dynamization) for nearest neighbor search.
+
+use super::kd::KdTree;
+use super::vp::VpTree;
+use super::{Metric, NearestNeighbors, Neighborhood};
+
+use std::iter::{Extend, Flatten, FromIterator};
+
+/// A dynamic wrapper for a static nearest neighbor search data structure.
+///
+/// This type applies [dynamization](https://en.wikipedia.org/wiki/Dynamization) to an arbitrary
+/// nearest neighbor search structure `T`, allowing new items to be added dynamically.
+#[derive(Debug)]
+pub struct Forest<T>(Vec<Option<T>>);
+
+impl<T, U> Forest<U>
+where
+    U: FromIterator<T> + IntoIterator<Item = T>,
+{
+    /// Create a new empty forest.
+    pub fn new() -> Self {
+        Self(Vec::new())
+    }
+
+    /// Add a new item to the forest.
+    pub fn push(&mut self, item: T) {
+        let mut items = vec![item];
+
+        for slot in &mut self.0 {
+            match slot.take() {
+                // Collect the items from any trees we encounter...
+                Some(tree) => {
+                    items.extend(tree);
+                }
+                // ... and put them all in the first empty slot
+                None => {
+                    *slot = Some(items.into_iter().collect());
+                    return;
+                }
+            }
+        }
+
+        self.0.push(Some(items.into_iter().collect()));
+    }
+
+    /// Get the number of items in the forest.
+    pub fn len(&self) -> usize {
+        let mut len = 0;
+        for (i, slot) in self.0.iter().enumerate() {
+            if slot.is_some() {
+                len |= 1 << i;
+            }
+        }
+        len
+    }
+}
+
+impl<T, U> Extend<T> for Forest<U>
+where
+    U: FromIterator<T> + IntoIterator<Item = T>,
+{
+    fn extend<I: IntoIterator<Item = T>>(&mut self, items: I) {
+        for item in items {
+            self.push(item);
+        }
+    }
+}
+
+impl<T, U> FromIterator<T> for Forest<U>
+where
+    U: FromIterator<T> + IntoIterator<Item = T>,
+{
+    fn from_iter<I: IntoIterator<Item = T>>(items: I) -> Self {
+        let mut forest = Self::new();
+        forest.extend(items);
+        forest
+    }
+}
+
+type IntoIterImpl<T> = Flatten<Flatten<std::vec::IntoIter<Option<T>>>>;
+
+/// An iterator that moves items out of a forest.
+pub struct IntoIter<T: IntoIterator>(IntoIterImpl<T>);
+
+impl<T: IntoIterator> Iterator for IntoIter<T> {
+    type Item = T::Item;
+
+    fn next(&mut self) -> Option<Self::Item> {
+        self.0.next()
+    }
+}
+
+impl<T: IntoIterator> IntoIterator for Forest<T> {
+    type Item = T::Item;
+    type IntoIter = IntoIter<T>;
+
+    fn into_iter(self) -> Self::IntoIter {
+        IntoIter(self.0.into_iter().flatten().flatten())
+    }
+}
+
+impl<T, U, V> NearestNeighbors<T, U> for Forest<V>
+where
+    U: Metric<T>,
+    V: NearestNeighbors<T, U>,
+{
+    fn search<'a, 'b, N>(&'a self, neighborhood: N) -> N
+    where
+        T: 'a,
+        U: 'b,
+        N: Neighborhood<&'a T, &'b U>,
+    {
+        self.0
+            .iter()
+            .flatten()
+            .fold(neighborhood, |n, t| t.search(n))
+    }
+}
+
+/// A forest of k-d trees.
+pub type KdForest<T> = Forest<KdTree<T>>;
+
+/// A forest of vantage-point trees.
+pub type VpForest<T> = Forest<VpTree<T>>;
+
+#[cfg(test)]
+mod tests {
+    use super::*;
+
+    use crate::metric::tests::test_nearest_neighbors;
+    use crate::metric::ExhaustiveSearch;
+
+    #[test]
+    fn test_exhaustive_forest() {
+        test_nearest_neighbors(Forest::<ExhaustiveSearch<_>>::from_iter);
+    }
+
+    #[test]
+    fn test_forest_forest() {
+        test_nearest_neighbors(Forest::<Forest<ExhaustiveSearch<_>>>::from_iter);
+    }
+
+    #[test]
+    fn test_kd_forest() {
+        test_nearest_neighbors(KdForest::from_iter);
+    }
+
+    #[test]
+    fn test_vp_forest() {
+        test_nearest_neighbors(VpForest::from_iter);
+    }
+}
-- 
cgit v1.2.3


From b4a39a3f22fac361f6a535d281eee5586078281b Mon Sep 17 00:00:00 2001
From: Tavian Barnes <tavianator@tavianator.com>
Date: Thu, 23 Apr 2020 14:58:47 -0400
Subject: metric/forest: Optimize bulk insertion

---
 src/metric/forest.rs | 47 +++++++++++++++++++++++++++--------------------
 1 file changed, 27 insertions(+), 20 deletions(-)

(limited to 'src/metric')

diff --git a/src/metric/forest.rs b/src/metric/forest.rs
index f23c451..29b6f55 100644
--- a/src/metric/forest.rs
+++ b/src/metric/forest.rs
@@ -4,7 +4,7 @@ use super::kd::KdTree;
 use super::vp::VpTree;
 use super::{Metric, NearestNeighbors, Neighborhood};
 
-use std::iter::{Extend, Flatten, FromIterator};
+use std::iter::{self, Extend, Flatten, FromIterator};
 
 /// A dynamic wrapper for a static nearest neighbor search data structure.
 ///
@@ -24,23 +24,7 @@ where
 
     /// Add a new item to the forest.
     pub fn push(&mut self, item: T) {
-        let mut items = vec![item];
-
-        for slot in &mut self.0 {
-            match slot.take() {
-                // Collect the items from any trees we encounter...
-                Some(tree) => {
-                    items.extend(tree);
-                }
-                // ... and put them all in the first empty slot
-                None => {
-                    *slot = Some(items.into_iter().collect());
-                    return;
-                }
-            }
-        }
-
-        self.0.push(Some(items.into_iter().collect()));
+        self.extend(iter::once(item));
     }
 
     /// Get the number of items in the forest.
@@ -60,9 +44,32 @@ where
     U: FromIterator<T> + IntoIterator<Item = T>,
 {
     fn extend<I: IntoIterator<Item = T>>(&mut self, items: I) {
-        for item in items {
-            self.push(item);
+        let mut vec: Vec<_> = items.into_iter().collect();
+        let new_len = self.len() + vec.len();
+
+        for i in 0.. {
+            let bit = 1 << i;
+
+            if bit > new_len {
+                break;
+            }
+
+            if i >= self.0.len() {
+                self.0.push(None);
+            }
+
+            if new_len & bit == 0 {
+                if let Some(tree) = self.0[i].take() {
+                    vec.extend(tree);
+                }
+            } else if self.0[i].is_none() {
+                let offset = vec.len() - bit;
+                self.0[i] = Some(vec.drain(offset..).collect());
+            }
         }
+
+        debug_assert!(vec.is_empty());
+        debug_assert!(self.len() == new_len);
     }
 }
 
-- 
cgit v1.2.3


From e53ace3e69ed6bacedb7de345df10d3e575a291e Mon Sep 17 00:00:00 2001
From: Tavian Barnes <tavianator@tavianator.com>
Date: Sun, 19 Apr 2020 16:58:52 -0400
Subject: metric/soft: Implement soft deletes

---
 src/metric.rs      |   1 +
 src/metric/soft.rs | 282 +++++++++++++++++++++++++++++++++++++++++++++++++++++
 2 files changed, 283 insertions(+)
 create mode 100644 src/metric/soft.rs

(limited to 'src/metric')

diff --git a/src/metric.rs b/src/metric.rs
index 7a5f5f7..b46c8da 100644
--- a/src/metric.rs
+++ b/src/metric.rs
@@ -2,6 +2,7 @@
 
 pub mod forest;
 pub mod kd;
+pub mod soft;
 pub mod vp;
 
 use ordered_float::OrderedFloat;
diff --git a/src/metric/soft.rs b/src/metric/soft.rs
new file mode 100644
index 0000000..0d7dcdb
--- /dev/null
+++ b/src/metric/soft.rs
@@ -0,0 +1,282 @@
+//! [Soft deletion](https://en.wiktionary.org/wiki/soft_deletion) for nearest neighbor search.
+
+use super::forest::{KdForest, VpForest};
+use super::kd::KdTree;
+use super::vp::VpTree;
+use super::{Metric, NearestNeighbors, Neighborhood};
+
+use std::iter;
+use std::iter::FromIterator;
+use std::mem;
+
+/// A trait for objects that can be soft-deleted.
+pub trait SoftDelete {
+    /// Check whether this item is deleted.
+    fn is_deleted(&self) -> bool;
+}
+
+/// Blanket [SoftDelete] implementation for references.
+impl<'a, T: SoftDelete> SoftDelete for &'a T {
+    fn is_deleted(&self) -> bool {
+        (*self).is_deleted()
+    }
+}
+
+/// [Neighborhood] wrapper that ignores soft-deleted items.
+#[derive(Debug)]
+struct SoftNeighborhood<N>(N);
+
+impl<T, U, N> Neighborhood<T, U> for SoftNeighborhood<N>
+where
+    T: SoftDelete,
+    U: Metric<T>,
+    N: Neighborhood<T, U>,
+{
+    fn target(&self) -> U {
+        self.0.target()
+    }
+
+    fn contains(&self, distance: f64) -> bool {
+        self.0.contains(distance)
+    }
+
+    fn contains_distance(&self, distance: U::Distance) -> bool {
+        self.0.contains_distance(distance)
+    }
+
+    fn consider(&mut self, item: T) -> U::Distance {
+        if item.is_deleted() {
+            self.target().distance(&item)
+        } else {
+            self.0.consider(item)
+        }
+    }
+}
+
+/// A [NearestNeighbors] implementation that supports [soft deletes](https://en.wiktionary.org/wiki/soft_deletion).
+#[derive(Debug)]
+pub struct SoftSearch<T>(T);
+
+impl<T, U> SoftSearch<U>
+where
+    T: SoftDelete,
+    U: FromIterator<T> + IntoIterator<Item = T>,
+{
+    /// Create a new empty soft index.
+    pub fn new() -> Self {
+        Self(iter::empty().collect())
+    }
+
+    /// Push a new item into this index.
+    pub fn push(&mut self, item: T)
+    where
+        U: Extend<T>,
+    {
+        self.0.extend(iter::once(item));
+    }
+
+    /// Rebuild this index, discarding deleted items.
+    pub fn rebuild(&mut self) {
+        let items = mem::replace(&mut self.0, iter::empty().collect());
+        self.0 = items.into_iter().filter(|e| !e.is_deleted()).collect();
+    }
+}
+
+impl<T, U: Extend<T>> Extend<T> for SoftSearch<U> {
+    fn extend<I: IntoIterator<Item = T>>(&mut self, iter: I) {
+        self.0.extend(iter);
+    }
+}
+
+impl<T, U: FromIterator<T>> FromIterator<T> for SoftSearch<U> {
+    fn from_iter<I: IntoIterator<Item = T>>(iter: I) -> Self {
+        Self(U::from_iter(iter))
+    }
+}
+
+impl<T: IntoIterator> IntoIterator for SoftSearch<T> {
+    type Item = T::Item;
+    type IntoIter = T::IntoIter;
+
+    fn into_iter(self) -> Self::IntoIter {
+        self.0.into_iter()
+    }
+}
+
+impl<T, U, V> NearestNeighbors<T, U> for SoftSearch<V>
+where
+    T: SoftDelete,
+    U: Metric<T>,
+    V: NearestNeighbors<T, U>,
+{
+    fn search<'a, 'b, N>(&'a self, neighborhood: N) -> N
+    where
+        T: 'a,
+        U: 'b,
+        N: Neighborhood<&'a T, &'b U>,
+    {
+        self.0.search(SoftNeighborhood(neighborhood)).0
+    }
+}
+
+/// A k-d forest that supports soft deletes.
+pub type SoftKdForest<T> = SoftSearch<KdForest<T>>;
+
+/// A k-d tree that supports soft deletes.
+pub type SoftKdTree<T> = SoftSearch<KdTree<T>>;
+
+/// A VP forest that supports soft deletes.
+pub type SoftVpForest<T> = SoftSearch<VpForest<T>>;
+
+/// A VP tree that supports soft deletes.
+pub type SoftVpTree<T> = SoftSearch<VpTree<T>>;
+
+#[cfg(test)]
+mod tests {
+    use super::*;
+
+    use crate::metric::kd::Cartesian;
+    use crate::metric::tests::Point;
+    use crate::metric::Neighbor;
+
+    #[derive(Debug, PartialEq)]
+    struct SoftPoint {
+        point: Point,
+        deleted: bool,
+    }
+
+    impl SoftPoint {
+        fn new(x: f64, y: f64, z: f64) -> Self {
+            Self {
+                point: Point([x, y, z]),
+                deleted: false,
+            }
+        }
+
+        fn deleted(x: f64, y: f64, z: f64) -> Self {
+            Self {
+                point: Point([x, y, z]),
+                deleted: true,
+            }
+        }
+    }
+
+    impl SoftDelete for SoftPoint {
+        fn is_deleted(&self) -> bool {
+            self.deleted
+        }
+    }
+
+    impl Metric for SoftPoint {
+        type Distance = <Point as Metric>::Distance;
+
+        fn distance(&self, other: &Self) -> Self::Distance {
+            self.point.distance(&other.point)
+        }
+    }
+
+    impl Metric<[f64]> for SoftPoint {
+        type Distance = <Point as Metric>::Distance;
+
+        fn distance(&self, other: &[f64]) -> Self::Distance {
+            self.point.distance(other)
+        }
+    }
+
+    impl Cartesian for SoftPoint {
+        fn dimensions(&self) -> usize {
+            self.point.dimensions()
+        }
+
+        fn coordinate(&self, i: usize) -> f64 {
+            self.point.coordinate(i)
+        }
+    }
+
+    impl Metric<SoftPoint> for Point {
+        type Distance = <Point as Metric>::Distance;
+
+        fn distance(&self, other: &SoftPoint) -> Self::Distance {
+            self.distance(&other.point)
+        }
+    }
+
+    fn test_index<T>(index: &T)
+    where
+        T: NearestNeighbors<SoftPoint, Point>,
+    {
+        let target = Point([0.0, 0.0, 0.0]);
+
+        assert_eq!(
+            index.nearest(&target),
+            Some(Neighbor::new(&SoftPoint::new(1.0, 2.0, 2.0), 3.0))
+        );
+
+        assert_eq!(index.nearest_within(&target, 2.0), None);
+        assert_eq!(
+            index.nearest_within(&target, 4.0),
+            Some(Neighbor::new(&SoftPoint::new(1.0, 2.0, 2.0), 3.0))
+        );
+
+        assert_eq!(
+            index.k_nearest(&target, 3),
+            vec![
+                Neighbor::new(&SoftPoint::new(1.0, 2.0, 2.0), 3.0),
+                Neighbor::new(&SoftPoint::new(3.0, 4.0, 0.0), 5.0),
+                Neighbor::new(&SoftPoint::new(2.0, 3.0, 6.0), 7.0),
+            ]
+        );
+
+        assert_eq!(
+            index.k_nearest_within(&target, 3, 6.0),
+            vec![
+                Neighbor::new(&SoftPoint::new(1.0, 2.0, 2.0), 3.0),
+                Neighbor::new(&SoftPoint::new(3.0, 4.0, 0.0), 5.0),
+            ]
+        );
+        assert_eq!(
+            index.k_nearest_within(&target, 3, 8.0),
+            vec![
+                Neighbor::new(&SoftPoint::new(1.0, 2.0, 2.0), 3.0),
+                Neighbor::new(&SoftPoint::new(3.0, 4.0, 0.0), 5.0),
+                Neighbor::new(&SoftPoint::new(2.0, 3.0, 6.0), 7.0),
+            ]
+        );
+    }
+
+    fn test_soft_index<T>(index: &mut SoftSearch<T>)
+    where
+        T: Extend<SoftPoint>,
+        T: FromIterator<SoftPoint>,
+        T: IntoIterator<Item = SoftPoint>,
+        T: NearestNeighbors<SoftPoint, Point>,
+    {
+        let points = vec![
+            SoftPoint::deleted(0.0, 0.0, 0.0),
+            SoftPoint::new(3.0, 4.0, 0.0),
+            SoftPoint::new(5.0, 0.0, 12.0),
+            SoftPoint::new(0.0, 8.0, 15.0),
+            SoftPoint::new(1.0, 2.0, 2.0),
+            SoftPoint::new(2.0, 3.0, 6.0),
+            SoftPoint::new(4.0, 4.0, 7.0),
+        ];
+
+        for point in points {
+            index.push(point);
+        }
+        test_index(index);
+
+        index.rebuild();
+        test_index(index);
+    }
+
+    #[test]
+    fn test_soft_kd_forest() {
+        test_soft_index(&mut SoftKdForest::new());
+    }
+
+    #[test]
+    fn test_soft_vp_forest() {
+        test_soft_index(&mut SoftVpForest::new());
+    }
+}
-- 
cgit v1.2.3


From 462daaffd5ec720ed80a2e7b1f445a73cabf5833 Mon Sep 17 00:00:00 2001
From: Tavian Barnes <tavianator@tavianator.com>
Date: Sun, 19 Apr 2020 16:59:24 -0400
Subject: metric/approx: Implement approximate nearest neighbor search

---
 src/metric.rs        |   1 +
 src/metric/approx.rs | 131 +++++++++++++++++++++++++++++++++++++++++++++++++++
 2 files changed, 132 insertions(+)
 create mode 100644 src/metric/approx.rs

(limited to 'src/metric')

diff --git a/src/metric.rs b/src/metric.rs
index b46c8da..268aefd 100644
--- a/src/metric.rs
+++ b/src/metric.rs
@@ -1,5 +1,6 @@
 //! [Metric spaces](https://en.wikipedia.org/wiki/Metric_space).
 
+pub mod approx;
 pub mod forest;
 pub mod kd;
 pub mod soft;
diff --git a/src/metric/approx.rs b/src/metric/approx.rs
new file mode 100644
index 0000000..c23f9c7
--- /dev/null
+++ b/src/metric/approx.rs
@@ -0,0 +1,131 @@
+//! [Approximate nearest neighbor search](https://en.wikipedia.org/wiki/Nearest_neighbor_search#Approximate_nearest_neighbor).
+
+use super::{Metric, NearestNeighbors, Neighborhood};
+
+/// An approximate [Neighborhood], for approximate nearest neighbor searches.
+#[derive(Debug)]
+struct ApproximateNeighborhood<N> {
+    inner: N,
+    ratio: f64,
+    limit: usize,
+}
+
+impl<N> ApproximateNeighborhood<N> {
+    fn new(inner: N, ratio: f64, limit: usize) -> Self {
+        Self {
+            inner,
+            ratio,
+            limit,
+        }
+    }
+}
+
+impl<T, U, N> Neighborhood<T, U> for ApproximateNeighborhood<N>
+where
+    U: Metric<T>,
+    N: Neighborhood<T, U>,
+{
+    fn target(&self) -> U {
+        self.inner.target()
+    }
+
+    fn contains(&self, distance: f64) -> bool {
+        if self.limit > 0 {
+            self.inner.contains(self.ratio * distance)
+        } else {
+            false
+        }
+    }
+
+    fn contains_distance(&self, distance: U::Distance) -> bool {
+        self.contains(self.ratio * distance.into())
+    }
+
+    fn consider(&mut self, item: T) -> U::Distance {
+        self.limit = self.limit.saturating_sub(1);
+        self.inner.consider(item)
+    }
+}
+
+/// An [approximate nearest neighbor search](https://en.wikipedia.org/wiki/Nearest_neighbor_search#Approximate_nearest_neighbor)
+/// index.
+///
+/// This wrapper converts an exact nearest neighbor search algorithm into an approximate one by
+/// modifying the behavior of [Neighborhood::contains].  The approximation is controlled by two
+/// parameters:
+///
+/// * `ratio`: The [nearest neighbor distance ratio](https://en.wikipedia.org/wiki/Nearest_neighbor_search#Nearest_neighbor_distance_ratio),
+///   which controls how much closer new candidates must be to be considered.  For example, a ratio
+///   of 2.0 means that a neighbor must be less than half of the current threshold away to be
+///   considered.  A ratio of 1.0 means an exact search.
+///
+/// * `limit`: A limit on the number of candidates to consider.  Typical nearest neighbor algorithms
+///   find a close match quickly, so setting a limit bounds the worst-case search time while keeping
+///   good accuracy.
+#[derive(Debug)]
+pub struct ApproximateSearch<T> {
+    inner: T,
+    ratio: f64,
+    limit: usize,
+}
+
+impl<T> ApproximateSearch<T> {
+    /// Create a new ApproximateSearch index.
+    ///
+    /// * `inner`: The [NearestNeighbors] implementation to wrap.
+    /// * `ratio`: The nearest neighbor distance ratio.
+    /// * `limit`: The maximum number of results to consider.
+    pub fn new(inner: T, ratio: f64, limit: usize) -> Self {
+        Self {
+            inner,
+            ratio,
+            limit,
+        }
+    }
+}
+
+impl<T, U, V> NearestNeighbors<T, U> for ApproximateSearch<V>
+where
+    U: Metric<T>,
+    V: NearestNeighbors<T, U>,
+{
+    fn search<'a, 'b, N>(&'a self, neighborhood: N) -> N
+    where
+        T: 'a,
+        U: 'b,
+        N: Neighborhood<&'a T, &'b U>,
+    {
+        self.inner
+            .search(ApproximateNeighborhood::new(
+                neighborhood,
+                self.ratio,
+                self.limit,
+            ))
+            .inner
+    }
+}
+
+#[cfg(test)]
+mod tests {
+    use super::*;
+
+    use crate::metric::kd::KdTree;
+    use crate::metric::tests::test_nearest_neighbors;
+    use crate::metric::vp::VpTree;
+
+    use std::iter::FromIterator;
+
+    #[test]
+    fn test_approx_kd_tree() {
+        test_nearest_neighbors(|iter| {
+            ApproximateSearch::new(KdTree::from_iter(iter), 1.0, std::usize::MAX)
+        });
+    }
+
+    #[test]
+    fn test_approx_vp_tree() {
+        test_nearest_neighbors(|iter| {
+            ApproximateSearch::new(VpTree::from_iter(iter), 1.0, std::usize::MAX)
+        });
+    }
+}
-- 
cgit v1.2.3