From 748157da3a8929a08582227a1b2e6870cb501738 Mon Sep 17 00:00:00 2001
From: "google-labs-jules[bot]"
 <161369871+google-labs-jules[bot]@users.noreply.github.com>
Date: Tue, 9 Jun 2026 09:45:40 +0000
Subject: [PATCH] Optimize spatial scans in clustering with inline squared
 distance comparisons

Replaced exact Euclidean distance calculations (`libm::sqrt`) in `auto_k_selection`
and `DBSCAN::region_query` with inline squared distance comparisons and inner
loop early exits. Also fixed a component counting logic bug in `auto_k_selection`.

Co-authored-by: teerthsharma <78080953+teerthsharma@users.noreply.github.com>
---
 .jules/bolt.md                          |  4 +++
 crates/aether-core/src/ml/clustering.rs | 40 +++++++++++++++++++------
 2 files changed, 35 insertions(+), 9 deletions(-)

diff --git a/.jules/bolt.md b/.jules/bolt.md
index 2ac6922..02b2903 100644
--- a/.jules/bolt.md
+++ b/.jules/bolt.md
@@ -5,3 +5,7 @@
 ## 2026-05-01 - Avoid High-Level Tensor Ops in Scalar Reductions
 **Learning:** High-level `Tensor` operations like `sub()` and `mul()` trigger intermediate heap allocations for shape and stride metadata. When computing scalar reductions (like MSE, distances, or loss functions), using these operations introduces severe memory overhead inside hot loops. Attempting to use `.min()` length truncation as a safeguard is an anti-pattern as it masks shape mismatch errors.
 **Action:** For scalar reductions, assert shape equality (`assert_eq!(a.shape, b.shape)`) and perform a single-pass iteration directly over the underlying borrowed data arrays (`a.data.borrow()`) to eliminate intermediate allocations and safely compute the result.
+
+## 2026-06-09 - Inline Squared Distance Comparisons in Clustering
+**Learning:** Exact distance calculations using libm::sqrt in spatial scan algorithms like DBSCAN and component counting are a bottleneck. Calculating exact Euclidean distance does unnecessary work when comparing against a constant threshold, and is further impacted by std/libm overhead.
+**Action:** In clustering modules, optimize spatial scans by using inline squared distance comparisons with an inner early exit condition (e.g., `!(sum < eps_sq)`). Always explicitly reject invalid (NaN/negative) thresholds before the loop.
diff --git a/crates/aether-core/src/ml/clustering.rs b/crates/aether-core/src/ml/clustering.rs
index 7e043ad..32fc8e0 100644
--- a/crates/aether-core/src/ml/clustering.rs
+++ b/crates/aether-core/src/ml/clustering.rs
@@ -13,7 +13,6 @@
 // ═══════════════════════════════════════════════════════════════════════════════
 //
 
-
 #![allow(dead_code)]
 
 use libm::{fabs, sqrt};
@@ -274,10 +273,15 @@ impl<const D: usize> KMeans<D> {
 
 /// Automatically determine optimal K using topological analysis
 pub fn auto_k_selection<const D: usize>(data: &[[f64; D]], n: usize, epsilon: f64) -> usize {
+    let n = n.min(MAX_POINTS);
+    if !(epsilon > 0.0) {
+        return n.clamp(1, MAX_POINTS);
+    }
+    let eps_sq = epsilon * epsilon;
+
     // Build epsilon-neighborhood graph and count connected components (β₀)
-    let mut components = n.min(MAX_POINTS);
+    let mut components = 0;
     let mut visited = [false; MAX_POINTS];
-    let n = n.min(MAX_POINTS);
 
     for start in 0..n {
         if visited[start] {
@@ -285,7 +289,7 @@ pub fn auto_k_selection<const D: usize>(data: &[[f64; D]], n: usize, epsilon: f6
         }
 
         // BFS from this point
-        components -= 1;
+        components += 1;
         let mut stack = [0usize; 64];
         let mut top = 1;
         stack[0] = start;
@@ -302,16 +306,21 @@ pub fn auto_k_selection<const D: usize>(data: &[[f64; D]], n: usize, epsilon: f6
             // Add neighbors
             for i in 0..n {
                 if !visited[i] && i != current {
-                    let dist = distance(&data[current], &data[i]);
-                    if dist < epsilon && top < 64 {
+                    let mut sum = 0.0;
+                    for d in 0..D {
+                        let diff = data[current][d] - data[i][d];
+                        sum += diff * diff;
+                        if !(sum < eps_sq) {
+                            break;
+                        }
+                    }
+                    if sum < eps_sq && top < 64 {
                         stack[top] = i;
                         top += 1;
                     }
                 }
             }
         }
-
-        components += 1;
     }
 
     // β₀ = number of connected components = suggested K
@@ -450,8 +459,21 @@ impl<const D: usize> DBSCAN<D> {
     ) -> heapless::Vec<usize, MAX_POINTS> {
         let mut neighbors = heapless::Vec::new();
 
+        if !(self.epsilon > 0.0) {
+            return neighbors;
+        }
+        let eps_sq = self.epsilon * self.epsilon;
+
         for j in 0..n {
-            if distance(&data[i], &data[j]) <= self.epsilon {
+            let mut sum = 0.0;
+            for d in 0..D {
+                let diff = data[i][d] - data[j][d];
+                sum += diff * diff;
+                if !(sum <= eps_sq) {
+                    break;
+                }
+            }
+            if sum <= eps_sq {
                 let _ = neighbors.push(j);
             }
         }