# Exploring Spatial Indexing for Accelerated Feature Retrieval in HPC

Margaret Lawson  The University of Illinois at Urbana-Champaign Sandia
National LaboratoriesAlbuquerqueNew Mexico<EMAIL_ADDRESS>,  William
Gropp  The University of Illinois at Urbana-Champaign<EMAIL_ADDRESS>and
Jay Lofstead  Sandia National LaboratoriesAlbuquerqueNew Mexico
<EMAIL_ADDRESS>
(2021)

###### Abstract.

Despite the critical role that range queries play in analysis and
visualization for HPC applications, there has been no comprehensive analysis
of indices that are designed to accelerate range queries and the extent to
which they are viable in an HPC setting. In this state of the practice paper
we present the first such evaluation, examining 20 open-source C and C++
libraries that support range queries. Contributions of this paper include
answering the following questions: which of the implementations are viable in
an HPC setting, how do these libraries compare in terms of build time, query
time, memory usage, and scalability, what are other trade-offs between these
implementations, is there a single overall best solution, and when does a
brute force solution offer the best performance? We also share key insights
learned during this process that can assist both HPC application scientists
and spatial index developers.

geometric range searching, spatial indexing, k-d tree, R-tree, octree

††copyright: iw3c2w3††journalyear: 2021

##  1\. Introduction

During analysis and visualization, HPC application scientists often need to
extract particular spatial subsets of their data. Scientists may need to
extract these subsets for sub-sampling (Lofstead et al., 2011), to perform
areal interpolation (Papadomanolakis et al., 2006a; Hinterberger et al.,
1994), or to perform spatial correlations (Lu et al., 2017). Scientists may
also extract spatial subsets that correspond to regions or features of
interest. This can allow scientists to concentrate their analysis and
visualization on more interesting data or to use different analysis routines
depending on the type of region. Examples of regions that may be extracted
include application boundaries and surfaces (Papadomanolakis et al., 2006b),
regions of the brain in neuroscience imagery (Beyer et al., 2013; Tauheed et
al., 2012), counties or other regional boundaries in satellite imagery (Beynon
et al., 2000), regions with different flow types in a fluid dynamics
simulation (Barone et al., 2021), and regions with different chemical
properties in a combustion simulation (Gosink et al., 2007). Examples of
spatial features of interest that may be extracted include supernovae in
astronomy (Betoule et al., 2013), blood vessels, lesions, and areas of
inflammation in pathology imagery (Aji et al., 2012; Liang et al., 2016,
2017), tissue regions with particular cell types in microscopy datasets
(Beynon et al., 2000), shocks in a fluid simulation (Beynon et al., 2000),
flame fronts, extinction events and vortices in a combustion simulation
(Gosink et al., 2007; Gosink et al., 2010), and tropical storms in a climate
simulation (Chiu et al., 2015). These regions and features may be identified
using techniques such as image segmentation (Pecha, 2016; Pradipkumar and
Rajamenakshi, 2016), machine-learning based approaches (Graff et al., 2014;
Liu et al., 2016; Tzeng and Ma, 2005; Cubuk et al., 2015), or other analysis
techniques (Bremer et al., 2010; Stockinger et al., 2005; Post et al., 2003;
Wang et al., 2008; Livnat et al., 1996; Yu et al., 2010). Extracting spatial
subsets is thus a critical component of many different types of analysis and
visualization routines for a wide range of HPC applications.

Extracting a particular spatial subset is trivial for regular meshes (meshes
in which the cells are congruent parallelotopes) where a single, simple
function can be used to map from a spatial coordinate to a data array offset.
In the simplest case (a Cartesian mesh), the mesh coordinates will be the same
as the array indices. In contrast, for unstructured meshes (meshes without
implied neighborhood connectivity), this operation is very computationally
intensive. Without any form of index, scientists will have to perform a linear
(brute force) search over all mesh points or elements for each region they
want to extract. With petascale simulations already using meshes that can
range from tens of billions (Rasquin et al., 2014; Wang et al., 2017; Houzeaux
et al., 2011; Zhao et al., 2020) to trillions (Godenschwager et al., 2013) of
elements and with the impending arrival of exascale machines (of Energy,
2019), this kind of linear search, which may need to be performed thousands or
millions of times during analysis and visualization, represents a significant
inefficiency that can delay the path to scientific discovery. A lot of
research has been done to develop spatial indices that can accelerate this
type of search (known as a range query) to help scientists quickly identify
what mesh points or elements fall within a given region. However, there has
been no comprehensive analysis of spatial index implementations that are
designed to accelerate range queries and the extent to which they are viable
in an HPC setting. To be viable in an HPC setting, a spatial indexing library
must have good performance for building the index and performing queries, good
scalability, and moderate memory overheads. It is important for the index to
have moderate memory overheads since the index will provide the best
performance when stored in memory, many applications require significant
amounts of memory (for MPI, for the mesh, for variable data, etc.) and nodes
have limited amounts of memory.

In this state of the practice paper, we present the first comprehensive review
of spatial index implementations that support range queries. We also provide
scientists with the information they need to decide when to use a spatial
index for range querying and which spatial index implementation will be best
suited to their use case. This information can help scientists accelerate
range queries by several orders of magnitude and thereby greatly accelerate
their analysis routines and the discovery process. Contributions of this paper
include presenting a thorough analysis of 20 free, open-source C and C++
libraries that support range queries with expected sub-linear query times,
sharing key insights learned during this process, and answering the following
questions:

  1. (1)

Which of the implementations are viable in an HPC setting?

  2. (2)

How do these libraries compare in terms of build time, query time, and memory
usage at different scales?

  3. (3)

What are other factors in deciding which library to use?

  4. (4)

Is there a single overall best solution?

  5. (5)

When does a brute force solution offer the best performance?

The rest of this paper is organized as follows. Section 2 discusses related
work. Section 3 provides an overview of the evaluated libraries. Section 4
presents the evaluation and results. Finally, Section 5 provides discussion,
offers additional insights, and presents areas for future work.

##  2\. Related Work

Range queries have been extensively studied in many different branches of
computer science including privacy and security (Zhu et al., 2015; Xu et al.,
2018), mobile computing (Cong and Jensen, 2016), networks and other graph-
based areas (Demirbas and Ferhatosmanoglu, 2003; Soheili et al., 2005; Wan et
al., 2019), cloud computing (Yu et al., 2015; Tang et al., 2016; Eldawy and
Mokbel, 2015; Zhong et al., 2012), and databases (Papadias et al., 2003;
Shekhar et al., 1999; Egenhofer, 1994; Orenstein, 1986; Güting, 1994). Below
we provide an overview of two categories of related work that we find to be
most related to the work presented here: theoretical research done on
geometric range searching and search structures and research that either
directly implements structures that support range queries or evaluates these
structures.

###  2.1. Theory: Geometric Range Searching and Geometric Search Structures

Efficient geometric range searching has been the focus of a significant amount
of research in the theoretical branch of computer science (Matoušek, 1994;
Agarwal, 1996; Willard, 1985; Chazelle, 1990; Hadjieleftheriou et al., 2005;
Jagadish, 1990; Pagel et al., 1993; Bentley and Maurer, 1980; Overmars, 1988).
This work has invented a number of geometric search structures, and here we
present an overview of the three that are most commonly supported by the
libraries we evaluate in this paper: octrees, k-d trees, and R-trees.

####  2.1.1. Octrees

An octree (Jackins and Tanimoto, 1980) is a tree in which each internal node
has exactly eight children and which is not necessarily balanced. An octree
can either be used to recursively partition a 3D space into eight octants
(region-based octrees) or to partition a set of 3D points based on their
coordinates (point-based octrees). Internally, octrees typically store point
data only in the leaf nodes for region-based octrees and typically store point
data both at leaf nodes and internal nodes for point-based octrees (Formaggia,
1999).

####  2.1.2. K-d trees

K-d trees (Bentley, 1975) are binary trees that store k-dimensional data and
represent a recursive subdivision of this data using $(k-1)$-dimensional
hyperplanes. Each (non-leaf) level in the tree corresponds to one of the
k-dimensions, and each internal node at that level represents a splitting
hyperplane for that dimension. K-trees are sensitive to the order in which
points are inserted and are not necessarily balanced. Typically the data will
be split before the median point along the longest side of the node (the
sliding midpoint rule), but other splitting strategies are possible.

####  2.1.3. R-trees

R-trees (Guttman, 1984) are balanced trees that contain a hierarchy of
d-dimensional boxes, where each child node is contained in the box represented
by its parent. Depending on the implementation, these partitions may overlap.
The leaves either store a d-dimensional point or the d-dimensional minimum
bounding box of the objects stored in the leaf. These objects may be a set of
points or other shapes. R-trees have a minimum and maximum branching factor
for all internal nodes (apart from the root), and there are algorithms
designed to build optimal R-trees using bulk (static) construction
(Leutenegger et al., 1997; García R et al., 1998). R*-trees are an R-tree
variant that is designed to minimize coverage and spatial overlap between
internal nodes.

###  2.2. Implementations and Evaluations of Structures Supporting Range
Queries

Several recent research efforts have implemented data structures that support
range queries. However, to our knowledge, none of these implementations are
open-source. This includes the development of a query processing technique for
unstructured tetrahedral meshes in PostgreSQL (Papadomanolakis et al., 2006a),
intersection queries in PostGIS (de Oliveira et al., 2020), point, range and
polygon queries in Microsoft’s SQL Server (Perlman et al., 2010), and parallel
range queries using GPUs (Nguyen and Rhodes, 2016; Nouri and Tu, 2018; Kim et
al., 2013; Prasad et al., 2015; Luo et al., 2012; You et al., 2013). One
project performed a comparison of nearest neighbor search implementations and
included some evaluation of range queries (Elseberg et al., 2012). However,
this work only considers six of the libraries we evaluate here, collects range
query results only for artificial datasets containing 60,000 points, and was
performed almost a decade ago.

###  2.3. Summary

Our work differs from related work in that it focuses on spatial index
implementations that are open-source and which support range queries. Most
related work is focused on structures that are theoretical rather than
implemented, structures that do not offer range queries, or structures that
offer range queries but are not open-source. This work also differs from
related work in that it evaluates structures for use in HPC and therefore uses
significantly larger datasets than those used in related work and offers an
analysis of memory requirements, and strong and weak scaling.

##  3\. Spatial Indexing Libraries Overview

We evaluate 20 libraries that support range queries. To be included in this
study, a library has to: provide bindings for C or C++, offer a free, open-
source version, and offer spatial indexing that supports range querying (with
expected sub-linear query times). We limit our search to libraries with C or
C++ bindings since these are the languages most commonly used in HPC for
analysis and visualization tasks (Amaral et al., 2020; Laguna et al., 2019).
Only one library, LEDA (Mehlhorn and Näher, 1999), is excluded for failing to
meet the free, open-source requirement. LEDA offers a free version, but range
queries are only offered in its paid version, which costs 3000 euros for pure
research efforts. We leave for future work evaluating libraries that offer
other language bindings.

