Add PointCloudHashGrid: uniform-grid counterpart to PointCloudBVH

Docs/Sphinx/source/ImplemPointCloud.rst

@@ -0,0 +1,65 @@
+.. _Chap:ImplemPointCloud:
+
+Point clouds
+============
+
+EBGeometry provides two turnkey classes for nearest-neighbor and closest-point work over a point
+cloud: ``PointCloudBVH`` (:file:`Source/EBGeometry_PointCloudBVH.hpp`) and ``PointCloudHashGrid``
+(:file:`Source/EBGeometry_PointCloudHashGrid.hpp`). They answer the same queries and expose the
+**same public interface** -- the same ``Hit`` result type, the same ``closestPoint`` /
+``closestPoints`` / ``nearestNeighbor`` / ``nearestNeighbors`` / ``allNearestNeighbors`` query
+methods, the same ``O(N)`` brute-force reference queries, and the same ``position()`` /
+``metadata()`` accessors -- so the two are drop-in interchangeable. They differ only in the spatial
+acceleration structure they build: a hierarchical tree, or a uniform grid. See :ref:`Chap:PointCloud`
+for the conceptual picture and the trade-off between the two.
+
+Both are built directly from a raw cloud -- point positions plus a parallel array of user metadata --
+and both return the matched point's **cloud index** (its position in the input arrays) together with
+the squared distance; the user metadata is reachable through ``metadata()``. ``closestPoint`` /
+``closestPoints`` answer an arbitrary external query point, while ``nearestNeighbor`` /
+``nearestNeighbors`` (and the batch ``allNearestNeighbors``) answer a point already in the cloud,
+excluding it from its own result and seeding the search from the group it lives in -- a strictly
+cheaper search an external point cannot use (see :ref:`Chap:PointCloud`).
+
+Each accelerated query also has an ``O(N)`` brute-force counterpart -- ``closestPointBruteForce`` /
+``closestPointsBruteForce`` / ``nearestNeighborBruteForce`` / ``nearestNeighborsBruteForce`` -- that
+answers the same question by a full linear scan. These are reference implementations for testing and
+debugging (verify an accelerated result against ground truth, or A/B-test a suspected bug against an
+unaccelerated path) and are not meant for production queries.
+
+PointCloudBVH
+-------------
+
+``PointCloudBVH`` specializes the :ref:`Chap:ImplemBVH` machinery: it is a subclass of
+``BVH::PackedBVH`` with two things the general path does not offer -- a much cheaper **index-based
+build**, and **turnkey query methods** that hide ``pruneTraverse()`` entirely. It is built by
+partitioning an index permutation in place with a longest-axis midpoint split and packing the
+``PointAoSoA`` leaves inline (no intermediate primitive list, no ``shared_ptr``, no separate packing
+pass), which is several times faster to build than a full Surface-Area-Heuristic tree and, for
+near-uniform clouds, just as tight to query. Because it is a BVH, it can also be composed as a
+primitive inside an outer BVH or CSG tree.
+
+See the `PointCloudBVH doxygen page
+<doxygen/html/classEBGeometry_1_1PointCloudBVH.html>`__ for the full interface, and
+``Examples/ClosestPointBVH`` / ``Examples/NearestNeighborBVH`` for worked usage.
+
+PointCloudHashGrid
+------------------
+
+``PointCloudHashGrid`` circumvents the tree entirely and stores the cloud in a **uniform grid**.
+Points are counting-sorted into a dense array of fixed-size cells (a CSR bucket array keyed by
+integer cell coordinates) -- an ``O(N)`` build with no recursive partitioning and no tree nodes. A
+query is an **expanding-shell** search outward from the query point's cell (Chebyshev radius
+0, 1, 2, ...); it stops, exactly and without ever missing a neighbor, as soon as the k-th best
+distance found is closer than any unvisited cell can hold. With the default cell size (~1 point/cell)
+that is almost always one or two shells.
+
+For a near-uniform cloud the grid both builds and queries faster than the BVH; for a strongly
+clustered or multi-scale cloud a single global cell size is a poor fit and ``PointCloudBVH`` is the
+better choice. The grid is also bounded-domain (dense cells sized to the bounding box, ``O(N)``
+memory for a compact cloud) and serves only point queries; unlike ``PointCloudBVH`` it cannot be
+composed as a primitive inside an outer BVH/CSG.
+
+See the `PointCloudHashGrid doxygen page
+<doxygen/html/classEBGeometry_1_1PointCloudHashGrid.html>`__ for the full interface, and
+``Examples/NearestNeighborHashGrid`` for a worked comparison against the BVH example.

Docs/Sphinx/source/PointCloud.rst

@@ -1,32 +1,56 @@
-.. _Chap:ImplemPointCloud:
+.. _Chap:PointCloud:
 
 Point clouds
 ============
 
-For nearest-neighbor and closest-point work over particle clouds, EBGeometry provides
-``PointCloudBVH`` (:file:`Source/EBGeometry_PointCloudBVH.hpp`). It specializes the :ref:`Chap:ImplemBVH`
-machinery: it is a subclass of ``BVH::PackedBVH`` with two things the general path does not offer -- a
-much cheaper **index-based build**, and **turnkey query methods** that hide ``pruneTraverse()``
-entirely.
-
-It is built directly from a raw cloud -- particle positions plus a parallel array of user metadata --
-by partitioning an index permutation in place with a longest-axis midpoint split and packing the
-``PointAoSoA`` leaves inline (no intermediate primitive list, no ``shared_ptr``, no separate packing
-pass), which is several times faster to build than a full Surface-Area-Heuristic tree and, for
-near-uniform clouds, just as tight to query.
-
-Queries return the matched particle's cloud index (and squared distance); the user metadata is
-reachable through ``metadata()``. ``closestPoint`` / ``closestPoints`` answer an arbitrary external
-point, while ``nearestNeighbor`` / ``nearestNeighbors`` (and the batch ``allNearestNeighbors``) answer
-a particle already in the cloud and additionally seed the search bound from that particle's own leaf
--- a strictly cheaper search an external point cannot use.
-
-Each accelerated query also has an ``O(N)`` brute-force counterpart -- ``closestPointBruteForce`` /
-``closestPointsBruteForce`` / ``nearestNeighborBruteForce`` / ``nearestNeighborsBruteForce`` -- that
-answers the same question by a full linear scan. These are reference implementations for testing and
-debugging (verify an accelerated result against ground truth, or A/B-test a suspected tree/traversal
-bug against an unaccelerated path) and are not meant for production queries.
-
-See the `PointCloudBVH doxygen page
-<doxygen/html/classEBGeometry_1_1PointCloudBVH.html>`__ for the full interface, and
-``Examples/ClosestPointBVH`` / ``Examples/NearestNeighborBVH`` for worked usage.
+A *point cloud* is simply a set of points in space, with no connectivity, orientation, or
+inside/outside notion. Unlike a surface mesh, there is therefore no *signed* distance to a point
+cloud -- only unsigned distances between points. The operations of interest are proximity queries:
+
+* **Closest point.** Given an arbitrary query point, find the cloud point nearest to it.
+* **k nearest.** Find the :math:`k` closest cloud points to a query, in ascending order of distance.
+* **k-nearest-neighbor graph.** For *every* point in the cloud, find its :math:`k` nearest *other*
+  points -- the classic all-nearest-neighbors problem.
+
+Answering one query by scanning all :math:`N` points is :math:`\mathcal{O}(N)`, so the
+all-nearest-neighbors graph is :math:`\mathcal{O}(N^2)` -- infeasible for large clouds. As with mesh
+distance queries, the cost is reduced by a spatial acceleration structure that lets a query rule out
+the vast majority of points without ever measuring the distance to them. EBGeometry offers two such
+structures, built on two different ideas.
+
+Hierarchical partitioning
+--------------------------
+
+The first idea is to build a tree over the points, recursively subdividing them into spatially tight
+groups -- exactly the bounding volume hierarchy described in :ref:`Chap:BVH`, with the points (or
+small groups of them) as the primitives. A query then descends the tree, at each node visiting the
+nearer child first and *pruning* any subtree whose bounding volume is already farther than the best
+match found so far. Because the subdivision follows the points themselves, the tree adapts to the
+cloud's density -- it stays balanced whether the points are uniform, lie on a surface, or clump into
+clusters -- and a query touches only :math:`\mathcal{O}(\log N)` nodes on average.
+
+Uniform grid
+------------
+
+The second idea dispenses with the tree entirely and lays a single regular grid of fixed-size cells
+over the cloud, bucketing each point into the cell that contains it. A query starts in the cell
+containing the query point and searches outward one shell of cells at a time (Chebyshev radius
+0, 1, 2, ...). It can stop as soon as the best distance found so far is closer than the nearest edge
+of the next unvisited shell -- no closer point can lie beyond that shell, so the search terminates
+*exactly*, never missing a neighbor. If the cells are sized to hold about one point each, a query
+resolves in the first shell or two.
+
+The grid trades adaptivity for simplicity. A single global cell size fits a **near-uniform** cloud
+best -- there, both building the grid (a counting sort) and querying it are cheaper than a tree. For
+a **strongly clustered or multi-scale** cloud a single cell size is a poor compromise (too coarse
+where the cloud is dense, too fine where it is sparse), and the density-adaptive hierarchy is the
+better choice.
+
+Self queries
+------------
+
+The all-nearest-neighbors graph is a special case worth calling out: every query point is *itself* a
+member of the cloud. A point therefore already knows which cell or leaf it lives in, so its search
+can be *seeded* from that group -- giving a tight pruning bound immediately -- and must exclude the
+point itself (otherwise it would trivially find itself at distance zero). This makes a batch of self
+queries strictly cheaper than the same number of independent external queries.