Table 1 contains basic information about the libraries. In the evaluation, we
perform box queries since this is the type of query used most commonly in
analysis and visualization. For libraries that do not support box queries, we
use a sphere (or radius) that is large enough to encompass the box, and then
perform a filtering step to ensure exact matches are returned. This will
result in a slight performance penalty for these libraries.

Table 1. Basic info. about the libraries evaluated in this paper Library | Ver. |  | Tree
---
type
| Search
---
type
| Data
---
types
| Can set
---
leaf size
3DTK (Nüchter and Lingemann, 2011) | - | k-d tree | Box | Points | Yes
ALGLIB (Bochkanov and Bystritsky, 2011) | 3.17.0 | k-d tree | Box | Points | No
ANN (Arya and Mount, 1998) | 1.1.2 | k-d tree |  | KNN w.
---
radius
Points | Yes
ANN (Arya and Mount, 1998) | 1.1.2 | BD-tree |  | KNN w.
---
radius
Points | Yes
Boost (Gehrels et al., 2020) | 1.74.0 | R-tree | Box | Points, boxes | Yes
CGAL (Fabri and Pion, 2009) | 5.2.0 | k-d tree | Box | Points | Yes
CGAL (Fabri and Pion, 2009) | 5.2.0 | Box | Box | Points | No
CGAL (Fabri and Pion, 2009) | 5.2.0 | Segment | Box | Boxes | No
CGAL (Fabri and Pion, 2009) | 5.2.0 | R-tree | Box | Points, shapes | No
FLANN (Muja and Lowe, 2009) | 1.9.1 | k-d tree | Sphere | Points | Yes
KDTREE(1) (Carvalho, 2020) | - | k-d tree | Sphere | Points | No
KDTREE2 (Kennel, 2004) | - | k-d tree | Sphere | Points | Yes
KDTREE(3) (Anon, 2020) | - | k-d tree | Sphere | Points | Yes
KDTREE(4) (Tsiombikas, 2019) | 0.5.7 | k-d tree | Sphere | Points | No
libkdtree++ (Krafft, 2020) | - | k-d tree | Sphere | Points | No
libkdtree(2) (Dietrich, 2013) | - | k-d tree | Box | Points | No
libnabo (Elseberg et al., 2012) | 1.0.7 | k-d tree |  | KNN w.
---
radius
Points | Yes
libspatialindex (Hadjieleftheriou, 2015) | 1.9.3 | R*-tree | Box | Boxes | Yes
nanoflann (Blanco and Rai, 2014) | 1.3.2 | k-d tree | Sphere | Points | Yes
Octree (Behley et al., 2015) | - | octree | Sphere | Points | Yes
PCL (Rusu and Cousins, 2011) | 1.11.1 | k-d tree | Sphere | Points | No
PCL (Rusu and Cousins, 2011) | 1.11.1 | octree | Box | Points | Yes
PicoTree (Broere, 2021) | 0.5.2 | k-d tree | Box | Points | Yes
R-tree (Green and Douglas, 2004) | - | R-tree | Box | Boxes | Yes
Spatial (Bougerel, 2017) | 2.1.8 | k-d tree | Box | Points, boxes | No

Note that some of the libraries share the same name or similar names (e.g.,
KDTREE and KDTREE2, libkdtree++ and libkdtree). We have added a number (shown
in parentheses) after some of these library names to help disambiguate them.
For the rest of the paper, we will refer to these libraries using the name
appended with the number (without the parentheses).

Insight: few libraries can store mesh elements (shapes/boxes). Since many
variables are calculated per mesh element, this is a significant limitation.

For more information about these libraries, readers can look at the Git
repository for this project (Lawson, 2021).

##  4\. Evaluation and Results

This section presents an overview of the evaluation setup including the
hardware and software, datasets, and scales used in testing. We also present
the evaluation results for storage and querying of mesh points and mesh
elements at two different scales, and present an analysis of the scalability
of the libraries. HPC application data is typically collected for each mesh
node or mesh element. This is why we perform an evaluation for both mesh
points and elements.

###  4.1. Experimental Setup

All experiments are run on the Vortex machine at Sandia National Labs. Vortex
uses RHEL7, and has the following per node: 318 GB high-bandwidth DRAM and
dual socket IBM POWER9 CPUs with 22 cores/socket and 4 hardware threads per
core (176 total). We use GCC 10.2.0 to compile all of the libraries apart from
PCL, which we can only build with GCC 8.3.1 on Vortex. We use four different
evaluation setups, which are described in Table 2. All experiments use the
same mesh, which was created for a typical (Domino et al., 2019) turbulent
low-Mach study that includes solution verification. The turbulent low-Mach
study uses meshes that range from around 150 million to 2 billion nodes and
elements. The mesh we use is unstructured and has 152.7 million nodes and
152.1 million hexahedral elements. The coordinates are 8 byte doubles. It
should be noted that there is no standard in HPC for the number of mesh points
or elements assigned per process during domain decomposition. The number will
depend on a number of factors such as the computational intensity of the
application, the memory per node, and the number of compute hours the
scientist has access to. However, using between $10,000$ and $100,000$ nodes
or elements per process during a simulation is common (Gahvari et al., 2011;
Hsu et al., 2011; Quintanas-Corominas et al., 2018; Gabriel et al., 2010; Feki
and Gabriel, 2008; Zaspel and Griebel, 2011; Colli Franzone et al., 2018;
Garcia-Gasulla et al., 2020; Landge et al., 2014). Then, during analysis and
visualization, scientists will often use a much smaller number of processors
(10% or fewer) (Marrinan et al., 2019; Peña-Monferrer et al., 2021; Cummings
et al., 2008) thus resulting in upwards of $100,000$ to $1,000,000$ mesh nodes
or elements per process. This is further evidenced by the fact that analysis
and visualization clusters often have far fewer cores and nodes than compute
clusters (Simakov et al., 2018; Salim et al., 2019). We therefore perform
“small” scale evaluation with approximately $1,000,000$ mesh nodes or elements
per process to reflect this common range, and “large” scale evaluation with
approximately $10,000,000$ mesh nodes or elements per process to reflect more
extreme use cases. We use OpenMPI 4.1.0 with one hardware thread per (MPI)
process. Apart from the weak scaling results, all evaluations use a single
node. But, as the weak scaling results demonstrate in Section 4.3.2, these
results should hold at any scale if the number of processes per node is held
constant. This is because each process performs entirely independent work
without any need for message passing or coordination, and the indices are also
kept entirely in-memory (eliminating the possibility of file system
contention).

Table 2. Experiment setups | Mesh points | Mesh elements
---|---|---
160 procs. (small scale) | 1,227,411 / proc. | 1,087,807 / proc.
16 procs. (large scale) | 10,163,37 / proc. | 10,878,078 / proc.

###  4.2. Performance Evaluations

Each performance evaluation consists of two basic parts. First, each process
creates an in-memory instance of the data structure being evaluated and
inserts the mesh coordinates or elements for its subdomain (assigned portion
of the simulation space). Second, each process performs a set of queries at
four different sizes, which we refer to as extra-small, small, medium and
large. These queries return approximately 0.001%, 0.1% 1%, and 10% of the
data. Each process performs 10,000 queries that are extra-small, 10,000 that
are small, 10,000 that are medium, and 1,000 that are large. This allows us to
evaluate how the spatial indices perform for retrieving features or regions of
interest that range from very small to very large. This testing workload is
designed to reflect the real-world scenario described and cited in the
introduction: when a scientist performing analysis or visualization with a
complex mesh type needs to extract data for particular regions of interest. In
the large scale evaluation, we evaluate all libraries that perform within a
factor of 10$\times$ of the best performing library for each of the query
sizes. This helps us to identify which libraries perform “best” across
different query sizes and at different scales. For all libraries that allow
the user to set a maximum leaf size, we test leaf sizes of 20 and 200. For
space reasons, we only present the best result for each library. The full set
of results can be found in the Git repository for this project (Lawson, 2021).
Several libraries also offer different insertion algorithms or splitting
strategies (for k-d or BD trees). However, we did not find any significant
performance differences between these algorithms and have therefore omitted
these results to save space. One exception is that we use Boost’s STR packing
algorithm for all of the Boost tests since it results in dramatically improved
results for data that can be bulk loaded (like mesh data). For libraries that
offer both static and dynamic trees, we use static trees since the mesh we use
in evaluation is static. We perform three runs for each configuration and
average the results. After completing the timing results, we also collect
memory results using the Massif heap profiler (val, 2020) from Valgrind 3.16.1
(Nethercote and Seward, 2007). For these memory results, we present data only
on the construction of the tree since memory increases during querying are
transient. We exclude memory allocated as a result of reading in the mesh
nodes or elements from the mesh file. We run these tests on a single process
that has the amount of data closest to the average (within 3% of the global
average).

In addition to the libraries, we evaluate a brute force solution. In the build
phase, the solution copies the 2D vector of mesh points or elements using the
assignment operator. In reality, no such copy would be required, we merely
include this result as a point of comparison. For the queries, the brute force
solution performs a linear search checking for each point or element if it
intersects the query range. Evaluating the brute force solution allows us to
determine how much of a performance advantage the spatial indices offer and to
evaluate whether this advantage is sufficient to justify the memory
requirements of the index.

In each of the result tables below, for the build throughput (insertions per
second), and the query throughput (queries per second), we classify a
library’s performance as fast (green cell with bolded text) if is within
10$\times$ of the best performing library for build time or for the given
query size. We classify the performance as moderate (yellow with plain text)
if it is within 100$\times$. All other results are classified as slow (red
with underlined text). For memory usage, we classify any configuration that
uses 2$\times$ or less of the raw memory size as having low memory
requirements (green with bolded text), any configuration that uses less than
5$\times$ the raw memory size as having moderate memory requirements (yellow
with plain text), and all other configurations as having high memory
requirements (red with underlined text). In all cases, the raw data size can
be seen by looking at the memory requirements for the brute force solution,
which makes a copy of the mesh points or elements. In the tables, memory
(abbreviated mem.) indicates the memory used by the tree, while peak memory
indicates the maximum memory used while creating the tree.

Insight: some libraries temporarily use large amounts of additional memory
making them a poor choice for significantly memory constrained environments.

####  4.2.1. Small Scale Mesh Points

Each query of extra-small, small, medium, and large size retrieves an average
of 0.00104%, 0.100%, 1.03% and 8.85% of the process’s assigned mesh points.
Even at this small scale, a few libraries run into errors and are thus not
included in the results. First, although FLANN and PCL can use GPUs, with 4
NVIDIA Tesla V100 GPUs with a total of 64 GB RAM, they run into out of memory
errors. In addition, since they are the only two libraries that support GPUs
it is not possible to perform a true cross-library comparison of this
capability. Next, the CGAL range tree experiences out of memory errors.
Finally, none of the ANN configurations complete within 24 hours. We therefore
include results only for the query sizes that complete in these 24 hours.

Insight: only two of the libraries can utilize GPUs and they quickly exhaust
the GPUs memory.