Docs/Sphinx/source/index.rst

@@ -84,6 +84,7 @@ Concepts
    GeometryRepresentations.rst
    DCEL.rst
    BVH.rst
+   PointCloud.rst
    Octree.rst
    CSG.rst
    Triangles.rst
@@ -101,7 +102,7 @@ Implementation
    ImplemCSG.rst
    ImplemDCEL.rst
    ImplemBVH.rst
-   PointCloud.rst
+   ImplemPointCloud.rst
    ImplemOctree.rst
    Parsers.rst
    SIMDClasses.rst

EBGeometry.hpp

@@ -28,6 +28,7 @@
 #include "Source/EBGeometry_Parser.hpp"
 #include "Source/EBGeometry_PointAoSoA.hpp"
 #include "Source/EBGeometry_PointCloudBVH.hpp"
+#include "Source/EBGeometry_PointCloudHashGrid.hpp"
 #include "Source/EBGeometry_PointSoA.hpp"
 #include "Source/EBGeometry_Polygon2D.hpp"
 #include "Source/EBGeometry_Random.hpp"

Examples/CMakeLists.txt

@@ -34,6 +34,7 @@ ebgeometry_add_example(ClosestPointBVH)
 ebgeometry_add_example(CSGUnion)
 ebgeometry_add_example(MeshSDF)
 ebgeometry_add_example(NearestNeighborBVH)
+ebgeometry_add_example(NearestNeighborHashGrid)
 ebgeometry_add_example(NestedBVH)
 ebgeometry_add_example(OctreeBoundingVolume)
 ebgeometry_add_example(PackedSpheres)

Examples/NearestNeighborHashGrid/CMakeLists.txt