Table 3. Small scale results for mesh points (160 processes, $\approx$1 million points/process). Even at small scale, some libraries perform orders of magnitude better on writing, querying, memory usage. Green is good, yellow is moderate, red is poor222A complete explanation of the classification scheme can be found in Section 4.2.. Configuration |  | Leaf
---
size
| Insertions/sec
---
(millions)
| X-small
---
queries/sec
| Small
---
queries/sec
| Medium
---
queries/sec
| Large
---
queries/sec
| Mem.
---
(MB)
| Peak mem.
---
(MB)
Brute force | - | 5.51 |  29.3 |  29.3 | 29 | 27.7 | 29.2 | 29.2
3DTK | 20 |  0.00383 |  10.2 |  9.32 | 7.41 | 3.97 |  612 |  780
ALGLIB | 8 | 1.11 | 36100 | 2360 | 415 | 69.6 |  152 |  384
ANN k-d tree | 200 | 1.16 |  35.4 | - | - | - | 54.6 | 54.6
ANN BD-tree | 200 | 0.704 |  36.1 | - | - | - | 44.9 | 54.6
Boost | 200 | 2.87 | 47900 | 3660 | 560 | 61.1 | 79.1 |  185
CGAL k-d tree | 200 | 33.2 | 9690 | 4060 | 625 | 62.9 | 77.2 | 111
FLANN | 200 | 1.98 |  78.6 | 55.5 | 25.4 | 7.64 | 54.8 | 54.8
KDTREE1 | 1 | 0.388 |  99.5 |  4.9 |  0.901 |  0.18 | 127 |  292
KDTREE2 | 200 | 2.83 | 10400 | 571 | 97.6 | 19.2 | 45.8 | 45.8
KDTREE3 | 200 | 1.65 | 5370 | 228 | 35.5 | 6.33 | 59.6 | 127
KDTREE4 | 1 |  0.297 | 563 | 65.3 | 15.3 | 3.55 | 115 | 115
libkdtree++ | 1 |  0.11 | 1120 | 95.8 | 20.2 | 4.49 |  272 |  272
libkdtree2 | 1 |  0.012 | 24000 | 994 | 109 | 12.5 |  285 |  285
libnabo tree | 200 | 1.5 |  69.1 | 64.5 | 33 | 8.63 | 73.4 | 94.8
libspatialindex | 20 |  0.00283 | 4290 | 330 | 40.9 | 4.97 | 104 |  196
nanoflann | 200 | 0.43 |  83.5 | 59 | 27.4 | 8.45 | 20.7 | 20.7
nanoflann (duplicated data) | 200 | 1.09 |  89 | 61.5 | 28.3 | 8.65 | 69.4 | 79.1
Octree | 200 | 4.92 | 8180 | 526 | 98.2 | 20.9 | 29.5 | 29.5
PCL k-d tree | 15 | 0.708 | 5600 | 256 | 40.4 |  0.965 | 51.5 | 60.5
PCL octree | 20 | 4.61 |  140 | 112 | 77.1 | 4.35 | 51.5 | 60.5
PicoTree | 20 | 2.12 | 23500 | 1620 | 311 | 57.5 | 105 | 111
R-tree | 20 |  0.183 | 42300 | 2720 | 447 | 68.9 | 123 | 123
Spatial | 1 | 0.451 | 24000 | 330 | 32.3 | 3.6 | 108 |  212

The results can be seen in Table 2. A few things should be noted about these
results. First, as indicated in the table, the libnabo results are for the
tree (heap) with $O(\log n)$ query times rather than the $O(n)$ linear vector
heap. Second, we use KDTREE2’s option to rearrange the data for better
performance. We found this offers substantially better query performance with
only a slight decrease in build throughput. Third, nanoflann allows you to
create an adaptor for an existing vector so that the storage is not
duplicated. We show the results for both no duplicated data (“nanoflann”) and
duplicate data (“nanoflann duplicated data”). Finally, three libraries,
KDTREE2, libkdtree2, and Octree, have artificially reduced memory requirements
because they only support storing coordinate data as floats rather than
doubles (and thus offer lower precision).

From the results we can see that the CGAL k-d tree has by far the highest
insertion throughput, followed by Octree and the PCL octree. However, even one
of the slower libraries, R-tree with a leaf size of 20, builds in
approximately 7 seconds so it should be kept in mind that (apart from 3DTK and
libspatialindex) all of these trees build quite quickly. This highlights the
performance benefits of being able to parallelize construction (and querying)
across the processes through domain decomposition, rather than attempting to
use a single index for the entire mesh (although this does increase the total
memory consumption). Shockingly, CGAL is able to create the k-d tree
significantly faster than the brute force solution copies the 2D vector of
mesh coordinates. There are a few libraries that perform well across all query
sizes: ALGLIB, Boost, the CGAL k-d tree, KDTREE2, libkdtree2, Octree, PicoTree
and R-tree.

Insight: query performance is strongly dictated by the number of matching mesh
points and when a large number of points are returned, most libraries fail to
beat the brute force solution.

The poor performance on the large queries is a reflection of the fact that all
of these libraries are tree based and as such can offer a significant
performance advantage when a large portion of the search space can be pruned
and a significant performance penalty when large portions of the tree must be
traversed. It should also be noted that all of the libraries use far more
memory than the theoretical lower bound of linear storage. Almost all of the
libraries use between 2-5$\times$ the space needed for just the raw coordinate
points (29.2 MB), and several use significantly more than this. Only nanoflann
with no duplicated storage is able to use anywhere close to linear memory, but
it performs more than 10$\times$ slower than the best performing library for
all query sizes.

####  4.2.2. Small Scale Mesh Elements

Each query of extra-small, small, medium, and large size retrieves an average
of 0.00213%, 0.212%, 1.68% and 12.68% of the process’s assigned mesh elements.

Table 4. Small scale results for hexahedral mesh elements (160 processes, $\approx$1 million elements/process) for the libraries that support 2D or 3D shapes. We find moderately to significantly reduced performance compared to point storage (Table 3), better average performance compared to brute force, and reduced memory overheads. Green is good, yellow is moderate, red is poor444A complete explanation of the classification scheme can be found in Section 4.2.. Configuration |  | Leaf
---
size
| Insertions/sec
---
(millions)
| X-small
---
queries/sec
| Small
---
queries/sec
| Medium
---
queries/sec
| Large
---
queries/sec
| Mem.
---
(MB)
| Peak mem.
---
(MB)
Brute force | - | 2.23 |  15.9 |  15.8 | 15.7 | 14.9 | 47.0 | 47.0
Boost | 20 | 2.11 | 36100 | 2250 | 270 | 33.6 | 61.8 | 160
CGAL R-tree | 1 | 15.3 | 4170 | 1420 | 228 | 34.8 | 75.3 | 96.3
libspatialindex | 200 |  0.051 | 2140 | 213 | 32.8 | 4.64 | 72.4 | 227
R-tree | 200 | 0.18 | 33200 | 2460 | 419 | 67.6 | 102 | 102
Spatial | 1 | 0.393 | 1010 | 233 | 40.2 | 5.78 | 110 | 198

The results of this evaluation are shown in Table LABEL:tab:results-small-
bboxes. One library that we test, CGAL’s segment tree, runs into out of memory
errors, and is therefore not included in the results. In these results we
again find that a few of the libraries’ data structures perform well across
the board: Boost, the CGAL R-tree, and the R-tree library. In this case, the
libraries overall use (relatively) far less memory, with most using less than
2$\times$ the memory needed for the raw element coordinate data (47.0MB).
Unlike with the point results, here all of the best performing libraries
outperform brute force at query sizes (including the large queries) since the
brute force solution experiences an approximately 2$\times$ slowdown as it now
has to consider two points per vector element (the lower and upper corner of
the bounding box) rather than just one.

Table 5. Large scale results for mesh points (16 processes, $\approx$10 million points/process) for libraries that perform best at small scale. Despite the fact that all of these libraries perform well at small scale, significant differences emerge at large scale. Green is good, yellow is moderate, red is poor4. Configuration |  | Leaf
---
size
| Insertions/sec
---
(millions)
| X-small
---
queries/sec
| Small
---
queries/sec
| Medium
---
queries/sec
| Large
---
queries/sec
| Mem.
---
(MB)
| Peak mem.
---
(MB)
Brute force | 20 | 4.48 |  8.08 |  8.08 | 8.04 | 7.71 | 240 | 240
ALGLIB | 8 | 1.03 | 29300 | 1270 | 170 | 21.7 | 1180 |  3170
Boost | 200 | 2.44 | 19400 | 939 | 121 | 15 | 568 |  1500
CGAL k-d tree | 200 | 29.8 | 1660 | 1140 | 143 | 17.6 | 537 | 886
KDTREE2 | 200 | 2.12 | 1610 | 63.1 | 12.3 | 2.65 | 374 | 374
libkdtree2 | 1 |  0.00285 | 7720 | 163 | 16.1 | 1.76 |  2340 |  2340
Octree | 200 | 3.28 | 1440 | 73.1 | 14.7 | 3.21 | 234 | 234
PicoTree | 20 | 1.53 | 16200 | 815 | 132 | 23.5 | 852 | 886
R-tree | 20 |  0.166 | 27000 | 1290 | 175 | 22.5 | 1020 | 1020

####  4.2.3. Large Scale Mesh Points

Each query of extra-small, small, medium, and large size retrieves an average
of 0.00117%, 0.0996%, 1.05% and 9.84% of the process’s assigned mesh points.
Results are shown in Table 5. Although at small scale all of these libraries
have good performance for all query sizes, here we see the relative
performance decline for the KDTREE2, libkdtree2 and Octree libraries, and for
the CGAL k-d tree with the extra-small queries. The rest still perform well
across the board.

Insight: at large scale the libraries outperform the brute force solution at
all query scales. These results can help estimate when the brute force
solution offers the best performance.

Table 6. Large scale results for hexahedral mesh elements (16 processes, $\approx$11 million elements/process) for libraries that perform best at small scale. CGAL offers the best write performance, R-tree generally offers the best read performance, and Boost does well across the board. Green is good, yellow is moderate, red is poor. Configuration |  | Leaf
---
size
| Insertions/sec
---
(millions)
| X-small
---
queries/sec
| Small
---
queries/sec
| Medium
---
queries/sec
| Large
---
queries/sec
| Mem.
---
(MB)
| Peak mem.
---
(MB)
Brute Force | - | 2.43 |  4.34 |  4.34 | 4.3 | 3.97 | 473 | 473
Boost | 20 | 1.78 | 21200 | 773 | 80.7 | 8.95 | 697 | 1950
CGAL R-tree | 1 | 15.3 | 424 | 530 | 70.5 | 8.87 | 1490 | 1490
R-tree | 200 | 0.163 | 17400 | 968 | 142 | 19.2 | 1060 | 1060

####  4.2.4. Large Scale Mesh Elements

The queries of extra-small, small, medium, and large size retrieve an average
of 0.00222%, 0.120%, 1.17% and 10.5% of the mesh elements assigned to each
process. Results are shown in Table 6. Once again, we see that CGAL (this time
with an R-tree) struggles with the extra-small queries but has exceptionally
fast build throughput, and Boost performs well across the board. Here Boost is
the only library that achieves good results in terms of memory requirements
given that the raw data takes 472 MB. Once again the R-tree library has
significantly slower insertion throughput than the alternative libraries, but
offers around a 2$\times$ speedup for the medium and large queries.

###  4.3. Scalability Evaluations

In this section, we present an evaluation of the strong and weak scaling for
the libraries.

####  4.3.1. Strong Scaling

We evaluate the strong scaling by comparing the performance as we use
10$\times$ fewer process (16 instead of 160) for the same mesh. As a result,
the data per process increases from small to large (by a factor of
8.28$\times$ for the point tests and 10.0$\times$ for the element tests). The
results are shown in Table 6 and Table 8. The values shown are the factor
difference (the large result divided by the small result). For the performance
results, we classify the scaling as good (green with bolded text) if it is
better than logarithmic $(\frac{1}{\log_{2}8.28}=0.328$ for the points and
$\frac{1}{\log_{2}10.0}=0.301$ for the elements). We classify the scaling as
moderate (yellow with plain text) if it is worse than logarithmic but better
than half of logarithmic. All other results are classified as poor (red with
underlined text). For the memory results, we classify the result as good
(green) if there is a sub-linear increase in memory usage, moderate (yellow)
if the increase is greater than linear but less than a 1.5$\times$ increase,
and poor (red) if there is a greater than 1.5$\times$ increase.

Table 7. Strong scaling results for mesh points. Few libraries have good query scalability. Green is good, yellow is moderate, red is poor66footnotemark: 6. Configuration |  | Leaf
---
size
Insertions/sec |  | X-small
---
queries/sec
| Small
---
queries/sec
| Medium
---
queries/sec
| Large
---
queries/sec
| Mem.
---
(MB)
| Peak mem.
---
(MB)
Brute force | - | 0.81 | 0.28 | 0.28 | 0.28 | 0.28 | 8.22 | 8.22
ALGLIB | 8 | 0.93 | 0.81 | 0.54 | 0.41 | 0.31 | 7.76 | 8.26
Boost | 200 | 0.85 | 0.41 | 0.26 | 0.22 | 0.25 | 7.18 | 8.11
CGAL k-d tree | 200 | 0.90 | 0.17 | 0.28 | 0.23 | 0.28 | 6.96 | 7.98
KDTREE2 | 200 | 0.75 |  0.15 |  0.11 |  0.13 |  0.14 | 8.16 | 8.16
libkdtree2 | 1 | 0.24 | 0.32 |  0.16 |  0.15 |  0.14 | 8.21 | 8.21
Octree | 200 | 0.67 | 0.18 |  0.14 |  0.15 |  0.15 | 7.95 | 7.95
PicoTree | 20 | 0.72 | 0.69 | 0.50 | 0.42 | 0.41 | 8.11 | 7.98
R-tree | 20 | 0.91 | 0.64 | 0.47 | 0.39 | 0.33 | 8.29 | 8.29

#####  4.3.1.1 Mesh Points

As shown in Table 6, apart from libkdtree2, the libraries scale very well in
terms of insertion throughput and they almost uniformly have a sub-linear
increase in memory usage. However, KDTREE2, libkdtree2 and Octree experience
significantly worse than logarithmic scaling for most or all query sizes. Only
PicoTree and R-tree experience good scaling across the board, with ALGLIB
experiencing good scaling at most query sizes. CGAL’s k-d tree experiences
moderate scaling for most query sizes.

Table 8. Strong scaling results for hexahedral mesh elements. The libraries generally have good scalability apart from large queries and memory overheads. Green is good, yellow is moderate, red is poor88footnotemark: 8. Configuration |  | Leaf
---
size
Insertions/sec |  | X-small
---
queries/sec
| Small
---
queries/sec
| Medium
---
queries/sec
| Large
---
queries/sec
| Mem.
---
(MB)
| Peak mem.
---
(MB)
Brute force | - | 1.09 | 0.27 | 0.27 | 0.27 | 0.27 | 10.0 | 10.0
Boost | 20 | 0.84 | 0.59 | 0.34 | 0.30 | 0.27 | 11.3 | 12.2
CGAL R-tree | 20 | 1.00 |  0.10 | 0.37 | 0.31 | 0.25 |  19.8 |  15.5
R-tree | 200 | 0.91 | 0.52 | 0.39 | 0.34 | 0.28 | 10.4 | 10.4

#####  4.3.1.2 Mesh Elements

From Table 8 we can see that all libraries achieve good scaling for insertion
throughput and just miss the cutoff for good scaling for the large query size.
For the remaining query sizes, all achieve good scaling apart from CGAL’s
R-tree with the extra-small query size. Interestingly, they all experience a
greater than linear increase in memory usage (in contrast to their sub-linear
scaling for point storage), and the CGAL R-tree requires almost 20$\times$ the
memory despite storing only 10$\times$ as many elements.

####  4.3.2. Weak Scaling

To evaluate the weak scaling, we use the same setup used for the performance
evaluations with one modification: we only perform 1000 extra-small, small,
and medium queries and 100 large queries (1/10th as many queries of each
category). This number of queries is sufficient to allow us to compare
differences between the libraries. Where more than the baseline number of
processes (160 and 16) are used, several processes write and read the same
mesh data. This allows us to ensure that the performance differences are the
result of the weak scaling changes rather than changes in the mesh data. For
space reasons, we only present the results for a representative subset of the
libraries. For each of the four experiment setups, we present scaling results
for one library that performed (more or less) well across the board and one
that performed moderately well. This demonstrates how a library should be
expected to scale depending on what performance category it falls into. We
consider two different weak scaling scenarios: when the number of processes
per node is held constant and when the number is varied. The results are shown
in Tables LABEL:tab:weak-scaling-fixed-procs-per-node and 10. For both
scenarios, we compare the results to the performance for the baseline setup
used in the performance evaluations (a single node, with 160 processes for a
small amount of data per process or 16 processes for a large amount of data).
We classify the scaling as good (indicated in green with bold text) if the
performance is at most 10% worse than the baseline, moderate (indicated in
yellow with plain text) if it is $(10\%-25\%]$ worse than the baseline, and
poor (indicated in red with underlined text) if the performance is over 25%
worse. In each of the tables, the results indicate the percent change from the
baseline with a negative number indicating the operation takes less time than
the baseline and a positive number indicating the operation takes longer than
the baseline.

Table 9. Weak scaling results for a representative subset of the libraries using a fixed number of processes per node. As expected, the libraries all exhibit approximately linear weak scaling. Green is good, yellow is moderate555A complete explanation of the classification scheme can be found in Section 4.3.2.. Configuration |  | Leaf
---
Size
Data |  | Data per
---
proc
Procs | Nodes |  | Procs per
---
node
Insertions |  | X-small
---
queries
| Small
---
queries
| Medium
---
queries
| Large
---
queries
CGAL k-d tree | 20 | Points | Small | 4800 | 30 | 160 | 4.7% | -6.4% | -1.1% | 2.0% | 2.2%
KDTREE | 1 | Points | Small | 4800 | 30 | 160 | 0.4% | -4.6% | 0.7% | -0.6% | -2.0%
Boost | 200 | Elems | Small | 4800 | 30 | 160 | -11.0% | -2.9% | 4.8% | -2.7% | -7.9%
Spatial | 20 | Elems | Small | 4800 | 30 | 160 | -2.1% | -48.2% | -2.6% | -0.1% | -2.2%
ALGLIB | 8 | Points | Large | 480 | 30 | 16 | 0.5% | 7.9% | 3.4% | 5.5% | 0.3%
KDTREE2 | 200 | Points | Large | 480 | 30 | 16 | 0.0% | 1.4% | 2.7% | 2.7% | 2.3%
Boost | 200 | Elems | Large | 480 | 30 | 16 | 11.6% | 19.4% | 9.8% | 8.7% | 1.5%
R-tree | 200 | Elems | Large | 480 | 30 | 16 | 4.7% | 0.9% | -6.3% | -1.7% | -5.5%

#####  4.3.2.1 Fixed Number of Processes per Node

When the number of processes per node is held constant, we would expect to
find linear weak scaling (constant performance when the work per process is
held constant) since each process performs independent work and does not use
shared resources such as the network or file system. To evaluate this, we look
at how the performance changes when using 15 or 30 nodes (compared to 1 node).
The 15 and 30 node results show very close agreement (within 5%) so we omit
the 15 node results for space reasons. The results are shown in Table
LABEL:tab:weak-scaling-fixed-procs-per-node.

Insight: the libraries achieve approximately linear weak scaling in all cases.
This indicates our results should hold at any scale if the processes per node
and work per process are held constant.

Table 10. Weak scaling results for a representative subset of the libraries using a variable number of processes per node . Some libraries are more sensitive than others to the number of processes per node. Green is good, yellow is moderate, red is poor. Configuration |  | Leaf
---
Size
Data |  | Data per
---
proc
Nodes | Procs |  | Procs per
---
node
Insertions |  | X-small
---
queries
| Small
---
queries
| Medium
---
queries
| Large
---
queries
CGAL k-d tree | 20 | Points | Small | 2 | 160 | 80 | 21.7% | 1.3% | 6.6% | 3.8% | -16.4%
CGAL k-d tree | 20 | Points | Small | 8 | 160 | 20 |  38.9% | -48.5% | -19.7% | -28.5% | -79.2%
KDTREE | 1 | Points | Small | 2 | 160 | 80 | 7.7% | 0.4% | 1.6% | -0.3% | -1.9%
KDTREE | 1 | Points | Small | 8 | 160 | 20 | 4.7% | -5.9% | -4.5% | -5.7% | -6.7%
Boost | 200 | Elems | Small | 2 | 160 | 80 | 9.1% | 15.4% |  33.7% | 19.6% | 11.6%
Boost | 200 | Elems | Small | 8 | 160 | 20 | -17.5% | -55.1% | -49.4% | -67.0% | -43.5%
Spatial | 20 | Elems | Small | 2 | 160 | 80 | 12.5% | 1.0% | 0.9% | -0.7% | -2.7%
Spatial | 20 | Elems | Small | 8 | 160 | 20 | -12.2% | -56.9% | -19.6% | -18.0% | -20.5%
ALGLIB | 8 | Points | Large | 1 | 32 | 32 | -0.9% | -1.9% | 1.3% | -4.4% | -14.7%
ALGLIB | 8 | Points | Large | 1 | 80 | 80 | 17.9% |  26.3% | 22.2% | 16.0% | 2.5%
KDTREE2 | 200 | Points | Large | 1 | 32 | 32 | -2.7% | -2.4% | 0.1% | -0.2% | -6.8%
KDTREE2 | 200 | Points | Large | 1 | 80 | 80 | 4.9% | 0.8% | 6.2% | 6.2% | -2.5%
Boost | 200 | Elems | Large | 1 | 32 | 32 | 14.7% | 20.3% | 6.1% | -0.5% | -6.5%
Boost | 200 | Elems | Large | 1 | 64 | 64 |  25.8% |  40.7% | 24.0% | 14.1% | -0.9%
R-tree | 200 | Elems | Large | 1 | 32 | 32 | 7.1% | -7.7% | -18.4% | -13.4% | -19.5%
R-tree | 200 | Elems | Large | 1 | 64 | 64 | 7.2% | 10.8% | 5.8% | 12.6% | 5.4%