@@ -0,0 +1,62 @@
+# SPDX-FileCopyrightText: 2026 Robert Marskar <robert.marskar@sintef.no>
+#
+# SPDX-License-Identifier: GPL-3.0-or-later
+
+# Standalone build file for the NearestNeighborHashGrid example. It can be built on its own
+# (e.g. copied out as a template for a new project) or as part of the top-level
+# EBGeometry build.
+cmake_minimum_required(VERSION 3.16)
+project(NearestNeighborHashGrid CXX)
+
+# Locate the header-only library. EBGEOMETRY_HOME is the directory that contains
+# EBGeometry.hpp (the repository root); it defaults to the in-tree location two
+# levels up. When built from the top-level EBGeometry project the
+# EBGeometry::EBGeometry target already exists and is used as-is.
+if(NOT TARGET EBGeometry::EBGeometry)
+  set(EBGEOMETRY_HOME "${CMAKE_CURRENT_SOURCE_DIR}/../.." CACHE PATH
+      "Path to the directory containing EBGeometry.hpp")
+  add_library(EBGeometry::EBGeometry INTERFACE IMPORTED)
+  target_include_directories(EBGeometry::EBGeometry INTERFACE "${EBGEOMETRY_HOME}")
+  target_compile_features(EBGeometry::EBGeometry INTERFACE cxx_std_17)
+endif()
+
+add_executable(NearestNeighborHashGrid main.cpp)
+
+# Give the binary a .ex extension so it's covered by the repo-wide *.ex .gitignore rule,
+# without needing a per-example ignore entry. Always place it directly in this
+# source directory (not wherever -B/--preset pointed the build tree), so all three
+# documented build methods (direct compile, GNU make, CMake) produce the binary in
+# the same place: ./NearestNeighborHashGrid.ex.
+set_target_properties(NearestNeighborHashGrid PROPERTIES
+  SUFFIX ".ex"
+  RUNTIME_OUTPUT_DIRECTORY "${CMAKE_CURRENT_SOURCE_DIR}")
+target_link_libraries(NearestNeighborHashGrid PRIVATE EBGeometry::EBGeometry)
+target_compile_options(NearestNeighborHashGrid PRIVATE -Wall -Wextra -Wno-unused-parameter)
+
+# Examples are meant to be fast to run, not to match the enclosing project's build
+# type (which may be an unoptimised Debug build, e.g. under the top-level project's
+# debug preset). Always build optimized for the local machine's ISA.
+target_compile_options(NearestNeighborHashGrid PRIVATE -O3 -march=native)
+
+# ---------------------------------------------------------------------------
+# Floating-point precision for this example.
+#
+# Cache variable, overridable at configure time with
+# -DEBGEOMETRY_PRECISION=float. It is handed to the compiler as
+# EBGEOMETRY_PRECISION and picked up in main.cpp by `using T = EBGEOMETRY_PRECISION`.
+# ---------------------------------------------------------------------------
+set(EBGEOMETRY_PRECISION double CACHE STRING "Floating-point precision (float or double)")
+set_property(CACHE EBGEOMETRY_PRECISION PROPERTY STRINGS float double)
+target_compile_definitions(NearestNeighborHashGrid PRIVATE EBGEOMETRY_PRECISION=${EBGEOMETRY_PRECISION})
+
+# ---------------------------------------------------------------------------
+# EBGEOMETRY_EXPECT() runtime assertions for this example.
+#
+# OFF by default. Enable at configure time with -DEBGEOMETRY_ENABLE_ASSERTIONS=ON.
+# When built from the top-level EBGeometry project, this option already exists
+# (defined by the root CMakeLists.txt) and its value is inherited as-is.
+# ---------------------------------------------------------------------------
+option(EBGEOMETRY_ENABLE_ASSERTIONS "Enable EBGEOMETRY_EXPECT() runtime assertion checks" OFF)
+if(EBGEOMETRY_ENABLE_ASSERTIONS)
+  target_compile_definitions(NearestNeighborHashGrid PRIVATE EBGEOMETRY_ENABLE_ASSERTIONS)
+endif()

Examples/NearestNeighborHashGrid/GNUmakefile

@@ -0,0 +1,32 @@
+# GNU makefile for building this EBGeometry example on its own.
+#
+# EBGEOMETRY_HOME must point at the directory that contains EBGeometry.hpp (the
+# EBGeometry repository root); it defaults to the in-tree location two levels up.
+# PRECISION selects the floating-point type (float or double). ASSERTIONS
+# (ON/OFF) toggles EBGEOMETRY_EXPECT() runtime assertion checks. Override any
+# of these on the command line, e.g.:
+#
+#   make PRECISION=float ASSERTIONS=ON EBGEOMETRY_HOME=/path/to/EBGeometry
+
+EBGEOMETRY_HOME ?= ../..
+PRECISION       ?= double
+ASSERTIONS      ?= OFF
+
+CXX      ?= g++
+CXXFLAGS ?= -std=c++17 -O3 -march=native -Wall -Wextra
+
+ifeq ($(ASSERTIONS),ON)
+  CXXFLAGS += -DEBGEOMETRY_ENABLE_ASSERTIONS
+endif
+
+TARGET := NearestNeighborHashGrid.ex
+
+$(TARGET): main.cpp
+	$(CXX) $(CXXFLAGS) -I$(EBGEOMETRY_HOME) -DEBGEOMETRY_PRECISION=$(PRECISION) $< -o $@
+
+.PHONY: run clean
+run: $(TARGET)
+	./$(TARGET)
+
+clean:
+	$(RM) $(TARGET)

Examples/NearestNeighborHashGrid/README.md

@@ -0,0 +1,59 @@
+Examples/NearestNeighborHashGrid
+--------------------------------
+
+Nearest-neighbor search over a point cloud using the turnkey
+[`PointCloudHashGrid`](https://rmrsk.github.io/EBGeometry/doxygen/html/classEBGeometry_1_1PointCloudHashGrid.html)
+class -- the uniform-grid counterpart to
+[`Examples/NearestNeighborBVH`](../NearestNeighborBVH/README.md), with the **same interface**, so
+switching between the tree and the grid is a one-line type change.
+
+`allNearestNeighbors(k)` computes the `k` nearest neighbors of every point in one batched pass; a
+self query excludes the point itself, so the neighbors returned are the nearest *other* points. The
+grid suits **near-uniform** clouds; for clustered or multi-scale clouds the density-adaptive
+`PointCloudBVH` is the better choice.
+
+Building
+--------
+
+This example is standalone and can be built in three ways. Each needs the path
+to the EBGeometry root -- the directory that contains `EBGeometry.hpp` -- which
+is two levels up from this folder (`../..`) when building in place. See
+[Direct compilation](https://rmrsk.github.io/EBGeometry/BuildingDirectCompile.html),
+[Building with GNU Make](https://rmrsk.github.io/EBGeometry/BuildingGNUMake.html), and
+[Building with CMake](https://rmrsk.github.io/EBGeometry/BuildingCMake.html) in the user
+documentation for more detail on each approach.
+
+**CMake**
+
+    cmake -S . -B build
+    cmake --build build
+
+The binary is `NearestNeighborHashGrid.ex`, in this same directory (same as the other two methods
+below). Build in single precision, against a library in a different location, or with
+`EBGEOMETRY_EXPECT()` runtime assertions enabled, with cache variables:
+
+    cmake -S . -B build -DEBGEOMETRY_PRECISION=float -DEBGEOMETRY_HOME=/path/to/EBGeometry -DEBGEOMETRY_ENABLE_ASSERTIONS=ON
+
+**GNU make**
+
+    make
+
+This produces `./NearestNeighborHashGrid.ex`. Override the defaults on the command line:
+
+    make PRECISION=float EBGEOMETRY_HOME=/path/to/EBGeometry ASSERTIONS=ON
+
+**Directly with a compiler**
+
+    g++ -std=c++17 -O3 -march=native -I../.. main.cpp -o NearestNeighborHashGrid.ex
+
+Add `-DEBGEOMETRY_PRECISION=float` for single precision, or `-DEBGEOMETRY_ENABLE_ASSERTIONS` to enable
+`EBGEOMETRY_EXPECT()` runtime assertion checks.
+
+Running
+-------
+
+    ./NearestNeighborHashGrid.ex
+
+Takes no arguments. It prints the build and per-point query times and one worked `nearestNeighbor()`
+result. The neighbor graph is checked against a brute-force scan when built with
+`-DEBGEOMETRY_ENABLE_ASSERTIONS`.