#####  4.3.2.2 Variable Number of Processes per Node

Next, we evaluate weak scaling when a different number of processes per node
is used. When a small amount of data per process is stored, we are already
maxing out the processes per node (160). We therefore evaluate when 20, 40,
60, and 80 processes per node are used. For space reasons we only present the
results for 20 and 80 processes per node. In this case, we use the same number
of total processes (160) but vary the number of nodes. For the case when a
large amount of data per process is stored, we only use 16 nodes per process
(far below the maximum), and therefore keep the number of nodes constant (one)
but increase the number of processes per node. For the mesh point case, we
present results for when 32 and 80 processes per node are used. For the mesh
elements case we present results for when 32 and 64 processes per node are
used because with 80 processes we exhaust the node’s memory. The results,
shown in Table 10, are much more varied. The libraries generally scale well
when the number of processes per node is only 2$\times$ different from
baseline. However, when the number of processes per node is 4, 5 or 8$\times$
different the results are sometimes greatly affected. This can be seen when
using 64 or 80 processes per node with a large amount of data per process and
with 20 processes per node with a small amount of data per process. For
example, with a small number of points per process, the CGAL k-d tree (for
points) and Boost library (for elements) both perform dramatically better when
only 20 processes per node are used. With a large amount of data per process,
we can see that ALGLIB and Boost both perform substantially worse when the
number of processes per node is increased. Thus, some of the libraries are
much more sensitive to the number of processes per node that are being
utilized.

44footnotetext: A complete explanation of the classification scheme can be
found in Section 4.3.1.

.

###  4.4. Applicability of results

Although we only use a single mesh type in this evaluation (unstructured,
hexahedral), we can use these results to estimate the performance we would get
with other mesh types. As shown in Tables 2 through 6, query performance is
greatly affected by the query size (percent of stored mesh points or elements
that match the query), exhibiting an approximately linear relationship.
Therefore, any mesh with more densely clustered points or elements near the
range being queried should result in decreased performance. Another factor
that should affect performance is the shape of the mesh elements. All of the
libraries that support mesh elements (2D or 3D shapes) use axis aligned
bounding boxes (AABBs) to determine potential intersections. Thus, the more a
mesh element deviates from an AABB, the worse the performance because of the
increase in tree branches that must be explored and potential matches
identified, and the need to use more complex intersection tests to achieve
accurate results. Beyond these factors, we would not expect much of a
performance difference. Although many of the libraries use data structures
that are not necessarily balanced (e.g., octrees or k-d trees), almost all of
the indices are balanced. This is accomplished through a combination of bulk
loading and re-balance functions. Therefore, we would not expect to find tree-
height related performance differences, or edge cases producing dramatically
unbalanced trees. We leave as future work directly evaluating how the use of
different mesh types affects these results.

Insight: our query results can be used to extrapolate expected performance for
different mesh types.

##  5\. Discussion, Insights and Future Work

In this work, we set out to answer the following questions regarding the 20
free, open-source C/C++ libraries that support range queries:

  1. (1)

Which of the implementations are viable in an HPC setting?

  2. (2)

How do these libraries compare in terms of build time, query time, and memory
usage at different scales?

  3. (3)

What are other factors in deciding which library to use?

  4. (4)

Is there a single overall best solution?

  5. (5)

When does a brute force solution offer the best performance?

The second question has been thoroughly addressed in the evaluation section
(see Section 4). In this section we summarize answers to the remaining
questions using the evaluation results and offer a number of insights that can
assist both HPC application scientists and spatial index developers. We also
present areas for future work.

###  5.1. Which of the Implementations Are Viable in an HPC Setting?

As mentioned above, to be viable in an HPC setting, a spatial indexing library
must have fast performance for building the index and performing queries, good
scalability, and moderate memory usage. In this section, we discuss which
libraries best meet these criteria.

Overall, the results show that for point storage, ALGLIB, Boost, CGAL’s k-d
tree, PicoTree and R-tree achieve the best query results. KDTREE2, libkdtree2
and Octree perform well at small scale but achieve worse performance at large
scale. Of these five best performing libraries, R-tree experiences orders of
magnitude worse insertion throughput. Therefore, if the application will be
performing a small number of queries (thus resulting in the index construction
being a larger proportion of work), the R-tree library should not be used, and
CGAL’s k-dtree, which offers almost an order of magnitude improvement over the
next best library, should be strongly considered. In addition, both the R-tree
library and ALGLIB require around 4$\times$ more memory than is required for
the raw data points. Therefore, these are not viable solutions for
particularly memory constrained environments. However, for point storage the
R-tree library does offer some of the fastest query throughputs. In fact, if
many medium and large queries will be performed at large scale the R-tree
library offers the fastest solution providing nearly double the query
throughput for medium and large queries.

For storing mesh elements, Boost, CGAL’s R-tree and the R-tree library offer
the best performance. CGAL’s R-tree offers close to a 10$\times$ improvement
in insertion throughput over Boost and a 100$\times$ throughput improvement
over R-tree both at large and small scales, but generally offers the worst
performance of the three for all query sizes (and particularly poor
performance for extra-small queries). It is therefore a viable solution if a
relatively small number of queries will be performed. Boost has by far the
smallest memory usage (although it has the highest peak memory usage) and is
thus the best solution for moderately memory constrained environments that can
tolerate temporary memory pressure. Finally, the R-tree library almost always
offers the best query performance for any query size, performing around
2$\times$ better than both Boost and CGAL for the medium and large queries at
both small scale and large scale. It is therefore a viable solution when
slightly higher memory usage can be tolerated, and a large number of queries
will be performed.

It is worth emphasizing that all of the best performing libraries use
significantly more memory than the raw data size and therefore are not
suitable for severely memory constrained environments. If memory is severely
constrained, we would recommend nanoflann with a data adapter (which is very
simple to write), which results in sub-linear memory usage, or the Octree
library for point storage and the Boost library for element storage. For
severely space constrained environments, there is another simple option: store
floats rather than doubles. This will cut memory requirements approximately in
half. We expect that four bytes worth of precision will be sufficient for
range queries in most cases given that the spatial index is not taking the
place of the mesh coordinate storage, but rather is offering an additional
data structure on top (users therefore are not losing mesh precision, simply a
small amount of query precision). It should also be kept in mind that it is
the number of mesh points per process that determines the total structure size
rather than the dataset size. Most simulations store a large number of
variables (anywhere from 5-100) and a large number of timesteps (anywhere from
100-100,000) and the number of mesh coordinates or elements will therefore be
a relatively small fraction of the overall dataset. It is also worth
emphasizing that, given how quickly these structures can be created, they can
easily be generated on the fly and therefore will not require any long-term
storage.

###  5.2. What Are Other Considerations in Deciding Which Library to Use?

There are a large number of different factors to consider when deciding which
library to use. Due to space constraints, we will only discuss a few of the
most important additional factors here, but a more complete overview can be
found in the Git repository for this project (Lawson, 2021).

####  5.2.1. Ease of Use

All of the libraries are easy to use when compiling a project. CGAL, PicoTree
and R-tree are header only libraries. For ALGLIB you simply link your
executable to the necessary sources, and for Boost you link to the library,
which is widely available and easy to compile. However, the libraries differ
greatly in terms of their flexibility for data storage. ALGLIB can only store
numeric data (coordinates, integer tags, and N additional real values per
coordinate). Boost and CGAL offer significant flexibility in what data can be
stored, but require the data to be based on classes defined in the library.
For Boost this requires containing a boost::geometry::model::point/box as part
of what is inserted (for point and box storage). For CGAL, this is more
complex. Here is an overview of the definitions we used to create a k-d tree
that contains 3D points and query it using a bounding box.

⬇

typedef CGAL::Simple_cartesian<double> cgal_kernel;

typedef cgal_kernel::Point_3 cgal_pt;

typedef boost::tuple<cgal_pt,uint32_t> cgal_pt_w_index;

typedef CGAL::Search_traits_3<cgal_kernel> cgal_traits;

typedef CGAL::Search_traits_adapter<cgal_pt_w_index,

CGAL::Nth_of_tuple_property_map<0, cgal_p_w_index>,

cgal_traits> Traits;

typedef CGAL::Sliding_midpoint<Traits> cgal_split_rule;

typedef CGAL::Kd_tree<Traits> cgal_kd_tree;

typedef CGAL::Fuzzy_iso_box<Traits> cgal_box;

CGAL is nearly infinitely flexible, but the extensive documentation only
demonstrates a small fraction of how its capabilities work and can be
combined, resulting in a substantial challenge even for experienced coders.
PicoTree allows you to very easily create a custom class and is incredibly
simple to use overall. The R-tree library is the simplest to adapt to the
evaluated use-case since it allows you to indicate an element data type, and
thus does not require any modifications (or custom classes) for use in this
work.

####  5.2.2. Support and Documentation

Of the five best performing libraries, only ALGLIB offers commercial support.
ALGLIB’s commercial version also provides multithreading support and SIMD
usage. ALGLIB, Boost and CGAL all offer substantial documentation. PicoTree
and R-tree offer little to no documentation but this is somewhat offset by the
fact that they have much smaller and simpler code bases, making it much easier
to find relevant functions.

####  5.2.3. Static vs. Dynamic Trees

ALGLIB offers only static trees, Boost offers dynamic trees (with support for
bulk loading), CGAL’s k-d tree supports deletion but is not self balancing,
CGAL’s R-tree is static, PicoTree offers static trees only, and R-tree offers
dynamic trees. Therefore, only Boost, CGAL’s k-d tree or R-tree would be
suitable for dynamic environments.

####  5.2.4. Theoretical Performance Bounds

All of the five best performing libraries for storing points use either a k-d
tree (ALGLIB, CGAL’s k-d tree, PicoTree) or an R-tree (Boost, the R-tree
library). For element storage, all of the best performing libraries use an
R-tree (Boost, CGAL’s R-tree and the R-tree library). Both k-d trees and
R-trees can answer a range query in $O(n^{1-1/d}+k)$ time (Agarwal, 1996)
whereas octrees can answer range queries in $O(\log n)$ time (Jackins and
Tanimoto, 1980). Therefore, it is perhaps a bit surprising that none of the
best performing implementations are octrees. Given the wide range of
performance we see across the libraries for each tree type, we should expect
the library’s implementation to be a much greater factor in determining
performance than the theoretical bounds of the data structure.

###  5.3. Is There a Single Overall Best Solution?

As the evaluation section makes clear, there is no universally good library
that achieves fast build throughput, fast query throughput for all query sizes
and uses close to linear memory. However, Boost, CGAL and R-tree offer some of
the best performance at large and small scales and have the advantage of being
able to support both points and boxes. Of the three libraries, Boost offers
the best overall performance in terms of good build and query times (and has
the lowest memory usage), while CGAL offers the best build times and R-tree
offers the best query times. There is therefore no single best solution, but
rather the answer will depend on the problem scale, memory availability, the
number of queries to be performed and the fraction of stored data that is
expected to match the queries. In addition, as discussed in the next section,
there are circumstances under which the best library to use is no library.

###  5.4. Under What Circumstances Does a Brute Force Solution Offer the Best
Performance?

The results demonstrate that if many of the queries will retrieve 10% or more
of the stored mesh data, and if there is a small amount of data per process
(e.g., less than 10 million points or elements), then one should strongly
consider using a brute force solution. The libraries evaluated in this paper
offer the best performance when the search space can be substantially reduced
(pruning branches) and pay a large performance penalty for the $\log n$ tree
traversals when a large portion of the tree must be searched. In addition, as
discussed previously, all of the libraries require significantly more memory
than the raw data requires, and therefore the brute force solution should be
used if users are severely memory constrained.

###  5.5. Future Work

We leave as future work evaluating these libraries for dynamic use (e.g., if
adaptive mesh refinement is used) and for mesh elements that are more complex
shapes.  We also plan to assess the portability and software design choices of
the libraries, and to evaluate how the use of different architectures affects
the evaluation results. Finally, we plan to evaluate libraries for languages
other than C and C++.

## Acknowledgements

This paper describes objective technical results and analysis. Any subjective
views or opinions that might be expressed in the paper do not necessarily
represent the views of the U.S. Department of Energy or the United States
Government.

Sandia National Laboratories is a multimission laboratory managed and operated
by National Technology and Engineering Solutions of Sandia, LLC, a wholly
owned subsidiary of Honeywell International, Inc., for the U.S. Department of
Energy’s National Nuclear Security Administration under contract DE-NA0003525.

Margaret Lawson acknowledges support from the United States Department of
Energy through the Computational Sciences Graduate Fellowship (DOE CSGF) under
grant number: DE- SC0020347.

This work was supported in part by the State of Illinois.

## References

  * (1)
  * val (2020) 2020\.  _Massif: a heap profiler_.  Retrieved Jan 20, 2021 from https://www.valgrind.org/docs/manual/ms-manual.html
  * Agarwal (1996) Pankaj K Agarwal. 1996\.  _Range Searching_.  Technical Report. DUKE UNIV DURHAM NC DEPT OF COMPUTER SCIENCE.
  * Aji et al. (2012) Ablimit Aji, Fusheng Wang, and Joel H Saltz. 2012\.  Towards building a high performance spatial query system for large scale medical imaging data. In _Proceedings of the 20th international conference on advances in geographic information systems_. 309–318.
  * Amaral et al. (2020) Vasco Amaral, Beatriz Norberto, Miguel Goulão, Marco Aldinucci, Siegfried Benkner, Andrea Bracciali, Paulo Carreira, Edgars Celms, Luís Correia, Clemens Grelck, et al. 2020\.  Programming languages for data-Intensive HPC applications: A systematic mapping study.  _Parallel Comput._ 91 (2020), 102584.
  * Anon (2020) Anon. 2020.  _KdTree_.  Retrieved January 20, 2021 from https://github.com/G3tupup/KdTree
  * Arya and Mount (1998) Sunil Arya and DM Mount. 1998.  ANN: library for approximate nearest neighbor searching. In _Proceedings of IEEE CGC Workshop on Computational Geometry, Providence, RI_.
  * Barone et al. (2021) Matthew F Barone, Jaideep Ray, and Stefan Domino. 2021\.  Feature Selection, Clustering, and Prototype Placement for Turbulence Data Sets. In _AIAA Scitech 2021 Forum_. 1750.
  * Behley et al. (2015) Jens Behley, Volker Steinhage, and Armin B. Cremers. 2015\.  Efficient Radius Neighbor Seach in Three-dimensional Point Clouds. In _Proc. of the IEEE International Conference on Robotics and Automation (ICRA)_.
  * Bentley (1975) Jon Louis Bentley. 1975\.  Multidimensional binary search trees used for associative searching.  _Commun. ACM_ 18, 9 (1975), 509–517.
  * Bentley and Maurer (1980) Jon Louis Bentley and Hermann A. Maurer. 1980.  Efficient worst-case data structures for range searching.  _Acta Informatica_ 13, 2 (1980), 155–168.
  * Betoule et al. (2013) M Betoule, J Marriner, N Regnault, J-C Cuillandre, P Astier, J Guy, C Balland, P El Hage, D Hardin, R Kessler, et al. 2013\.  Improved photometric calibration of the SNLS and the SDSS supernova surveys.  _Astronomy & Astrophysics_ 552 (2013), A124.
  * Beyer et al. (2013) Johanna Beyer, Ali Al-Awami, Narayanan Kasthuri, Jeff W Lichtman, Hanspeter Pfister, and Markus Hadwiger. 2013.  ConnectomeExplorer: Query-guided visual analysis of large volumetric neuroscience data.  _IEEE transactions on visualization and computer graphics_ 19, 12 (2013), 2868–2877.
  * Beynon et al. (2000) Michael Beynon, Renato Ferreira, Tahsin Kurc, Alan Sussman, and Joel Saltz. 2000.  DataCutter: Middleware for filtering very large scientific datasets on archival storage systems. In _IEEE Symposium on Mass Storage Systems_. 119–134.
  * Blanco and Rai (2014) Jose Luis Blanco and Pranjal Kumar Rai. 2014.  _nanoflann: a C++ header-only fork of FLANN, a library for Nearest Neighbor (NN) wih KD-trees_.  Retrieved Jan 20, 2021 from https://github.com/jlblancoc/nanoflann
  * Bochkanov and Bystritsky (2011) Sergey Bochkanov and Vladimir Bystritsky. 2011.  ALGLIB-a cross-platform numerical analysis and data processing library.  _ALGLIB Project. Novgorod, Russia_ (2011).
  * Bougerel (2017) Sylvain Bougerel. 2017\.  _Spatial C++ Library_.  Retrieved Jan 20, 2021 from http://spatial.sourceforge.net/index.html
  * Bremer et al. (2010) Peer-Timo Bremer, Gunther Weber, Julien Tierny, Valerio Pascucci, Marc Day, and John Bell. 2010\.  Interactive exploration and analysis of large-scale simulations using topology-based data segmentation.  _IEEE Transactions on Visualization and Computer Graphics_ 17, 9 (2010), 1307–1324.
  * Broere (2021) Jonathan Broere. 2021\.  _PicoTree_.  Retrieved Jan 20, 2021 from https://github.com/Jaybro/pico_tree
  * Carvalho (2020) J. Frederico Carvalho. 2020\.  _KDTree_.  Retrieved Jan 20, 2021 from https://github.com/crvs/KDTree
  * Chazelle (1990) Bernard Chazelle. 1990\.  Lower bounds for orthogonal range searching: I. the reporting case.  _Journal of the ACM (JACM)_ 37, 2 (1990), 200–212.
  * Chiu et al. (2015) Hsuan-Te Chiu, Jerry Chou, Venkat Vishwanath, and Kesheng Wu. 2015.  In-memory query system for scientific datasets. In _2015 IEEE 21st International Conference on Parallel and Distributed Systems (ICPADS)_. IEEE, 362–371.
  * Colli Franzone et al. (2018) Piero Colli Franzone, Luca F Pavarino, and Simone Scacchi. 2018\.  A numerical study of scalable cardiac electro-mechanical solvers on HPC architectures.  _Frontiers in physiology_ 9 (2018), 268.
  * Cong and Jensen (2016) Gao Cong and Christian S Jensen. 2016.  Querying geo-textual data: Spatial keyword queries and beyond. In _Proceedings of the 2016 International Conference on Management of Data_. 2207–2212.
  * Cubuk et al. (2015) Ekin D Cubuk, Samuel Stern Schoenholz, Jennifer M Rieser, Brad Dean Malone, Joerg Rottler, Douglas J Durian, Efthimios Kaxiras, and Andrea J Liu. 2015\.  Identifying structural flow defects in disordered solids using machine-learning methods.  _Physical review letters_ 114, 10 (2015), 108001\.
  * Cummings et al. (2008) J Cummings, A Pankin, N Podhosrzki, G Park, S Ku, R Barreto, S Klasky, CS Chang, H Strauss, L Sugiyama, et al. 2008\.  Plasma edge kinetic-MHD modeling in tokamaks using Kepler workflow for code coupling, data management and visualization.  _Communications in Computational Physics_ 4, 3 (2008), 675–702.
  * de Oliveira et al. (2020) Matheus A de Oliveira, Marcelo de M Menezes, Salles VG de Magalhaes, Bruno F Coelho, and Viçosa Campus da UFV. 2020\.  Evaluating the usage of exact queries on 3D spatial databases.  (2020), 162–167.
  * Demirbas and Ferhatosmanoglu (2003) Murat Demirbas and Hakan Ferhatosmanoglu. 2003.  Peer-to-peer spatial queries in sensor networks. In _Proceedings Third International Conference on Peer-to-Peer Computing (P2P2003)_. IEEE, 32–39.
  * Dietrich (2013) Jörg Dietrich. 2013\.  _libkdtree_.  Retrieved Jan 20, 2021 from https://github.com/joergdietrich/libkdtree
  * Domino et al. (2019) Stefan P. Domino, Philip Sakievich, and Matthew Barone. 2019\.  An assessment of atypical mesh topologies for low-Mach large-eddy simulation.  _Computers & Fluids_ 179 (2019), 655–669.  https://doi.org/10.1016/j.compfluid.2018.12.002
  * Egenhofer (1994) Max J. Egenhofer. 1994\.  Spatial SQL: A query and presentation language.  _IEEE Transactions on knowledge and data engineering_ 6, 1 (1994), 86–95.
  * Eldawy and Mokbel (2015) Ahmed Eldawy and Mohamed F Mokbel. 2015.  Spatialhadoop: A mapreduce framework for spatial data. In _2015 IEEE 31st international conference on Data Engineering_. IEEE, 1352–1363.
  * Elseberg et al. (2012) Jan Elseberg, Stéphane Magnenat, Roland Siegwart, and Andreas Nüchter. 2012.  Comparison of nearest-neighbor-search strategies and implementations for efficient shape registration.  _Journal of Software Engineering for Robotics_ 3, 1 (2012), 2–12.
  * Fabri and Pion (2009) Andreas Fabri and Sylvain Pion. 2009.  CGAL: The Computational Geometry Algorithms Library. In _Proceedings of the 17th ACM SIGSPATIAL International Conference on Advances in Geographic Information Systems_ (Seattle, Washington) _(GIS ’09)_. Association for Computing Machinery, New York, NY, USA, 538–539.  https://doi.org/10.1145/1653771.1653865
  * Feki and Gabriel (2008) Saber Feki and Edgar Gabriel. 2008.  Incorporating historic knowledge into a communication library for self-optimizing high performance computing applications. In _2008 Second IEEE International Conference on Self-Adaptive and Self-Organizing Systems_. IEEE, 265–274.
  * Formaggia (1999) Luca Formaggia. 1999\.  Data structures for unstructured mesh generation.  (1999).
  * Gabriel et al. (2010) Edgar Gabriel, Saber Feki, Katharina Benkert, and Michael M Resch. 2010.  Towards performance portability through runtime adaptation for high-performance computing applications.  _Concurrency and Computation: Practice and Experience_ 22, 16 (2010), 2230–2246.
  * Gahvari et al. (2011) Hormozd Gahvari, Allison H Baker, Martin Schulz, Ulrike Meier Yang, Kirk E Jordan, and William Gropp. 2011\.  Modeling the performance of an algebraic multigrid cycle on HPC platforms. In _Proceedings of the international conference on Supercomputing_. 172–181.
  * Garcia-Gasulla et al. (2020) Marta Garcia-Gasulla, Filippo Mantovani, Marc Josep-Fabrego, Beatriz Eguzkitza, and Guillaume Houzeaux. 2020.  Runtime mechanisms to survive new HPC architectures: a use case in human respiratory simulations.  _The International Journal of High Performance Computing Applications_ 34, 1 (2020), 42–56.
  * García R et al. (1998) Yván J García R, Mario A López, and Scott T Leutenegger. 1998.  A greedy algorithm for bulk loading R-trees. In _Proceedings of the 6th ACM international symposium on Advances in geographic information systems_. 163–164.
  * Gehrels et al. (2020) Barend Gehrels, Bruno Lalande, Mateusz Loskot, Adam Wulkiewicz, Menelaos Karavelas, and Vissarion Fisikopoulos. 2020.  _Geometry_.  Retrieved Jan 20, 2021 from https://www.boost.org/doc/libs/1_75_0/libs/geometry/doc/html/index.html
  * Godenschwager et al. (2013) Christian Godenschwager, Florian Schornbaum, Martin Bauer, Harald Köstler, and Ulrich Rüde. 2013.  A framework for hybrid parallel flow simulations with a trillion cells in complex geometries. In _Proceedings of the International Conference on High Performance Computing, Networking, Storage and Analysis_. 1–12.
  * Gosink et al. (2007) Luke Gosink, John Anderson, Wes Bethel, and Kenneth Joy. 2007.  Variable interactions in query-driven visualization.  _IEEE Transactions on Visualization and Computer Graphics_ 13, 6 (2007), 1400–1407.
  * Gosink et al. (2010) Luke J Gosink, Christoph Garth, John C Anderson, E Wes Bethel, and Kenneth I Joy. 2010.  An application of multivariate statistical analysis for query-driven visualization.  _IEEE Transactions on Visualization and Computer Graphics_ 17, 3 (2010), 264–275.
  * Graff et al. (2014) Philip Graff, Farhan Feroz, Michael P Hobson, and Anthony Lasenby. 2014.  SKYNET: an efficient and robust neural network training tool for machine learning in astronomy.  _Monthly Notices of the Royal Astronomical Society_ 441, 2 (2014), 1741–1759.
  * Green and Douglas (2004) Melinda Green and Greg Douglas. 2004.  _R-Trees_.
  * Güting (1994) Ralf Hartmut Güting. 1994\.  An introduction to spatial database systems.  _the VLDB Journal_ 3, 4 (1994), 357–399.
  * Guttman (1984) Antonin Guttman. 1984\.  R-trees: A dynamic index structure for spatial searching. In _Proceedings of the 1984 ACM SIGMOD international conference on Management of data_. 47–57.
  * Hadjieleftheriou (2015) Marios Hadjieleftheriou. 2015\.  LibSpatialIndex.
  * Hadjieleftheriou et al. (2005) Marios Hadjieleftheriou, George Kollios, Petko Bakalov, and Vassilis J Tsotras. 2005.  Complex spatio-temporal pattern queries. In _VLDB_ , Vol. 5. 877–888.
  * Hinterberger et al. (1994) Hans Hinterberger, Kathrin Anne Meier, and Hans Gilgen. 1994\.  Spatial data reallocation based on multidimensional range queries. a contribution to data management for the earth sciences. In _Seventh International Working Conference on Scientific and Statistical Database Management_. IEEE, 228–239.
  * Houzeaux et al. (2011) G Houzeaux, R de la Cruz, and M Vázquez. 2011.  Parallel uniform mesh subdivision in alya.  _Partnership for Advanced Computing in Europe, Available on-line at www. prace-ri. eu_ (2011).
  * Hsu et al. (2011) Ming-Chen Hsu, Ido Akkerman, and Yuri Bazilevs. 2011\.  High-performance computing of wind turbine aerodynamics using isogeometric analysis.  _Computers & Fluids_ 49, 1 (2011), 93–100.
  * Jackins and Tanimoto (1980) Chris L Jackins and Steven L Tanimoto. 1980.  Oct-trees and their use in representing three-dimensional objects.  _Computer Graphics and Image Processing_ 14, 3 (1980), 249–270.
  * Jagadish (1990) H Vi Jagadish. 1990\.  Spatial search with polyhedra. In _Proceedings. Sixth International Conference on Data Engineering_. IEEE Computer Society, 311–312.
  * Kennel (2004) Matthew B Kennel. 2004\.  KDTREE 2: Fortran 95 and C++ software to efficiently search for near neighbors in a multi-dimensional Euclidean space.  _arXiv:physics/0408067_ (2004).
  * Kim et al. (2013) Jinwoong Kim, Sul-Gi Kim, and Beomseok Nam. 2013.  Parallel multi-dimensional range query processing with R-trees on GPU.  _J. Parallel and Distrib. Comput._ 73, 8 (2013), 1195–1207.
  * Krafft (2020) Martin F. Krafft. 2020\.  _libkdtree_.  Retrieved Jan 20, 2021 from https://github.com/nvmd/libkdtree
  * Laguna et al. (2019) Ignacio Laguna, Ryan Marshall, Kathryn Mohror, Martin Ruefenacht, Anthony Skjellum, and Nawrin Sultana. 2019.  A large-scale study of MPI usage in open-source HPC applications. In _Proceedings of the International Conference for High Performance Computing, Networking, Storage and Analysis_. 1–14.
  * Landge et al. (2014) Aaditya G Landge, Valerio Pascucci, Attila Gyulassy, Janine C Bennett, Hemanth Kolla, Jacqueline Chen, and Peer-Timo Bremer. 2014.  In-situ feature extraction of large scale combustion simulations using segmented merge trees. In _SC’14: Proceedings of the International Conference for High Performance Computing, Networking, Storage and Analysis_. IEEE, 1020–1031.
  * Lawson (2021) Margaret Lawson. 2021\.  _Benchmarking Suite for Range Tree Libraries_.  https://github.com/mlawsonca/benchmarking_suite_range_searching_libraries
  * Leutenegger et al. (1997) Scott T Leutenegger, Mario A Lopez, and Jeffrey Edgington. 1997\.  STR: A simple and efficient algorithm for R-tree packing. In _Proceedings 13th International Conference on Data Engineering_. IEEE, 497–506.
  * Liang et al. (2016) Yanhui Liang, Hoang Vo, Ablimit Aji, Jun Kong, and Fusheng Wang. 2016.  Scalable 3d spatial queries for analytical pathology imaging with mapreduce. In _Proceedings of the 24th ACM SIGSPATIAL International Conference on Advances in Geographic Information Systems_. 1–4.
  * Liang et al. (2017) Yanhui Liang, Hoang Vo, Jun Kong, and Fusheng Wang. 2017\.  ispeed: an efficient in-memory based spatial query system for large-scale 3d data with complex structures. In _Proceedings of the 25th ACM SIGSPATIAL International Conference on Advances in Geographic Information Systems_. 1–10.
  * Liu et al. (2016) Yunjie Liu, Evan Racah, Joaquin Correa, Amir Khosrowshahi, David Lavers, Kenneth Kunkel, Michael Wehner, William Collins, et al. 2016\.  Application of deep convolutional neural networks for detecting extreme weather in climate datasets.  _arXiv preprint arXiv:1605.01156_ (2016).
  * Livnat et al. (1996) Yarden Livnat, Han-Wei Shen, and Christopher R Johnson. 1996\.  A near optimal isosurface extraction algorithm using the span space.  _IEEE transactions on Visualization and Computer Graphics_ 2, 1 (1996), 73–84.
  * Lofstead et al. (2011) Jay Lofstead, Milo Polte, Garth Gibson, Scott Klasky, Karsten Schwan, Ron Oldfield, Matthew Wolf, and Qing Liu. 2011\.  Six degrees of scientific data: reading patterns for extreme scale science IO. In _Proceedings of the 20th international symposium on High performance distributed computing_ _(HPDC ’11)_. ACM, 49–60.  http://doi.acm.org/10.1145/1996130.1996139
  * Lu et al. (2017) Hao Lu, Sudip K Seal, Wei Guo, and John Poplawsky. 2017\.  Spherical Region Queries on Multicore Architectures. In _Proceedings of the Seventh Workshop on Irregular Applications: Architectures and Algorithms_. 1–4.
  * Luo et al. (2012) Lijuan Luo, Martin DF Wong, and Lance Leong. 2012.  Parallel implementation of R-trees on the GPU. In _17th Asia and South Pacific Design Automation Conference_. IEEE, 353–358.
  * Marrinan et al. (2019) Thomas Marrinan, Greg Eisenhauer, Matthew Wolf, Joseph A Insley, Silvio Rizzi, and Michael E Papka. 2019.  Parallel streaming between heterogeneous HPC resources for real-time analysis.  _Journal of Computational Science_ 31 (2019), 163–171.
  * Matoušek (1994) Jiří Matoušek. 1994\.  Geometric range searching.  _ACM Computing Surveys (CSUR)_ 26, 4 (1994), 422–461.
  * Mehlhorn and Näher (1999) Kurt Mehlhorn and Stefan Näher. 1999.  _LEDA: A platform for combinatorial and geometric computing_.  Cambridge university press.
  * Muja and Lowe (2009) Marius Muja and David G Lowe. 2009.  Fast approximate nearest neighbors with automatic algorithm configuration.  _VISAPP (1)_ 2, 331-340 (2009), 2.
  * Nethercote and Seward (2007) Nicholas Nethercote and Julian Seward. 2007.  Valgrind: a framework for heavyweight dynamic binary instrumentation.  _ACM Sigplan notices_ 42, 6 (2007), 89–100.
  * Nguyen and Rhodes (2016) Cuong Nguyen and Philip J Rhodes. 2016.  Accelerating range queries for large-scale unstructured meshes. In _2016 IEEE International Conference on Big Data (Big Data)_. IEEE, 502–511.
  * Nouri and Tu (2018) Zhila Nouri and Yi-Cheng Tu. 2018.  GPU-based parallel indexing for concurrent spatial query processing. In _Proceedings of the 30th International Conference on Scientific and Statistical Database Management_. 1–12.
  * Nüchter and Lingemann (2011) Andreas Nüchter and K Lingemann. 2011.  _3DTK–The 3D Toolkit_.  Retrieved Jan 20, 2021 from http://slam6d.sourceforge.net
  * of Energy (2019) U.S. Department of Energy. 2019\.  _U.S. Department of Energy and Intel to Build First Exascale Supercomputer_.  U.S. Department of Energy.  https://www.energy.gov/articles/us-department-energy-and-intel-build-first-exascale-supercomputer
  * Orenstein (1986) Jack A Orenstein. 1986\.  Spatial query processing in an object-oriented database system. In _Proceedings of the 1986 ACM SIGMOD international conference on Management of data_. 326–336.
  * Overmars (1988) Mark H Overmars. 1988\.  Efficient data structures for range searching on a grid.  _Journal of Algorithms_ 9, 2 (1988), 254–275.
  * Pagel et al. (1993) Bernd-Uwe Pagel, Hans-Werner Six, Heinrich Toben, and Peter Widmayer. 1993.  Towards an analysis of range query performance in spatial data structures. In _Proceedings of the twelfth ACM SIGACT-SIGMOD-SIGART symposium on Principles of database systems_. 214–221.
  * Papadias et al. (2003) Dimitris Papadias, Jun Zhang, Nikos Mamoulis, and Yufei Tao. 2003\.  Query processing in spatial network databases. In _Proceedings 2003 VLDB Conference_. Elsevier, 802–813.
  * Papadomanolakis et al. (2006a) Stratos Papadomanolakis, Anastassia Ailamaki, Julio C Lopez, Tiankai Tu, David R O’Hallaron, and Gerd Heber. 2006a.  Efficient query processing on unstructured tetrahedral meshes. In _Proceedings of the 2006 ACM SIGMOD international conference on Management of data_. 551–562.
  * Papadomanolakis et al. (2006b) Stratos Papadomanolakis, Anastassia Ailamaki, Julio C Lopez, Tiankai Tu, David R O’Hallaron, and Gerd Heber. 2006b.  Efficient query processing on unstructured tetrahedral meshes. In _Proceedings of the 2006 ACM SIGMOD international conference on Management of data_. ACM, 551–562.
  * Pecha (2016) Marek Pecha. 2016\.  Image segmentation techniques in the HPC environment and their applications.  (2016).
  * Peña-Monferrer et al. (2021) Carlos Peña-Monferrer, Robert Manson-Sawko, and Vadim Elisseev. 2021.  HPC-cloud native framework for concurrent simulation, analysis and visualization of CFD workflows.  _Future Generation Computer Systems_ 123 (2021), 14–23.
  * Perlman et al. (2010) Eric Perlman, Randal Burns, Michael Kazhdan, Rebecca R Murphy, William P Ball, and Nina Amenta. 2010.  Organization of data in non-convex spatial domains. In _International Conference on Scientific and Statistical Database Management_. Springer, 342–359.
  * Post et al. (2003) Frits H Post, Benjamin Vrolijk, Helwig Hauser, Robert S Laramee, and Helmut Doleisch. 2003\.  The state of the art in flow visualisation: Feature extraction and tracking. In _Computer Graphics Forum_ , Vol. 22. Wiley Online Library, 775–792.
  * Pradipkumar and Rajamenakshi (2016) Kodrani Kajal Pradipkumar and R Rajamenakshi. 2016.  Segmentation of large scale medical images using HPC: Classification of methods and challenges.  _International Journal of Advance Engineering and Research (IJAERS)_ 3, 1 (2016), 56–61.
  * Prasad et al. (2015) Sushil K Prasad, Michael McDermott, Xi He, and Satish Puri. 2015\.  GPU-based Parallel R-tree Construction and Querying. In _2015 IEEE International Parallel and Distributed Processing Symposium Workshop_. IEEE, 618–627.
  * Quintanas-Corominas et al. (2018) Adrià Quintanas-Corominas, Pere Maimí, Eva Casoni, Albert Turon, Joan Andreu Mayugo, Gerard Guillamet, and Mariano Vázquez. 2018.  A 3D transversally isotropic constitutive model for advanced composites implemented in a high performance computing code.  _European Journal of Mechanics-A/Solids_ 71 (2018), 278–291.
  * Rasquin et al. (2014) Michel Rasquin, Cameron Smith, Kedar Chitale, Seegyoung Seol, Ben Matthews, J Martin, Onkar Sahni, Raymond Loy, Mark S Shephard, and Kenneth E Jansen. 2014\.  Scalable fully implicit finite element flow solver with application to high-fidelity flow control simulations on a realistic wing design.  _Computing in Science and Engineering_ 16, 6 (2014), 13–21.
  * Rusu and Cousins (2011) Radu Bogdan Rusu and Steve Cousins. 2011.  3d is here: Point cloud library (pcl). In _2011 IEEE international conference on robotics and automation_. IEEE, 1–4.
  * Salim et al. (2019) Michael A Salim, Thomas D Uram, J Taylor Childers, Prasanna Balaprakash, Venkatram Vishwanath, and Michael E Papka. 2019.  Balsam: Automated scheduling and execution of dynamic, data-intensive hpc workflows.  _arXiv preprint arXiv:1909.08704_ (2019).
  * Shekhar et al. (1999) Shashi Shekhar, Sanjay Chawla, Sivakumar Ravada, Andrew Fetterer, Xuan Liu, and Chang-tien Lu. 1999.  Spatial databases-accomplishments and research needs.  _IEEE transactions on knowledge and data engineering_ 11, 1 (1999), 45–55.
  * Simakov et al. (2018) Nikolay A Simakov, Joseph P White, Robert L DeLeon, Steven M Gallo, Matthew D Jones, Jeffrey T Palmer, Benjamin Plessinger, and Thomas R Furlani. 2018.  A Workload Analysis of NSF’s Innovative HPC Resources Using XDMoD.  _arXiv preprint arXiv:1801.04306_ (2018).
  * Soheili et al. (2005) Amir Soheili, Vana Kalogeraki, and Dimitrios Gunopulos. 2005\.  Spatial queries in sensor networks. In _Proceedings of the 13th annual ACM international workshop on Geographic information systems_. 61–70.
  * Stockinger et al. (2005) Kurt Stockinger, John Shalf, Kesheng Wu, and E Wes Bethel. 2005\.  _Query-driven visualization of large data sets_.  IEEE.
  * Tang et al. (2016) Mingjie Tang, Yongyang Yu, Qutaibah M Malluhi, Mourad Ouzzani, and Walid G Aref. 2016.  Locationspark: A distributed in-memory data management system for big spatial data.  _Proceedings of the VLDB Endowment_ 9, 13 (2016), 1565–1568.
  * Tauheed et al. (2012) Farhan Tauheed, Laurynas Biveinis, Thomas Heinis, Felix Schurmann, Henry Markram, and Anastasia Ailamaki. 2012.  Accelerating range queries for brain simulations. In _2012 IEEE 28th International Conference on Data Engineering_. IEEE, 941–952.
  * Tsiombikas (2019) John Tsiombikas. 2019\.  _kdtree_.  Retrieved January 20, 2021 from https://github.com/jtsiomb/kdtree
  * Tzeng and Ma (2005) Fan-Yin Tzeng and Kwan-Liu Ma. 2005.  Intelligent feature extraction and tracking for visualizing large-scale 4d flow simulations. In _SC’05: Proceedings of the 2005 ACM/IEEE Conference on Supercomputing_. IEEE, 6–6.
  * Wan et al. (2019) Shaohua Wan, Yu Zhao, Tian Wang, Zonghua Gu, Qammer H Abbasi, and Kim-Kwang Raymond Choo. 2019.  Multi-dimensional data indexing and range query processing via Voronoi diagram for internet of things.  _Future Generation Computer Systems_ 91 (2019), 382–391.
  * Wang et al. (2008) Chaoli Wang, Hongfeng Yu, and Kwan-Liu Ma. 2008.  Importance-driven time-varying data visualization.  _IEEE Transactions on Visualization and Computer Graphics_ 14, 6 (2008), 1547–1554.
  * Wang et al. (2017) Xiao-qing Wang, Xian-long Jin, Da-zhi Kou, and Jia-hui Chen. 2017\.  A parallel approach for the generation of unstructured meshes with billions of elements on distributed-memory supercomputers.  _International Journal of Parallel Programming_ 45, 3 (2017), 680–710.
  * Willard (1985) Dan E Willard. 1985\.  New data structures for orthogonal range queries.  _SIAM J. Comput._ 14, 1 (1985), 232–253.
  * Xu et al. (2018) Guowen Xu, Hongwei Li, Yuanshun Dai, Kan Yang, and Xiaodong Lin. 2018.  Enabling efficient and geometric range query with access control over encrypted spatial data.  _IEEE Transactions on Information Forensics and Security_ 14, 4 (2018), 870–885.
  * You et al. (2013) Simin You, Jianting Zhang, and Le Gruenwald. 2013\.  Parallel spatial query processing on gpus using r-trees. In _Proceedings of the 2Nd ACM SIGSPATIAL international workshop on analytics for big geospatial data_. 23–31.
  * Yu et al. (2010) Hongfeng Yu, Chaoli Wang, Ray W Grout, Jacqueline H Chen, and Kwan-Liu Ma. 2010.  In situ visualization for large-scale combustion simulations.  _IEEE computer graphics and applications_ 30, 3 (2010), 45–57.
  * Yu et al. (2015) Jia Yu, Jinxuan Wu, and Mohamed Sarwat. 2015.  Geospark: A cluster computing framework for processing large-scale spatial data. In _Proceedings of the 23rd SIGSPATIAL international conference on advances in geographic information systems_. 1–4.
  * Zaspel and Griebel (2011) Peter Zaspel and Michael Griebel. 2011.  Massively parallel fluid simulations on amazon’s hpc cloud. In _2011 First International Symposium on Network Cloud Computing and Applications_. IEEE, 73–78.
  * Zhao et al. (2020) Zhong Zhao, Yang Zhang, Lei He, Xinghua Chang, and Laiping Zhang. 2020.  A large-scale parallel hybrid grid generation technique for realistic complex geometry.  _International Journal for Numerical Methods in Fluids_ 92, 10 (2020), 1235–1255.
  * Zhong et al. (2012) Yunqin Zhong, Jizhong Han, Tieying Zhang, Zhenhua Li, Jinyun Fang, and Guihai Chen. 2012\.  Towards parallel spatial query processing for big spatial data. In _2012 IEEE 26th International Parallel and Distributed Processing Symposium Workshops & PhD Forum_. IEEE, 2085–2094.
  * Zhu et al. (2015) Hui Zhu, Rongxing Lu, Cheng Huang, Le Chen, and Hui Li. 2015.  An efficient privacy-preserving location-based services query scheme in outsourced cloud.  _IEEE Transactions on Vehicular Technology_ 65, 9 (2015), 7729–7739.