Modular and Efficient Approximate Nearest Neighbors Search with Vectorium and kANNolo

Modular and Efficient Approximate Nearest Neighbors Search with Vectorium and kANNolo Silvio Martinicoa,b,∗ , Franco Maria Nardinib , Cosimo Rullib,1 and Rossano Venturinia a Department of Computer Science, University of Pisa, Largo Bruno Pontecorvo 3, 56127, Pisa, Italy b Institute of Information Science and Technologies “Alessandro Faedo” (ISTI), National Research Council of Italy (CNR), Via Giuseppe Moruzzi

arXiv:2611.03927v1 [cs.IR] 2 Nov 2026

1, 56127, Pisa, Italy

ARTICLE INFO

ABSTRACT

Keywords: approximate nearest neighbor search vector retrieval quantization graph-based indexing learned sparse retrieval

Approximate Nearest Neighbors (ANN) search over neural embeddings powers modern retrieval systems such as Retrieval-Augmented Generation. Embeddings come in three representations (dense, sparse and multi-vector), and efficient search combines an index structure, a similarity function and a quantization scheme. Existing libraries fix one combination, hard-code each one separately, or factor only one side of the problem, and none covers all three representations. We present V ECTORIUM and K ANN OLO , two open-source Rust libraries that split this design space into a data-management layer and an indexing layer. V ECTORIUM abstracts the vector, the distance function, the quantizer and the dataset independently, and provides 11 encoders, from 32-bit floats to sub-bit codes. On top of it, K ANN OLO implements Hierarchical Navigable Small-World (H NSW ) graph and Inverted Vector File (IVF) indexes that are generic over representation, encoder and metric, and a two-stage index that composes any first-stage index with any rescoring dataset. Since every combination is specialized at compile time, this modularity adds no run-time overhead, and the resulting indexes reach state-ofthe-art performance: our Product Quantization matches FAISS in accuracy at competitive speed, our R AB ITQ and H NSW match or outperform their reference implementations, and our IVF outperforms FAISS. Against six dense and five sparse systems on seven public collections, our indexes give the best memory–throughput trade-off on dense data in 9 of 12 collection–accuracy pairs and perform well on out-of-distribution queries. On sparse data, the systems built on V ECTORIUM lead on all three collections, with SEISMIC and SEISMICWAVE the fastest and KANNOLO’s indexes the most compact.

1. Introduction Vector search is a core component of many Artificial Intelligence (AI) applications. Retrieval-Augmented Generation (RAG), semantic search, recommendation and multimodal retrieval all rely on the same operation: embed a query, then find the items of a collection whose embeddings are closest to it under a similarity function. These systems run on millions or billions of high-dimensional vectors, and at this scale an exhaustive scan is impractical. Approximate Nearest Neighbors (ANN) search therefore trades a small loss in accuracy for orders of magnitude of throughput. Embeddings come in different representations. A dense dual encoder produces a compact real vector per document. A learned sparse encoder produces a high-dimensional vector with only a small number of nonzero coordinates. A multi-vector encoder produces one dense vector per token and scores a document with the MaxSim operator. The three representations offer different trade-offs between memory, effectiveness and efficiency. Which of them is preferable depends on the task or on the data modality: for example, sparse representations are often well suited to text, while dense ∗ Corresponding author

silvio.martinico@isti.cnr.it (S. Martinico); francomaria.nardini@isti.cnr.it (F.M. Nardini); cosimo@linkup.so (C. Rulli); rossano.venturini@unipi.it (R. Venturini) ORCID (s): 0009-0005-7280-6147 (S. Martinico); 0000-0003-3183-334X (F.M. Nardini); 0000-0003-0194-361X (C. Rulli); 0000-0002-9830-3936 (R.

Venturini) 1 Present address: Linkup Technologies, 28 avenue des Pépinières, 94260 Fresnes, France.

S. Martinico et al.: Preprint submitted to Elsevier

ones usually work better for images. The representations can also be combined to improve the quality of the results [1]. Researchers and practitioners therefore need to move freely between the three. They need to explore these trade-offs, design new retrieval techniques, build their RAG or agentic pipelines, and study how the representations interact. Serving any of the three efficiently requires two main ingredients: an index data structure for efficient search, and a quantization scheme that reduces the memory occupied by the vectors. Each ingredient involves further choices: the index structure (e.g. a proximity graph or an inverted index), the similarity function (e.g. Euclidean distance or inner product), the quantization scheme, and the optimizations applied on top of them (e.g. compressing the edges of a graph index). Representation, index, metric and compression are largely independent of one another. Together they span a large cross-product of configurations.

The gap. Libraries differ in how much of that cross-product

a user can reach. Earlier libraries either implement one index over one compression scheme [2, 3], which offers no flexibility, or cover more combinations by hard-coding each of them. FAISS [4, 5], to date the most widely adopted ANN library for dense data, exposes one class per combination of index structure and quantizer. Each class has its own construction and search code path (Section 2.5). The code therefore grows with every combination, which hurts extensibility and maintainability. Recent works have developed more modular frameworks, but each of them factors only one side of the problem. PARLAYANN [6] implements several graph Page 1 of 39

Modular and Efficient Approximate Nearest Neighbors Search with Vectorium and kANNolo

algorithms on shared parallel primitives, and ANNL IB [7] separates the graph algorithm from the data structure that stores the graph. Both are generic over the point type and the distance, but provide few concrete implementations of either: ANNL IB leaves both to the user, and PARLAYANN ships three dense point types with no quantizers. G LASS [8] factors the other side. It selects the quantized representations of dense vectors used to build, search and refine independently of one another, but its graph layout is fixed inside each index. Two limits therefore remain. The first is how far the factoring goes, which determines how easy a library is to understand, maintain and extend. No library abstracts each level of the data side (the vector, the distance function, the quantizer and the dataset) independently of the others and builds its indexes on top of those abstractions. With such a design, a new component would reach every combination without code specific to each of them. In addition, speed is often obtained with distance kernels hand-written for each instruction set, as in hnswlib, FAISS and G LASS. This code is duplicated across instruction sets and grows with every combination it serves. The second limit is the representation of the data. ANN libraries target one representation, in most cases the dense one. Support for two is rare and reduces performance on each [9]. No ANN library covers all three under a unified interface, so an optimization written for one representation does not reach the others.

Our contribution. We present V ECTORIUM and KAN-

NOLO, two Rust libraries that split this design space into two layers. V ECTORIUM (Section 4) is the data-management layer. It defines how a vector is represented in memory, how two vectors are scored under a metric, how a collection is compressed, and how it is stored and accessed. KANNOLO (Section 5) is the indexing layer. It implements standard ANN indexes, H NSW and IVF, each generic over the representations, encoders, metrics and datasets below it. It also composes them into new indexes through a two-stage index. Every component is generic, and the factoring adds no run-time overhead. The compiler specializes each index for the encoding it runs on, so the search loop contains no indirection. A new encoder is therefore available to every index as soon as it is written, and a new index works with every encoding already in the library. The resulting implementations reach state-of-the-art performance: where they implement the same algorithm as a hand-tuned library, they are competitive with it or faster (Sections 4 and 5). In exchange, compilation takes longer and binaries are larger. The distance computations follow the same principle. Each is a few lines of generic code that the compiler vectorizes for the target hardware, with explicit SIMD only for the packed encodings. They remain competitive with hand-optimized implementations (Section 4.4.3). Each component can thus be read, maintained and extended on its own. Finally, V ECTORIUM can be used independently of K ANN OLO . It also underlies the learned sparse retrieval engine SEISMIC [10, 11] and the multi-vector retrieval engine TACHIOM [12], and

S. Martinico et al.: Preprint submitted to Elsevier

it can support the development of any new index technique by providing its data-management layer. The paper is organized as follows. Section 2 covers the background and the related work on the three representations and on the index and compression techniques this work builds on. Section 3 fixes the notation, the evaluation measures, the collections and the measurement conventions. Sections 4 and 5 present V ECTORIUM and KANNOLO respectively, each with the experiments that support its design decisions. Section 6 evaluates the two libraries end to end, and Section 7 concludes.

1.1. Contributions The contributions of this work are the following: • We introduce V ECTORIUM , a data-management library for dense, sparse and multi-vector data. It offers 11 encoders behind one interface, from lossless f32 down to sub-bit codes. They include well-known quantizers that we reimplement, and our own encoders for learned sparse vectors, BlockedSparseEncoder and DotPacking8Encoder [13]. The reimplementations match or beat their reference implementations: our Product Quantizer matches FAISS in accuracy at competitive speed, and our R AB ITQ scan is 1.4–2.4× faster than the original authors’ implementation (Section 4). • We present KANNOLO, extended from its first version [14] (Section 2.5 lists the differences). It provides H NSW and IVF indexes generic over V ECTORIUM ’s representations and encodings, a graph layer shared by all graph indexes, and an experimental study of the optimizations a modern H NSW needs. Its H NSW and IVF match or outperform the reference implementations, hnswlib and FAISS (Section 5). • We present the R ERANKINDEX, KANNOLO’s twostage index, which composes any first-stage index with any rescoring dataset. It was proposed for a sparse-to-multi-vector pipeline [1]. We generalize it into a framework that combines any representations, encodings and indexes. We show that it expresses the two-stage designs behind many recent graph indexes, and that it improves on its first-stage graph index used alone. • An evaluation against six dense and five sparse stateof-the-art systems on seven publicly available collections (Section 6). Two of the sparse systems, SEISMIC and S EISMIC WAVE , are built on V ECTORIUM . On dense data, our indexes offer the best memory– throughput trade-off in 9 of 12 collection–accuracy target pairs. At 99% accuracy our two-stage H NSW reaches 2.5–3.3× the throughput of hnswlib with an index half the size. On the out-of-distribution queries, it is one of only two systems, with L ORANN, to reach 99% accuracy. On sparse data, the systems built on V ECTORIUM lead. On all three collections, SEISMIC Page 2 of 39

Modular and Efficient Approximate Nearest Neighbors Search with Vectorium and kANNolo

and SEISMICWAVE are the fastest systems at every target, and our indexes are the most compact. Our twostage index is the fastest graph index: it is 5.2–48.4× faster than PYANN with a third of its memory. • We release both libraries open source under the MIT license, KANNOLO1 and V ECTORIUM ,2 with Python bindings that give access to dense, sparse and multivector retrieval through a single interface.

2. Background and Related Work

log-saturated term weights, and was later refined for effectiveness and for efficiency [30, 31]. Learned sparse retrieval (LSR) has a second advantage besides that machinery: its representations are rooted in a vocabulary, so they are interpretable.

Multi-vector (late-interaction) embeddings. Pioneered

by COLBERT [32], multi-vector representations keep one dense vector per token. A document is scored by the MaxSim operator, the sum over query tokens of the best inner product against any document token, |𝑞| ∑

This section introduces the vector representations used in retrieval and the quantization schemes that compress them. It then states the approximate search problem and surveys the index families and adjacency-list compression methods that Sections 4 and 5 build on. It closes with the libraries and systems that cover this design space, and with the relation of V ECTORIUM and KANNOLO to the first version of KANNOLO and to our earlier publications.

Every query token interacts with every document token. This gives strong effectiveness but makes exhaustive scoring costly: a document costs |𝑥| vectors of storage and |𝑞| ⋅ |𝑥| inner products of work. COLBERT V2 [33] compresses the token vectors as centroid-plus-residual codes.

2.1. Vector Representations and Retrieval Models

2.2. Vector Quantization

To scale, retrieval systems rely on cheap data representations. These representations are built offline and evaluated efficiently at query time to gather a set of relevant results, which slower but more accurate methods such as cross-encoders can optionally re-rank. Vectors are the most widespread such representation, since they express a set of features compactly. Before representations were learned, they were engineered, from term weighting such as BM25 [15] in text retrieval to local descriptors such as SIFT [16] in computer vision. Neural representations, or embeddings, have superseded them. They fall into three families: dense, sparse and multi-vector. The last is also called late-interaction, after the MaxSim scoring it relies on.

Dense embeddings. A dual encoder built on a pretrained

transformer [17, 18] maps a text into a single real vector of fixed dimension. The vector typically has a few hundred to a few thousand components, all of them populated, and texts are ranked by inner product or cosine similarity. DPR [19] established the architecture for passage retrieval. D RAGON [20] reaches state-of-the-art zero-shot effectiveness through diverse data augmentation, and SNOW FLAKE [21] is a recent open family of general-purpose text encoders. Most ANN work targets this representation.

Learned sparse embeddings. A learned sparse encoder maps a text into a vector over the vocabulary in which only a few dozen to a few hundred coordinates are nonzero. The representation is therefore semantically informed and can still be indexed by the inverted-index machinery of lexical retrieval [22, 23, 24]. uniCOIL [25, 26] learns a weight per term over a document expanded by doc2query [27]. SPLADE [28, 29] adds sparsity-inducing regularization and 1 https://github.com/TusKANNy/kannolo 2 https://github.com/TusKANNy/vectorium

S. Martinico et al.: Preprint submitted to Elsevier

MaxSim(𝑞, 𝑥) =

𝑖=1

max ⟨𝑞𝑖 , 𝑥𝑗 ⟩.

1≤𝑗≤|𝑥|

Vector search usually relies on in-memory indexes. These indexes store the vectors themselves, so the footprint of the vectors is one of the main bottlenecks of the search. Quantization techniques reduce this footprint by mapping the continuous values of a vector onto a smaller discrete set. Scalar quantization casts every component into a narrower numeric type. For instance, it reduces 32-bit floating point numbers to the 16-bit half-precision f16 and bf16 formats, or to a fixed-point integer under a data-dependent range. The simplest uniform quantizers scale every component on a single global range. Research has since moved to more sophisticated techniques that keep as much information as possible in the few bits available. Locally-adaptive Vector Quantization (LVQ) [34] takes the range from the minimum and maximum component of each individual vector, so every vector is quantized against its own bounds. LVQ also adds a residual second level to the resulting codes, used for reranking. R AB ITQ [35] makes binary quantization, a single bit per component, viable. Each vector is centered by the collection’s mean and randomly rotated. This gives a residual 𝑟 = 𝑃 (𝑥 − 𝜇) whose information is spread evenly over the coordinates. The code stores the sign of each component of 𝑟, together with two scalars of metadata. From these, the distance is estimated without bias and with a theoretical error bound. The codes can therefore be used for ranking on their own, and not only as a filter. Extended R AB ITQ [36] generalizes the scheme to 𝑏 bits per component: a sign bit and 𝑏 − 1 magnitude bits. In scalar quantization the distance computation is usually an asymmetric distance computation (ADC), meaning that the query is left unquantized. Quantizing the query as well can reduce the computation to bit operations, as in R AB ITQ. Scalar quantization has its limits at low bit budgets. It spends at least one bit on each component, and the points it can represent lie on a regular grid, which the distribution of Page 3 of 39

Modular and Efficient Approximate Nearest Neighbors Search with Vectorium and kANNolo

the data may not follow. Vector quantization techniques learn a set of codes, the codebook, and encode each vector as a reference to its closest code. Product Quantization (PQ) [37] makes this viable at scale. It splits the 𝑑 dimensions into 𝑚 subspaces of dsub = 𝑑∕𝑚 dimensions and learns a codebook of 2nbits centroids per subspace by 𝑘-means. A vector is encoded as the 𝑚 indices of its nearest centroids, which take 𝑚⋅nbits bits regardless of 𝑑. Scoring is asymmetric: the query is not quantized. Every partial distance the scan can need is precomputed into a per-query table of 𝑚 × 2nbits entries, so a candidate costs 𝑚 lookups and 𝑚 − 1 sums. Refinements rotate the space before splitting it (OPQ [38]), quantize residuals iteratively [39, 40], drop the independence between subspaces [41], or weight the error by its effect on the inner product [42]. The fast-scan kernels of FAISS shrink a subspace to 16 centroids so that a sub-table fits in a vector register [5]. Quantization compresses the values of a vector. A sparse vector, however, is made of two streams: the weights and the integer identifiers of the components they belong to. The identifiers cannot be quantized in a lossy way, since a perturbed identifier names a different component and not a nearby one. Once the weights are quantized down to a byte, the identifiers account for the larger share of the footprint. They are therefore compressed losslessly, with the integer codes of inverted indexes. D OTPACKING8 [13] was introduced to cut this identifier overhead. It bit-packs the gaps between consecutive identifiers and computes the inner product directly on the compressed form, so a compressed scan runs close to the speed of an uncompressed one.

2.3. Approximate Nearest Neighbor Search Vector search is usually framed as 𝑘-Nearest-Neighbor (𝑘-NN) search. Given a collection  = {𝑥1 , … , 𝑥𝑛 } ⊂ ℝ𝑑 of 𝑛 vectors in dimension 𝑑, and given a query 𝑞 and an integer 𝑘, 𝑘-NN search asks for the 𝑘 elements of  that a similarity function ranks highest [43]. Under squaredEuclidean distance the answer minimizes ‖𝑞 − 𝑥‖2 . Under inner product the answer maximizes ⟨𝑞, 𝑥⟩, and the task is also called Maximum Inner Product Search (MIPS). However, 𝑘-NN commonly serves as an umbrella term for any task that asks for the best-matching vectors under a given similarity function. Answering such a query exactly costs a scan of the collection, and the classical space-partitioning structures of Section 2.3.1 do not avoid that cost in high dimension. Approximate Nearest Neighbor (ANN) search therefore trades a limited loss of quality for a large gain in throughput [44, 45, 46]. An index restricts the scan to a subset of , and indexes are compared by the quality/throughput trade-off they offer [47]. Many formulations of ANN exist, distance-bounded and error-bounded among them. The most widely adopted family of ANN techniques, which this paper studies, offers no guarantee at all. Its results are evaluated empirically, by the fraction of the exact top-𝑘 they recover at a given cost (Section 3.2).

S. Martinico et al.: Preprint submitted to Elsevier

2.3.1. Indexes for Dense Retrieval The index families below were developed for dense vectors. The classical ones divide the search space and avoid exhaustive distance computations by pruning irrelevant regions: hierarchical partitioning trees [48, 49, 50, 51, 52] and Locality-Sensitive Hashing (LSH) families [44, 53, 54, 45]. Both struggle with the curse of dimensionality. They degrade towards a linear scan, or they need large memory overheads to hold accuracy.

Inverted Vector File. One partition-based structure

nonetheless remains in wide use, at scale and on accelerators [4]. An Inverted Vector File (IVF) index partitions the collection into nlist clusters, the Voronoi cells [55] of a set of k-means centroids, and stores each cluster contiguously. A query is compared against the centroids, and only the nprobe most promising clusters are scanned. This is the geometric analogue of scanning the posting lists of a query’s terms in an inverted index [56]. IVF is attractive for operational as well as algorithmic reasons: construction is one k-means, and a cluster scan is a dense computation over contiguous memory. Two refinements are common. First, the centroid set can itself be indexed, by a graph or by a second level of clustering, once routing to it becomes expensive. Second, the vectors of a cluster can be stored as centroid-relative residuals 𝑟 = 𝑣 − 𝑐, which a lossy quantizer reconstructs more accurately. This is the IVFADC scheme introduced with Product Quantization [37]. ScaNN [42] trains the quantizer on the error of the score it estimates, where other quantizers use the reconstruction error. SOAR [57] assigns a vector to several clusters.

Proximity graphs. Elsewhere the field has shifted towards

proximity graph indexes, sometimes using a tree to initialize the traversal [58, 59, 60]. Empirical benchmarks place them at the top of the accuracy/throughput trade-off for inmemory search [61, 62, 47]. These indexes represent vectors as nodes and navigate towards a query by best-first traversal. They combine short-range edges for local resolution with long-range edges that let a traversal cross the space. Hierarchical structures such as H NSW [63, 64] and optimized single-layer graphs such as NSG [65, 66] and VAMANA [67] follow this design. Exact 𝑘-NN graph construction requires a prohibitive quadratic number of distance computations, so these indexes build their graphs approximately, from cheap initializations such as NN-Descent [68]. Modern indexes prune neighbors to reduce redundancy and promote longrange edges [69, 70, 64, 67]. They thus favor diversity of direction over proximity alone, and they bound node degree to balance traversal speed against the memory footprint of the adjacency lists. One of the most widely adopted members of the family is the Hierarchical Navigable Small World (H NSW ) graph [64]. It organizes the collection into a hierarchy of proximity graphs. The lowest level, the ground level, holds every vector. Each level above it holds a decreasing random sample Page 4 of 39

Modular and Efficient Approximate Nearest Neighbors Search with Vectorium and kANNolo

of the level below and provides long-range links that shorten navigation. A query descends greedily through the upper levels to locate its neighborhood cheaply. Beam search over the ground level then refines the result. Three parameters govern the trade-off. The degree bound 𝑀 limits how many outgoing edges a node keeps per level (2𝑀 on the ground level). The beam width used to find a new node’s neighbors at insertion is 𝑒𝑓construction , and the beam width at query time is 𝑒𝑓search . Neighbor selection uses a robust pruning heuristic that approximates a relativeneighborhood graph. An edge is kept only if no alreadyselected neighbor is closer to its endpoint than the endpoint is to the node. This diversifies a node’s neighborhood by direction and not by distance alone. H NSW ’s construction is incremental: nodes are inserted one at a time, each running the ordinary search to find its neighbors in the graph. The standard beam search of H NSW halts once the best remaining frontier candidate is worse than the worst element of the current results. Two families of replacement criteria exist. Relaxation-based criteria widen that comparison so the search continues a little past it. Adaptive Beam Search is one of them, with a provable bound on navigable graphs [71]. Saturation-based criteria such as Patience-in-Proximity [72] halt once the result set stops changing. A third line of work learns the decision from the state of the search [73]. Today, graph indexes still represent the state of the art in dense retrieval. However, they rely on scattered memory accesses and do little computation, so they are usually memory-bound. To alleviate this problem, a large part of the literature has moved towards two-stage pipelines [34, 74, 8, 35]. A set of candidates is first generated over quantized vectors and then refined against uncompressed or higherprecision representations. The graph is traversed quickly with the cheap representations, and the larger, more expensive ones are read only for a small subset of candidates. Another line of research goes further and stores heavily compressed representations of the neighbors directly inside the adjacency lists. This replication gives locality in the graph and trades memory for search speed. The two main examples of this line of research are NGT-QG [75, 76] by Yahoo Japan and SYMPHONYQG [77].

2.3.2. Indexes for Learned Sparse Retrieval There are two main families of indexes for learned sparse embeddings: one inherited from lexical retrieval, and one borrowed from the dense techniques above. The first is the inverted index. A learned sparse vector is a vector over the vocabulary, so it can be indexed by the machinery that BM25 and its predecessors built for term-weighted documents (Section 2.1). The index holds one inverted list per component, and a query touches only the lists of its own nonzero components. This structure allows pruning at query time without any geometric index over the collection. Block-Max Pruning [78] and Dynamic Superblock Pruning [79] bound the contribution of a block of postings and skip the block when the bound cannot enter S. Martinico et al.: Preprint submitted to Elsevier

the result set. Lightweight Superblock Pruning (L SP) [80] is a lighter-weight variant of the latter. SEISMIC [10, 11] partitions each inverted list into geometrically cohesive blocks with summary vectors, and evaluates only the blocks whose summary is promising. SINDI [81] is a recent alternative in the same setting. The other family is taken from the dense world. A proximity graph needs only a similarity function that it can evaluate, so the indexes of Section 2.3.1 apply to learned sparse embeddings without changes. PYANN [82], an approach based on H NSW , won the Sparse Track of the NeurIPS 2023 Big ANN challenge [83]. SEISMICWAVE combines the two worlds: it augments SEISMIC with a 𝜅-NN graph over documents [84] to expand the set of results and refine the search.

2.4. Adjacency List Compression In graph indexes the vectors account for only part of the footprint. Adjacency lists carry a significant memory overhead of their own, especially when the vectors are lowdimensional or heavily quantized [85]. Since the order of a list does not matter, it can be sorted and stored as the gaps between consecutive identifiers [86], which a variablelength integer code stores in few bytes. Combined with graph reordering [87, 88], STREAM VBYTE [89] compresses the lists by up to 2.3× with a query time equal to or lower than that of the uncompressed graph [90, 91].

2.5. Libraries and Systems The families above are composable: any index structure can be paired with any quantizer, over any of the three representations of Section 2.1. Libraries differ in how much of this cross-product they cover. FAISS [4, 5] is the most widely used and most extensively tested of them for dense data. It factors the space through a family of loosely related classes, each with its own construction and search code path: IndexIVFFlat and IndexIVFPQ, the scalar-quantizer and ByResidual variants of both, and their counterparts with a graph quantizer. Other libraries narrow the problem: each covers one algorithm thoroughly and leaves the cross-product aside. hnswlib [2], maintained by H NSW ’s authors under the nmslib organization [9], implements one algorithm over uncompressed dense vectors. DiskANN [3] implements VAMANA [67] with a focus on out-of-core search. PARLAYANN [6] provides deterministic parallel builds for several graph algorithms. ANNL IB [7] builds on the same primitives and separates the graph algorithm from the data structure that stores the graph. It composes filtered search, deletions and snapshots as modules over them, for dense vectors under a point type and distance function supplied by the user. Finally, the entries that populate public benchmarks [92, 82] are typically specialized to one algorithm and one representation. Vector databases and search engines form a different layer. Milvus [93], Qdrant [94], Vespa [95] and Elasticsearch [96] serve dense and learned sparse embeddings behind one interface, and Weaviate [97] serves dense embeddings. Qdrant, Vespa, Elasticsearch and Weaviate also Page 5 of 39

Modular and Efficient Approximate Nearest Neighbors Search with Vectorium and kANNolo

score late-interaction token vectors with MaxSim, mostly as a reranking stage. They are services built for deployment, with persistence, filtering and distribution, and each representation is served by a separate index type. The libraries above, and ours, form the layer where indexes are designed, composed and compared. K ANN OLO was introduced at ECIR 2025 [14] as a modular Rust library offering H NSW over both dense and sparse representations. This paper presents the system after that first release, with large structural differences. The data-management layer is now a separate library, V ECTORIUM (Section 4). It owns representations, distance kernels, encoders, storage, clustering and the index contract. K ANN OLO (Section 5) implements H NSW , IVF and the two-stage R ERANKINDEX generically over the representations and encodings that V ECTORIUM supplies. Beyond the ECIR version, this paper therefore adds the layer separation and its cost model, the encoder family and its quantization schemes, the multi-vector representation, the IVF family, and the graph backend abstraction with its optimizations and compressed adjacency layouts. It also adds the measurements that motivate each choice. Three components have been published separately and are summarized here: the R ERANKINDEX with its two pruning rules [1], generalized here to any pair of indexes and datasets; the D OTPACKING8 encoding and the JUMPD OT kernel [13], which V ECTORIUM integrates; and the graph compression scheme of [90], developed in KANNOLO.

3. Preliminaries and Experimental Setup This section fixes the notation, the measures, the collections and the measurement conventions used by Sections 4 and 5. These conventions apply throughout those two sections unless a study states otherwise.

3.1. Notation We keep the notation of Section 2.3: a collection  holds 𝑛 vectors of dimension 𝑑, each identified by its position 𝑖, and 𝑞 is a query. For a proximity graph, 𝑁 is its number of nodes (𝑁 = 𝑛 at the ground level), 𝑀 its maximum degree and 𝑒𝑓search the beam width of its search. 𝑚 is the number of subspaces of a Product Quantizer, 𝜆 the earlytermination relaxation factor, and 𝛽 the early-exit patience of the reranking stage.

3.2. Evaluation Measures Accuracy@𝑘. The percentage of the exact top-𝑘 that a

search recovers over a set of queries ,

𝑞

𝑞

1 ∑ | 𝑘 ∩ 𝑘 | , Accuracy@𝑘 ≔ 100 ⋅ || 𝑞∈ 𝑘 where 𝑞𝑘 is the set the index returns for query 𝑞 and 𝑘𝑞 is the exact top-𝑘 set.

Throughput. Queries per second (QPS), measured in wall-

clock time over the whole query set. Unless a study says S. Martinico et al.: Preprint submitted to Elsevier

Table 1 Development collections; 𝓁2 /IP denotes squared-Euclidean / inner product. ⋆ the average number of non-zero components is ≈ 120. Collection

Type

Sift1M dense Sift1B dense Dragon dense SnowFlake dense CoCondenser sparse

𝑑 Metric 128 128 768 1 024 ⋆ 30 522

𝓁2 𝓁2 IP IP IP

𝑛 Queries 1M 1B 8.8M 8.8M 8.8M

10 000 — 6 980 6 980 6 980

otherwise, queries are issued one at a time, and the reported figure is the reciprocal of the mean per-query latency.

Footprint. Bytes per vector, as computed by the SpaceUsage

accounting of V ECTORIUM . This accounting traverses a structure and sums the bytes it owns. When owned bytes cannot be compared with those of an external system, we report peak resident memory instead and state this in the text.

Build time. Wall-clock seconds for the complete construc-

tion, excluding serialization on both sides of any comparison.

3.3. Development Collections Table 1 lists the collections used for the design studies of Sections 4 and 5. They span two representations, three dense dimensionalities and two metrics, so each design decision is tested in the regimes where it can behave differently. SIFT1M [16] contains hand-crafted, integer-valued features compared under squared-Euclidean distance. Its subsampled billion-scale counterpart SIFT1B provides the 10k– 100M range of sizes used in the scaling study of Section 5.3.1. The three remaining collections are all built on the passage corpus of M S M ARCO-V1 [98] and differ only in the encoder used to embed it: D RAGON [20] and SNOW FLAKE [21] are dense encoders, and SPLADE COCONDENSER [28, 99] is a sparse encoder. After its first mention, each collection is referred to by the name of its encoder alone.

3.4. Hardware and Software All measurements run on a dual-socket Intel Xeon Silver 4314 (Ice Lake-SP, 2.40 GHz) with 16 cores and 32 hardware threads per socket, for 64 threads in total. Each socket has 24 MiB of shared L3 cache. The machine has two NUMA nodes and 503 GiB of memory, and runs Ubuntu 24.04 LTS. The processor supports the AVX-512 instructionset extensions, including AVX-512VPOPCNTDQ. The portability measurements of Appendix C.1 also run on this machine with a build restricted to AVX2. They also run on a Raspberry Pi 5 (Cortex-A76, Armv8.2 with 128-bit NEON, 16 GiB), built from the same source with the same nightly toolchain. V ECTORIUM and KANNOLO are built with the Rust nightly toolchain of 2026-06-19 (rustc 1.98.0-nightly), which is required for the portable-SIMD and algebraic floating-point intrinsics of Section 4.3. They are compiled Page 6 of 39

Modular and Efficient Approximate Nearest Neighbors Search with Vectorium and kANNolo

in release mode with target-cpu=native. FAISS appears as a baseline in two versions: 1.15.0 for the ProductQuantization comparison of Section 4.4.3, and 1.14.3 for the IVF comparison of Section 5.5. hnswlib appears at release 0.9.0 in the validation of Section 5.4.2. The end-to-end comparison of Section 6 uses release 0.8.0, the one pinned by V IBE.

3.5. Measurement Conventions Parallel scans. The exhaustive scans of Section 4 run on

all 64 threads with numactl --interleave=all, on a machine verified to be idle. When an external system is compared, both sides run under the same settings.

Micro-benchmarks. Micro-benchmarks

run singlethreaded, pinned to one core. We report the minimum over five repetitions. Hardware performance counters (retired instructions, cycles, last-level-cache misses) are restricted to the measured loop, so instruction counts and IPC measure the kernel alone and exclude the harness.

Search timing. We time one query at a time. The refer-

ence implementations we compare against parallelize across queries without shared state between them, and so do ours. A batched figure would therefore measure thread scaling and not the index itself. Every index class of KANNOLO’s Python interface exposes this as a batch search with a configurable number of threads.

Ground truth. The exact top-𝑘 under the collection’s own

metric, computed in f32. It is either shipped with the collection or computed once and reused across every configuration.

4. V ECTORIUM V ECTORIUM is the data management layer beneath every index in the system. It owns how a vector is represented in memory, how two vectors are scored under a given metric, how a collection is compressed, and how that collection is stored and accessed. One set of storage and scoring primitives is thus shared across systems and index families, and no index family reimplements it. V ECTORIUM also defines the Index contract that the indexes above it implement, and ships the exact flat index; Section 5.1 describes both together with KANNOLO’s indexes.

4.1. Overview and Design Principles V ECTORIUM is built from four core abstractions, each layered on top of the previous one. A vector view is data in flight. It is a lightweight, zero-copy view over a vector’s storage and describes how a given vector type is laid out in memory (Section 4.2). A Distance is a score together with its ordering. It keeps all generic frontier and heap code independent of the metric (Section 4.3). A VectorEncoder, paired with its QueryEvaluator, unifies compression and scoring. It maps input vectors to a stored encoding and defines how a query is scored against it, so all distance logic is kept in one S. Martinico et al.: Preprint submitted to Elsevier

place (Section 4.4). A Dataset is storage, parameterized by an encoder. It holds encoded vectors and hands out views of them (Section 4.5). Figure 1 shows how the four are stacked and how the three orthogonal axes of layout, metric and compression cut across them. Three principles govern their design. First, encoders and datasets are type parameters, so every combination of layout, metric and compression is a distinct monomorphized type. Distance kernels inline into the caller, no dynamic dispatch appears on the hot path, and a single generic implementation, such as the H NSW index of Section 5.4, serves every encoding. Second, storage layout (dense, sparse, packed, multi-vector), metric (squared-Euclidean or inner product) and compression scheme vary independently, and any meaningful combination is obtained by naming its type, with no new code. Third, storage is decoupled from scoring. A Dataset holds no distance logic, and all scoring is done by the encoder’s evaluator, so one dataset and one index search loop serve every compression scheme.

The cost of genericity. Monomorphization moves the cost

of the factoring from run time to build time. The Python bindings of KANNOLO 0.15.1 instantiate every combination they expose. A clean release build of these bindings takes 113 s on our machine, 83 s of it for KANNOLO itself, and produces a 5.3 MiB library. The beam search of Section 5.3 is compiled 164 times, once per combination of graph layout, encoder and metric. These copies make up 3.1% of the generated code. Each copy inlines its distance kernel into the search loop, so scoring makes no indirect call. Calling the kernel through a trait object would make each distance 1.07– 1.52× slower on cache-resident vectors and up to 1.07× slower under random access (harness of Table C.1). A Rust program compiles only the combinations it uses. The Python bindings expose a fixed set of index classes, so their number of copies grows with the product of the graph layouts, encoders and metrics they cover. A combination outside that set needs a new binding class and a rebuild. Table 2 summarizes the concrete encoders that V ECTORIUM ships and the position of each on the three axes. Its bits/dim column also maps the compression spectrum. The range goes from the lossless f32 cast at 32 bits/dim down to one bit per dimension for the binary code of R AB ITQ, and below one bit for a sufficiently short PQ code.

4.2. Vector Representations The lowest layer of V ECTORIUM fixes how a vector is held in memory while it is encoded, stored or scored, independently of any encoder or metric. Every representation is a zero-copy view over its backing storage, cheap enough to pass by value in the innermost distance loops. There are four. A dense view covers one contiguous, fixed-dimension run of values. A sparse view covers two parallel runs: the component indices and their nonzero values. A packed view covers an opaque run of words that only the encoder that produced it can interpret, such as the bit-packed blob of the D OTPACKING8 encoder (Section 4.4.6). A multi-vector view covers a flat buffer of token vectors that share one dimension. Page 7 of 39

Modular and Efficient Approximate Nearest Neighbors Search with Vectorium and kANNolo kANNolo

Index

FlatIndex · HNSW · IVF · RerankIndex fetch by id or range

vector view

Distance

score(view)

Vectorium Training

learned encoder parameters

sample

Dataset<E>

stores encoded vectors, hands out views; no distance logic

encode

VectorEncoder

layout × metric × compression input · query · encoded types resolved at compile time

builds

QueryEvaluator

built once per query, scores many candidates; all scoring lives here

calls per-query state

Distance kernels

₂

ℓ · IP · MaxSim, one per (metric, layout) pair

parameters

input vectors

build time query time

query

Figure 1: The layered architecture of Vectorium and its interface to the indexes of kANNolo.

Table 2 Capability matrix of Vectorium’s concrete encoders. 𝑑 is the vector dimension, 𝑚 the number of PQ subspaces (one u8 code each), 𝑏 = total_bits ∈ {2, 4, 8} the RaBitQ-ext code width. Encoder

Layout

Stored type

Metric

Bits/dim

PlainDenseQuantizer ScalarDenseQuantizer ProductQuantizer RabitqQuantizer RabitqExtQuantizer ScalarSparseQuantizer BlockedSparseEncoder DotPacking8Encoder ScalarMultiVecQuantizer MultiVecProductQuantizer MultiVecTwoLevelProductQuantizer

dense dense dense dense dense sparse packed sparse packed sparse multi-vector multi-vector multi-vector

f32 f16, bf16, 16/8-bit fixed point u8 codes u64 sign code u64 𝑏-bit fields f32, f16, bf16, 16/8-bit fixed point u16/f16 blocks packed gaps; fixed-point, f16 or scalar values f32, f16, bf16, 16/8-bit fixed point u8 codes centroid + u8

𝓁2 /IP 𝓁2 /IP 𝓁2 /IP 𝓁2 /IP 𝓁2 /IP 𝓁2 /IP IP IP MaxSim MaxSim MaxSim

32 16/8 8𝑚∕𝑑 1 + 64∕𝑑 𝑏 + 64∕𝑑 32/16/8 16 var. 32/16/8 8𝑚∕𝑑 (8𝑚+32)∕𝑑

Each view is generic over its value type, so one encoder body covers f32, f16, bf16 and fixed-point values alike. An encoder (Section 4.4.1) names three representations separately: the input pushed in at build time, typically plain f32; the query; and the encoded vector stored after compression. Because the three are distinct, a query can be scored against a compressed stored vector without either side matching the other’s representation. An asymmetric scoring scheme is thus expressed in the type system and does not depend on a convention.

4.3. Distance Computation The next layer turns two vector representations into a score. A single Distance trait fixes what a score is and how scores are ordered. A flat family of kernels computes the scores, each operating directly on the slices behind a view (Section 4.2). The encoders’ evaluators invoke these kernels (Section 4.4.1). A distance is a thin wrapper around a single f32. The Distance contract only requires it to be totally ordered, and SquaredEuclideanDistance and DotProduct implement it. DotProduct reverses its ordering, because a larger inner product denotes a better match, while for squared-Euclidean distance a smaller value does. “Smaller is better” therefore holds under either metric. Every generic consumer above S. Martinico et al.: Preprint submitted to Elsevier

this layer, including the candidate heaps, is written once and never branches on the metric in use.

A unified relaxation predicate. Early termination and pruning ask whether a candidate is within a relaxed bound, whose direction depends on the metric. The Distance trait hides the direction behind a single predicate. A minimize metric accepts a candidate when 𝑑 ≤ threshold+|threshold|⋅ 𝜆, a maximize metric when 𝑑 ≥ threshold − |threshold| ⋅ 𝜆, and 𝜆 = 0 recovers the exact comparison. The threshold enters in absolute value because an estimated distance can be negative, as the unbiased estimator of R AB ITQ returns for near-duplicate pairs. On such a bound, the multiplicative form threshold ⋅ (1 + 𝜆) would tighten the window when it should widen it. The kernel family. The scoring work is done by a small

set of free functions, one per (metric, layout) pair that occurs. Each comes in a checked form that validates its invariants and an unchecked form that assumes them. Sparse data is scored against a densified query: the evaluator expands the sparse query once into a dense buffer, into which each candidate’s stored components index directly. Sparse-versus-sparse comparisons, needed only during construction, use a sorted merge. The Page 8 of 39

Modular and Efficient Approximate Nearest Neighbors Search with Vectorium and kANNolo

sparse layouts compute squared-Euclidean distance as ‖𝑞‖2 + ‖𝑣‖2 − 2⟨𝑞, 𝑣⟩, so they visit only the intersection of the active components of the two operands and not their union. The same layer hosts the maxsim kernel of the multi-vector encoders (Section 4.4.7), which decodes each document token to f32 once and reuses it across all query tokens.

Portable auto-vectorization. Explicit SIMD in V ECTO-

RIUM is confined to the packed encodings, described at the

end of this paragraph. The kernel family, which scores every other encoding, contains none. Hand-written vector code would have to be duplicated for each instruction set and again for every (metric, layout) pair, and each copy would fix the register width and the reduction schedule in the source. Each kernel is a few lines of plain Rust. The compiler, with native CPU optimizations (Section 3.4), lowers it to the vector ISA of the target. Two choices make this possible. First, the reductions use Rust’s algebraic (fast-math) floating-point operators, which let the compiler reassociate and contract the multiply–add chains. A strictly ordered IEEE reduction fixes the sum to one scalar schedule, which can be neither widened into SIMD lanes nor fused. This relaxation applies only to scoring: the ordering and relaxation predicates above compare exactly. Second, the unchecked kernels present loops of known, equal trip count, with no bounds checks that would prevent widening. The packed encodings are scored differently. R AB ITQ (Section 4.4.4) scores by popcounts and shift-and-mask unpacking, and D OTPACKING8 (Section 4.4.6) decodes bit-packed gaps. These kernels use Rust’s portable SIMD, which fixes the lane width but is still independent of the ISA. Only the packed sparse encodings (Sections 4.4.5 and 4.4.6) call per-ISA intrinsics (AVX2, and NEON on Armv8). Without fast-math the kernels do not vectorize at all. The strict variant lowers to a scalar reduction over a single accumulator, and the algebraic one to packed fused multiply–adds over several independent accumulators. The algebraic f32 dot-product and squared-Euclidean kernels are 13× faster at 𝑑=128 and 18–20× faster at higher dimensions, and exhaustive top-10 results stay identical on SIFT and D RAGON. The same source vectorizes unchanged when rebuilt for AVX2 alone and on Armv8 with NEON (Appendix C.1). Expert-tuned intrinsics may still outperform compiled code; Section 4.4.3 measures this difference against FAISS.

Batched evaluation. Every kernel that can exploit

batching also comes in a batch-of-six form, exposed on the QueryEvaluator contract (Section 4.4.1). It scores one query against six vectors in a single fused pass. The gain comes from query reuse: each query element is loaded once and used six times. On cache-resident f32 vectors the batched kernel is 1.32–1.73× faster than six single-candidate calls. Under random access, as among graph neighbors, it is 1.24–1.59× faster for f32 and f16 (Appendix C.1). The contract provides the batched call as a default fan-out, and S. Martinico et al.: Preprint submitted to Elsevier

an encoder overrides it only where a single-pass kernel is faster.

4.4. Encoders and Quantization A VectorEncoder binds a compression scheme to a layout and a metric. It defines how an input vector is mapped to its stored encoding and how a query is scored against that encoding through the kernels of the previous layer. It is the single unit in which storage, compression and scoring meet. The subsubsections that follow state the contract and then present the concrete encoders along the compression axis, from lossless to increasingly lossy.

4.4.1. The Encoder and Evaluator Abstraction The VectorEncoder is the central contract of V ECTORIUM . An encoder supplies everything that a dataset or index needs about an encoding, so the dataset or index never needs to know which encoding it is. This comprises its distance (Section 4.3), its three representations (Section 4.2), the evaluator that computes scores, and the query-side parameters used to build that evaluator. All of them are resolved at compile time.

The evaluator. All scoring goes through a QueryEvaluator. It is created once from a query and then scores many candidates, so the once-per-query work (e.g. building PQ’s distance tables) is amortized over all the distances it computes. Evaluators are constructed in two ways: from a user query at search time, and from a stored dataset vector at build time. In the second case the encoded vector is first decoded back to f32. One scoring machinery thus serves both search and index construction. Query-side parameters. An evaluator may take a param-

eter that belongs to the query and not to the stored data. For example, for the 1-bit R AB ITQ encoder (Section 4.4.4) it sets how finely the query is quantized. Because this parameter is not part of the stored encoding, one stored dataset can serve every width concurrently from parallel threads.

Encoding, decoding and stored-vector pairs. Encoding

appends directly into the dataset’s own backing storage, with no intermediate allocation. Decoding reconstructs an f32 vector for the build-time workflows that need the original space. Graph construction (Section 5) also scores pairs of stored vectors against each other many times. The contract therefore includes a pairwise method, which encoders override to score their stored form directly, without decoding either operand.

4.4.2. Scalar Quantization Scalar quantization (Section 2.2) populates the compression axis by changing the numeric type in which values are stored. ScalarDenseQuantizer is a component-wise cast from an input value type to a stored value type, generic over the metric. ScalarSparseQuantizer is the same cast over the sparse layout. Neither has a codebook or a learned range.

Page 9 of 39

Modular and Efficient Approximate Nearest Neighbors Search with Vectorium and kANNolo Table 3 Precision against accuracy and footprint. Exhaustive FlatIndex top-10 search against the exact f32 ground truth; throughput parallel across queries on 64 threads (Section 3.5). Dataset

Prec. B/vec Accuracy@10 (%)

Sift1M (128-d, 𝓁2 )

f32 f16 bf16

512 256 256

100.00 819 100.00 1 401 100.00 1 109

Dragon (768-d, IP)

f32 f16 bf16

3 072 1 536 1 536

100.00 99.94 99.59

10 22 20

CoCondenser f32 (sparse, IP) f16

728 488 488

100.00 99.92 99.23

24 26 25

bf16

QPS

When the two value types coincide, the cast is the identity and the quantizer is lossless; the PlainDenseQuantizer and PlainSparseQuantizer aliases name this case. When the stored type is narrower, the same code path yields a reducedprecision encoding. The stored type may be any value type that can be converted from f32: f32 itself, the half-precision f16 and bf16, or the fixed-point FixedU8Q and FixedU16Q. Whatever the stored precision, the query is kept in f32. f16 is the default storage precision throughout the system, on the assumption that it costs no measurable accuracy. We measure accuracy, footprint and exhaustive search throughput as a function of stored precision. To do so, we run a full top-10 scan against the exact ground truth on three collections with three value distributions: SIFT1M (integers), D RAGON (real numbers) and COCONDENSER (real, positive numbers). Table 3 supports this choice. f16 ranks like f32: Accuracy@10 is unchanged on SIFT1M and drops by 0.06 points on D RAGON and 0.08 on the sparse set. These drops are far below any difference that an approximate index later introduces. f16 also halves the dense footprint and, because the scan is memory-bandwidth-bound, roughly doubles throughput. bf16 loses an order of magnitude more accuracy than f16 on D RAGON and COCONDENSER. Its 7 mantissa bits, against 10 for f16, resolve values more coarsely. On SIFT1M the two tie, because the small-integer SIFT values are representable on the coarser grid of bf16. Finally, the sparse row shows the asymmetry of a compression that touches only the values. It removes only about a third of the per-nonzero payload, because each entry still carries an uncompressed u16 component index.

4.4.3. Product Quantization Product Quantization (Section 2.2) is the first encoder in this section with a learned codebook. ProductQuantizer implements it with KSUB = 256 centroids per subspace, so a code is one u8 and a sub-table takes 1 KiB. Narrower codes would have to be bit-unpacked on every lookup. Wider ones would inflate both the codebook and the per-query table with S. Martinico et al.: Preprint submitted to Elsevier

no clear gain. The encoded form is therefore 𝑚 bytes, and the footprint is 8𝑚∕𝑑 bits per dimension. The subspace count is a compile-time parameter, so the code length is part of the encoder’s type. The per-vector stride of a PQ dataset is therefore a compile-time constant. The divisibility invariants that the kernels rely on (𝑚 mod 4 = 0 and 𝑑 mod 𝑚 = 0) are asserted once at construction and assumed in all later code. The 𝑚 codebooks are 𝑚 independent clustering problems, each a k-means over dsub-dimensional projections. Training therefore parallelizes across subspaces with no synchronization, as one task per subspace over the k-means primitive of Section 4.6. The evaluator builds the 𝑚 × 256 lookup table once per query. A small ProductQuantizerDistance trait, implemented once per metric, selects the table builder at compile time, so both metrics share one scoring path. Building the table costs 𝑚 ⋅ 256 ⋅ dsub = 256𝑑 multiply–adds, roughly the cost of 256 uncompressed full-dimension distances. Each subsequent candidate costs 𝑚 loads and 𝑚 − 1 adds. A PQ evaluator therefore recovers the cost of its table after a few hundred candidates, which is the regime in which index search uses it.

Centroid layout. K-means naturally produces centroids as

an array of structures (AoS), with the dsub values of centroid 𝑖 contiguous. With that layout, the table build needs one horizontal reduction per table entry. At a typical dsub of 8 or 16, this final horizontal sum costs more than the multiply–adds before it, and it is repeated 256 times per subspace. V ECTORIUM stores the codebooks transposed, as a flat structure of arrays (SoA) of shape [𝑚 × dsub × KSUB]. At dimension 𝑘, the values of consecutive centroids are then contiguous and load in one vector instruction. The loop nest is inverted, with dimensions in the outer loop and centroids in the inner loop. The accumulation becomes vertical: partial distances accumulate in parallel lanes, with no horizontal reduction anywhere. The kernel tiles the centroids into groups of 𝑇 =32 whose accumulators stay in registers across all dimension passes, so the table is written once, at the end. This is the same query reuse as in the batched kernels: here a query element is reused across 𝑇 centroid lanes, there across six candidates. As with the kernels of Section 4.3, the loop body is plain Rust and the compiler performs the vectorization. Without tiling, each of the dsub passes would read and write the running partial of every centroid, ≈2 ⋅ dsub ⋅ KSUB table accesses per subspace, against KSUB writes with tiling. Timed alone on the three deployed shapes (Appendix C.2), the shipped kernel builds the table 4.5–5.2× faster than the natural AoS build and 3.9–4.4× faster than an AoS build with fast-math. Almost all of the gain comes from the layout.

Scoring. A candidate costs 𝑚 table lookups, one per code

byte, which are summed. Their addresses depend on the data, so they cannot be vectorized. Scoring consumes the codes four at a time into four independent accumulators. PQ keeps the default batched call of Section 4.3. The query has already Page 10 of 39

Modular and Efficient Approximate Nearest Neighbors Search with Vectorium and kANNolo

been consumed into the table, so there is nothing to reuse, and six codes index a 1 KiB sub-table spanning sixteen cache lines at six unrelated positions.

1.25× faster than ours. From 32 MB to 850 MB the two are within 3.4%, and at 1.7 GB the AVX-512 kernel stays 1.13× ahead.

Validating against a reference implementation. We

4.4.4. R AB ITQ: Binary and Multi-Bit Codes R AB ITQ quantizes every component individually, so its minimum code width is one bit per component: 𝑑∕8+8 bytes once the metadata word is counted. V ECTORIUM ships it as two encoders that share their preprocessing, metadata layout and metric handling. RabitqQuantizer is the 1-bit case, and RabitqExtQuantizer [36] widens the code to total_bits ∈ {2, 4, 8} bits per component.

check our PQ encoder end to end against FAISS’s reference implementation, IndexPQ. Both sides train on the same onemillion-vector sample with the same parameters: KSUB= 256, 25 k-means iterations, one restart, non-spherical, exhaustive assignment. We disable the subsampling of FAISS’s max_points_per_centroid, which costs it 0.6 accuracy points. We link FAISS’s dynamic-dispatch build, which its documentation recommends for deployment, and both sides search one query at a time (Section 3.5). Code lengths are {0.5, 1, 2} bits per component, i.e. 𝑚 ∈ {8, 16, 32} on SIFT1M and 𝑚 ∈ {48, 96, 192} on D RAGON. Table 4 reports the result. Accuracy matches. No configuration differs by more than 0.43 points, which is about the spread of PQ’s own retraining (up to 0.40 points across three trainings, Table 6). The sign of the difference changes across the six configurations, and footprint is identical to three decimals. Search is faster for us when the codes fit in cache, and faster for this FAISS build when they do not. Our advantage is 1.62, 1.30 and 1.20× as the SIFT1M code array grows from 8 to 32 MB. On D RAGON the two are within a few percent at 425 and 850 MB, and at 1.7 GB FAISS is 1.17× faster. Our advantage holds while the code array fits in the 48 MiB of L3 that the two sockets of the machine provide together (Section 3.4), since the parallel scan runs on both. Once the scan is memory-bandwidth-bound, the measurement no longer reflects the lookup kernel. Training is faster in FAISS, by 1.13–1.57×, and the gap is widest where dsub is largest. Its assignment step scores all 256 centroids in one call, which fills a vector register exactly at dsub=16; ours iterates over the centroids one at a time. Encoding changes sharply with dsub. FAISS is 2.8–4.1× faster than us at dsub ∈ {4, 8}, but 16–17× slower at dsub=16. At that width its encoder switches to an algorithm that writes out the full 256𝑚-entry distance table of every vector before reducing any of them. Both results follow from how FAISS implements PQ. Its 8-bit scan is hand-written once per instruction set, in 1 146 lines across eleven files. Its table builds dispatch to functions specialized at dsub ∈ {1, 2, 4, 8}. This code also serves nbits = 16 and arbitrary widths, while we fix the width at 8. Its 736 lines of 8-bit SIMD code are directly comparable to ours, which consist of one 16-line function with no intrinsics plus 45 lines for the two table builds. The specialized code makes FAISS up to 17% faster on the largest D RAGON scan. However, each hard-coded threshold can select the wrong kernel, as the 16× slower encoder at dsub=16 shows. On this machine the per-instruction-set dispatch also selects a slower kernel. The AVX-512 scan loads its table through a microcoded gather instruction. It runs 1.2–2.0× slower than FAISS’s own scalar loop on SIFT1M, though up to 1.18× faster on the largest D RAGON codes. Forced to its scalar loop, FAISS scans the two smallest SIFT1M code arrays 1.23 and S. Martinico et al.: Preprint submitted to Elsevier

Rotation and stored metadata. V ECTORIUM uses an

𝑂(𝑑 log 𝑑) rotation of four rounds. Each round applies a random ±1 sign flip, a fast Hadamard transform over the largest power-of-two prefix, and a Kac’s-walk butterfly [100] that mixes the remaining tail. It is seeded by the caller for reproducibility. A 1-bit code stores the sign of each component of 𝑟 = 𝑃 (𝑥 − 𝜇), sixty-four to a u64 word. A 𝑏-bit code adds 𝑏 − 1 magnitude bits per component. The estimator needs two per-document constants: ‖𝑟‖2 , and ‖𝑟‖ already divided by the rescale factor that maps magnitudes onto code levels. In both encoders they are bit-cast into the two halves of one extra u64 word, which keeps a division off the per-candidate path. The rescale factor is estimated once at training time by averaging over random unit vectors. The encoding does not depend on the metric, which is a query-time property. Under squared-Euclidean distance the query is centered by the same mean. Under inner product it is only rotated, and the term ⟨𝜇, 𝑞⟩ is added back in closed form, with one dot product per query.

Code layout. All 𝑏 bits of a component are stored in one

aligned bit-field of a single component-major slab of 𝑏 ⋅ 𝑑∕8 + 8 bytes. A scan therefore reads each query element only once. It spends one widening conversion and one fused multiply–add per 16 components, the f32 lane count of a 512bit vector register. That count is the same at every width; only the number of shift-and-mask extractions changes. For this reason the scan cost grows much more slowly with the code width than in a table-lookup scheme (Table 6). Only the aligned widths {2, 4, 8} are supported, which together with 1 bit cover the useful range densely.

Query quantization. The 1-bit encoder scalar-quantizes

the query to query_bits ∈ {1, 2, 4, 8} bits per component and scores symmetrically, as popcounts over bit planes. query_bits is passed at search time, so one stored index serves every setting concurrently. We sweep it over one stored index on each of SIFT1M, D RAGON and SNOW FLAKE. Accuracy saturates at 4 bits. Going from 1 to 2 bits gains up to +9.3 points of Accuracy@10 (from 53.9 to 63.2 on D RAGON). Going from 4 to 8 gains at most +0.36 points and costs up to 20% of the throughput (148 against 118 QPS on D RAGON). V ECTORIUM therefore defaults to 4, and the comparisons below use this setting. The extended encoder does not quantize the query at all. Page 11 of 39

Modular and Efficient Approximate Nearest Neighbors Search with Vectorium and kANNolo Table 4 Validation of ProductQuantizer against Faiss IndexPQ (dynamic-dispatch build, plain ADC path) at nbits = 8 on both sides. Exhaustive top-10 scan over all queries, conventions of Section 3.5 on both sides, minimum of three repetitions. Footprint is identical throughout and omitted; bold marks a difference above 5%. Accuracy@10 (%)

Train (s)

Encode (s)

QPS

Collection

𝑚

ours

Faiss

ours

Faiss

ours

Faiss

ours

Faiss

Sift1M (𝑑=128, 𝓁2 )

8 16 32

32.1 55.3 71.6

31.7 55.6 71.5

9.4 13.6 24.3

6.0 11.0 21.4

0.20 0.31 0.55

3.31 0.11 0.13

7 587 4 125 2 155

4 676 3 172 1 802

Dragon (𝑑=768, IP)

48 96 192

53.8 68.6 83.2

53.8 69.0 83.3

55.3 81.1 144.0

37.2 65.0 126.9

10.4 16.6 29.3

172.3 5.9 7.2

140.8 64.7 26.7

136.0 65.9 31.2

This gave better accuracy and better throughput, because a quantized query costs scan time in proportion to its width.

Validating the implementation against the reference. We test our R AB ITQ implementation against the authors’ own. The reference library ships IVF and H NSW indexes with an encoder tailored to them. To isolate quantizer, rotator and estimator, we add a flat scan to it. Both sides then run the library’s production scan path, and only the encoder differs. The comparison runs on SIFT1M and D RAGON, with all queries and top-10, under the parallel-scan conventions of Section 3.5 on both sides. At total_bits = 1 both sides quantize the query to 4 bits, the width that the reference fixes at compile time. At wider codes neither side quantizes the query. Every row therefore compares the same estimator on the same query representation. The reference seeds its rotation from the system entropy source, so its figures are the mean of three repetitions. Their spread, up to 1.1 points on SIFT1M and 0.32 on D RAGON, limits how precisely the differences below can be interpreted. Table 5 reports the outcome. Accuracy matches. The largest difference at any width is 0.80 points on SIFT1M and 0.29 on D RAGON, with no systematic sign, and neither exceeds the run-to-run spread of the reference on that collection. Our footprint is smaller at every width, by 8 bytes/vector at 1 bit and 16 at wider codes. The reference stores more per-vector metadata beside each code. This includes an error bound that lets its graph index prune on a cheap lower bound before computing a full distance; we do not have this first pass. Our throughput is higher, by 1.44– 1.74× on SIFT1M and 1.44–2.37× on D RAGON, and it comes from the layout, width and accumulator choices described above. Build time is not compared, because centering on one global mean, as our encoders do, had to be added to the reference in our harness.

Comparing R AB ITQ and Product Quantization. Both encoders are available from one bit per component upwards, so a user must choose between them in this range. The R AB ITQ papers compare it with PQ and its variants at equal code length, and find its distance estimates more accurate and faster to compute [35, 36]. We compare the two as the S. Martinico et al.: Preprint submitted to Elsevier

encoder of an exhaustive scan on three dense collections of increasing dimensionality (Table 1), with codes of the same size. We report accuracy, throughput and build time from harnesses that differ only in the encoder. Footprints are matched by setting 𝑚 = total_bits ⋅ 𝑑∕8. This leaves each R AB ITQ configuration 8 bytes/vector larger, because of its metadata word. PQ trains on a one-million-vector sample, and the 1-bit R AB ITQ rows use a 4-bit query. We follow the parallel-scan conventions of Section 3.5 with one repetition, except for PQ’s accuracy, which is the mean of three trainings. Table 6 reports the outcome. The better encoder depends on both the data distribution and 𝑑, and the balance moves monotonically towards R AB ITQ as 𝑑 grows. At the lowest bit rate, the accuracy difference Δ = PQ − R AB ITQ is +29.5 points at 𝑑 = 128, +1.9 at 𝑑 = 768 and −2.3 at 𝑑 = 1 024. Within each collection it roughly shrinks as the code widens. All differences but one, D RAGON at 4 bits, are larger than the largest spread of PQ’s three trainings. On SIFT1M, whose hand-crafted features are strongly clustered, PQ is more accurate at every width. On D RAGON it leads only up to 4 bits, and on SNOW FLAKE R AB ITQ leads in both accuracy and throughput at every footprint. Throughput shows the opposite trend to accuracy, and the gap widens with the code width. Across the SIFT1M matrix the throughput of R AB ITQ spans 1.4×, while that of PQ spans 9.8×. At the widest code R AB ITQ is 3.6, 8.9 and 10.8× faster on the three collections. PQ is the faster encoder only on SIFT1M below ≈40 bytes/vector. At matched footprint, 𝑚 scales with 𝑑. The per-document work of PQ is 𝑚 dependent lookups into a table of 𝑚 ⋅ 256 floats. This table takes 16 KiB and fits in L1 at SIFT1M’s 𝑚=16, but takes 1 MiB at SNOW FLAKE’s 𝑚=1 024. The per-document work of R AB ITQ is a fixed number of vectorized multiply–adds over contiguous bytes. Build times differ by orders of magnitude. R AB ITQ has no training phase, while PQ must fit 𝑚 codebooks before it can encode anything: training takes 493 of the 634 seconds on SNOW FLAKE at 𝑚=1 024. The ratio runs from 62× to 1 062× and generally widens with 𝑚. Fewer k-means iterations or a smaller training sample reduce PQ’s build time by at most one order of magnitude before accuracy degrades appreciably, while the gap spans two to three orders. Page 12 of 39

Modular and Efficient Approximate Nearest Neighbors Search with Vectorium and kANNolo Table 5 Validation against the authors’ reference implementation of RaBitQ. Exhaustive top-10 scan, all queries, conventions of Section 3.5 on both sides. Best of three repetitions, the reference’s accuracy the mean of three since its rotation is seeded non-deterministically. Grey figures are our factor of advantage per column. Accuracy@10 (%) 1 bit

2 bit

QPS

4 bit

8 bit

1 bit

Sift1M (𝑑=128, 𝓁2 ) ours 25.7 55.5 84.5 ref. 24.9 55.5 83.9

97.5 97.0

1 988 1.64× 1 664 1.44× 1 608 1.74× 1 443 1.74× 1 209 1 157 922 831

Dragon (𝑑=768, IP) ours 66.6 82.3 93.9 ref. 66.9 82.5 94.1

99.3 99.3

151 1.44× 105

2 bit

Bytes/vec

4 bit

60 1.85× 32

8 bit

52 2.37× 22

30 1.86× 16

1 bit

2 bit

4 bit

8 bit

24 32

40 56

72 88

136 152

104 112

200 216

392 408

776 792

Table 6 RaBitQ against Product Quantization at matched footprint. Exhaustive FlatIndex top-10 search, all queries, Accuracy@10 (%) against the exact ground truth, footprint from SpaceUsage, conventions of Section 3.5. Build is the complete construction cost. PQ’s training is not deterministic: its accuracy is the mean of three trainings, which differ by at most 0.40 points. QPS is from one repetition, so it can be lower than the best-of-three figures of Table 5. Bytes/vec

Accuracy@10 (%)

QPS

Build (s)

RaBitQ

PQ

RaBitQ

PQ

RaBitQ

PQ

RaBitQ

PQ

1 bit / 𝑚=16 2 bit / 𝑚=32 4 bit / 𝑚=64 8 bit / 𝑚=128

24 40 72 136

16 32 64 128

25.7 55.5 84.5 97.5

55.2 71.5 87.8 99.6

1 988.1 1 661.6 1 606.1 1 442.3

3 966.3 2 060.2 889.1 402.8

0.031 0.063 0.058 0.060

14.3 25.2 40.0 63.7

Dragon 𝑑=768, IP

1 bit / 𝑚=96 2 bit / 𝑚=192 4 bit / 𝑚=384 8 bit / 𝑚=768

104 200 392 776

96 192 384 768

66.6 82.3 93.9 99.3

68.5 83.2 94.2 98.4

148.3 56.8 46.6 30.2

60.8 26.3 12.1 3.4

1.262 2.676 2.601 2.615

98.1 175.6 289.5 475.0

SnowFlake 𝑑=1 024, IP

1 bit / 𝑚=128 2 bit / 𝑚=256 4 bit / 𝑚=512 8 bit / 𝑚=1 024

136 264 520 1 032

128 256 512 1 024

82.1 90.4 96.6 99.7

79.8 88.0 95.4 98.6

143.1 44.6 29.3 17.3

46.4 19.7 9.0 1.6

1.842 3.756 3.637 3.682

132.3 233.8 386.2 634.0

Dataset

Width

Sift1M 𝑑=128, 𝓁2

4.4.5. Blocked Sparse Encoding BlockedSparseEncoder stores the same u16 components and f16 values as the scalar quantizer, in a layout designed for scoring. Each vector is split into blocks of 16 (component, value) pairs: 32 bytes of components followed by 32 bytes of values, which fill one 64-byte cache line. A query with at least 33 nonzeros is densified, as with the scalar quantizer. A shorter query is kept as a sorted list and scored with JUMPD OT [13]. The kernel walks the blocks of the candidate in order and skips every block whose last component precedes the current query component. It compares the query component against the 16 components of a block at once. Unlike the kernel family, this kernel uses explicit SIMD (AVX2, and NEON on Armv8). 4.4.6. Compressed Sparse Encoding: D OTPACKING 8 The scalar quantizer of Section 4.4.2 compresses only the values of a sparse vector. Each nonzero still carries a raw u16 component index, and on learned-sparse embeddings [28, 29, 99] this index is the largest remaining cost once the values are compressed. DotPacking8Encoder is the S. Martinico et al.: Preprint submitted to Elsevier

only encoder in the library that compresses it. The component stream is a sorted list of integer identifiers, stored as the gaps between consecutive ones. It is compressed with D OTPACKING8 [13], the bit-packing code for inner products of Section 2.2. Blocks of eight gaps are packed at the width of their largest gap. A byte shuffle decodes a block into a SIMD register, and a prefix sum over the gaps gives the component identifiers, which gather the query weights. A compressed vector is thus scored without being materialized. V ECTORIUM integrates the authors’ implementation. It stores the values in 8-bit fixed point, f16, or per-component 4-, 8- or 16-bit scalar codes, so both halves of the vector can be compressed. The gaps between consecutive component identifiers of a vector, and hence its compressed size, depend on how the coordinate axes are numbered. Any permutation 𝜋 of the axes may be applied to the stored components, provided the same 𝜋 is applied to every query. The encoder computes such a 𝜋 at training time, numbering components by decreasing document frequency, and stores it with its inverse. It relabels

Page 13 of 39

Modular and Efficient Approximate Nearest Neighbors Search with Vectorium and kANNolo

components through 𝜋 before packing and undoes the relabeling on decode, so the permutation is invisible to the rest of the system.

4.4.7. Multi-Vector (Late-Interaction) Encoders The encoders above compress one vector per document. Late-interaction retrieval [32, 33] (Section 2.1) represents a document as a set of token vectors and scores it by MaxSim. V ECTORIUM added this representation last, after the abstractions of the preceding subsections were fixed. It therefore tests directly whether these abstractions accept a new representation at the cost the design claims. Supporting it required four additions, one per layer, all behind an optional compiletime feature: the multi-vector view of Section 4.2, the maxsim kernel of Section 4.3, a multi-vector encoder contract, and the multi-vector dataset of Section 4.5. The Dataset contract itself, its access and prefetch primitives, persistence, and the Index and RerankIndex contracts of Section 5 were reused unmodified. A multi-vector collection is a Dataset whose encoded vector is a token buffer, and a MaxSim score is treated like a DotProduct.

Encoders. The same compression ladder is reproduced

for this representation. ScalarMultiVecQuantizer applies the cast of Section 4.4.2 to every token value. It iterates over document tokens in the outer loop, so a stored f16 token is widened once per document token and not once per query token. MultiVecProductQuantizer encodes each token with the Product Quantizer of Section 4.4.3, at 𝑚 bytes per token. For 𝑄 query tokens it builds a single table of shape [𝑚 × KSUB × 𝑄] with the query-token index innermost. A code byte thus selects the contribution of that subspace to every query token at once, and a document token is scored by 𝑚 vector accumulations in place of 𝑚 ⋅ 𝑄 scalar lookups. This applies the reasoning of the SoA centroid layout one level up. MultiVecTwoLevelProductQuantizer assigns each token to one of a set of k-means centroids and applies Product Quantization only to the residual. A token therefore costs 4 bytes of centroid identifier plus 𝑚 bytes of code. The coarse contribution, every document token’s centroid against every query token, is a small dense matrix product delegated to a GEMM. The residual contribution uses the accumulation above. Of all the encoders, this one tests the storage/scoring separation of Section 4.1 most. Its stored token is a byte payload whose “dimension”, 4 + 𝑚, only the encoder can read. The payload is laid out in blocks, all coarse identifiers and then all codes, so each scoring phase walks its own contiguous region. An optional variant normalizes each residual and appends the discarded norm as an f16 third block. This raises the width to 4 + 𝑚 + 2 and requires no change to the storage layer. The MaxSim kernel is monomorphized over the query-token count 𝑄 for the powers of two from 8 to 512. This turns its two 𝑄-element buffers into stack arrays and its inner loops into loops with a fixed trip count. A runtime-𝑄 variant covers other lengths; typical late-interaction queries, such as the 32 tokens of COLBERT V2, do not need it.

S. Martinico et al.: Preprint submitted to Elsevier

K ANN OLO ’s and V ECTORIUM ’s indexes could be used to search multi-vector representations. However, MaxSim is expensive enough to compute that neither is well suited to this representation. We therefore leave its indexing to engines specialized for it [101, 102], such as TACHIOM [12]. This engine builds its late-interaction index on the multi-vector views and dataset, the two-level encoder and its MaxSim kernel, and routes queries to its token centroids through KANNOLO’s H NSW . Late interaction is also used as the reranking step of a two-stage pipeline. This use is supported through a multi-vector Dataset inside R ERANKINDEX (Section 5) and is evaluated across all three encoders in [1]. There, a learned sparse first stage, searched with KANNOLO’s H NSW or with SEISMIC, is followed by COLBERT V2 MaxSim rescoring. This pipeline matches or exceeds the effectiveness of the token-level systems E MVB [102], WARP [103] and IGP [104], and is up to 24× faster than the fastest of them at comparable memory.

4.5. Datasets A Dataset holds a whole collection of encoded vectors. It is the pure storage of Section 4.1, parameterized by a VectorEncoder. It reports its size, returns a vector by identifier or by storage range, iterates over its vectors, and issues prefetch hints. It knows nothing else about the encoding. Four dataset types cover the four layouts of Section 4.2. Dense data is held in one contiguous buffer. Sparse data uses the compressed-sparse-row form: an offsets array plus parallel component and value arrays. Packed encodings are concatenated into one blob of words with the same offsets table; this is the storage form of D OTPACKING8 and of the blocked encoder. Multi-vector documents are stored as an offsets table over a flat token buffer. Datasets convert into one another in two ways. A zerocopy conversion keeps the stored bytes and only reinterprets the encoder. An example is switching a dense or sparse dataset between squared-Euclidean distance and inner product, which share a layout. A re-encoding conversion walks the data once to change its stored form, for example casting f32 to f16 or encoding a plain dataset into ProductQuantization codes. This conversion is also where an encoder receives the parameters it cannot learn from the data (Section 4.4.4). Every dataset persists and reloads without a rebuild, in a fixed-width little-endian format that is stable across hosts.

Identifier/range duality. Each stored vector can be

fetched by its identifier or by its span in the backing storage. For the sparse layout the difference matters. Fetching by identifier first reads the vector’s boundaries from the offsets table. This is a dependent load, which must complete before the vector itself can be requested. An index that walks a candidate list can resolve each range once, prefetch the vector body, and fetch it by range only when it scores it. The memory latency then overlaps with other work. For dense data the two methods coincide. The prefetch hint touches at most the first 256 bytes of the span. This covers a short or Page 14 of 39

Modular and Efficient Approximate Nearest Neighbors Search with Vectorium and kANNolo

compressed vector and leaves the tail of a long one to the hardware prefetcher. Under random access on COCONDENSER, range addressing is 1.46× faster than identifier addressing, and 2.07× with a prefetch issued ahead of use; in storage order the three strategies tie (Appendix C.3).

4.6. Clustering The layers above also depend on a primitive that is neither a representation, a metric nor an encoder: k-means. It learns the per-subspace codebooks of the Product Quantizer (Section 4.4.3), the coarse level of the two-level multi-vector quantizer (Section 4.4.7), and the partition from which an IVF index is built (Section 5). It comes in a plain and a spherical variant. The spherical variant 𝓁2 -normalizes the centroids after each update, so an inner-product collection is clustered by angle. KMeans is generic over the representation. It clusters any dataset that can average a set of its vectors and subtract a centroid from one of them; the dense and sparse datasets currently support both operations. It returns the centroids as a Dataset whose encoder the caller chooses. The metric of that encoder is the clustering metric, which is independent of the metric used later to search the collection. K-means therefore has no distance loop of its own. Assignment is a search over the centroid dataset with V ECTORIUM ’s distance kernels, either an exhaustive scan or a top-1 query against any Index built over the centroids. The training corpus and the centroids may use different storage types. A collection already stored as f16 is therefore clustered directly, with no widening pass, while means accumulate in f32. The update step groups assignments with a stable counting sort and sums each cluster in parallel, which fixes the summation order. With a fixed seed and exhaustive assignment, the centroids are bit-identical on a given target at any thread count. On sparse data the mean of a cluster is much denser than its members. After every update, each centroid can therefore be pruned to its largest-magnitude components, keeping a given factor times the average number of nonzeros in its cluster.

5. K ANNOLO K ANN OLO is the indexing layer above V ECTORIUM . It implements V ECTORIUM ’s index contract (Section 5.1) with two standard index families: a graph index (H NSW [64]) and an Inverted Vector File (IVF) index [4]. It also provides the search kernels, construction algorithms and adjacency layouts these indexes are built from, and the twostage R ERANKINDEX that composes them. Both families are generic over the representations, encoders, metrics and datasets, and name no encoding of their own, so one implementation of each serves dense and sparse data. Multivector data enters KANNOLO as the rescoring dataset of R ERANKINDEX (Section 4.4.7).

S. Martinico et al.: Preprint submitted to Elsevier

5.1. System Overview V ECTORIUM defines the Index contract, and KANNOLO is one library that implements it. Any other library built on V ECTORIUM can implement it as well. The contract covers search only: given a query and per-family search parameters, it returns the top 𝑘 results as (distance, id) pairs. Construction differs across families, so each family exposes its own. IVF, for example, also chooses the encoding its clusters are stored with (Section 5.5). V ECTORIUM also ships the default implementation, FlatIndex, which makes a dataset searchable by exhaustive scan in Θ(𝑛 log 𝑘) time. Because this index lives in the data layer, V ECTORIUM can work on its own with no dependency on the indexing layer, for example to assign points to centroids during k-means (Section 4.6).

RerankIndex. RerankIndex is the two-stage composition of

Section 2.3.1 written as a single type. It pairs a first-stage index with a companion reranking dataset. The two may use different vector encodings and distance metrics. For example, stage one can be an H NSW index over sparse representations that retrieves candidates cheaply, while stage two holds multi-vector embeddings and rescores the candidates with that representation, as in [1]. Alternatively, stage one can search a PQ-compressed copy of the same collection that stage two holds at full precision. The only requirement is that the two representations are aligned on the same vector identifiers, and the caller is responsible for it. The pipeline proceeds in two stages. First, the first-stage index retrieves 𝑘c candidates using its own search parameters. Then the reranking dataset scores every candidate against the query using the full representation, and the top 𝑘final are returned. Two optional optimizations reduce the cost of the second stage [1]: • Candidate Pruning (CP, 𝛼). The first-stage scores are read here as similarities, so a higher value is better. Before reranking, candidates whose first-stage score falls below (1 − 𝛼) times the score of the 𝑘final th candidate are discarded. This exploits the correlation between first-stage and second-stage scores to prune clearly suboptimal entries before the expensive reranking step. • Early-Exit Reranking (EE, 𝛽). The first-stage results are already sorted by first-stage score. During reranking, if 𝛽 consecutive candidates fail to enter the current top-𝑘final heap, the reranking loop terminates. The assumption is that the remaining candidates are unlikely to improve the ranking. A third option covers the case where the two stages score complementary parts of one representation, i.e. where the vector is (or is approximated by) the sum of the corresponding first- and second-stage vectors. The reranking dataset then holds a residual, and the final score is the sum of the two scores. This decomposition holds only under an additive metric, and the caller must ensure that it is meaningful. Page 15 of 39

Modular and Efficient Approximate Nearest Neighbors Search with Vectorium and kANNolo type Rescoring = DenseDataset<PlainDenseQuantizer<f16, DotProduct>>; type Codes = DenseDataset<RabitqQuantizer<DotProduct>>; let exact: Rescoring = DenseDataset::from_raw( values, n, PlainDenseQuantizer::new(dim)); let cfg = HNSWBuildConfiguration::default() .with_num_neighbors(16) .with_ef_construction(200); let graph: HNSW<Rescoring, Graph> = HNSW::build_index(exact.clone(), &cfg); let first: HNSW<Codes, Graph> = graph .permute_and_encode::<PlainNeighbors>() .convert_dataset_into_ref(RabitqConfig::default()); let index = RerankIndex::new(first, exact); let top10 = Index::search(&index, (q, q), 10, &RerankSearchParams::new(100, params, ()));

Figure 2: Building and searching the dense two-stage index of Section 6. Hnsw is built over f16 vectors and permuted. Its vectors are replaced by 1-bit RaBitQ codes, and candidates are rescored on the f16 copy, with 𝑘c = 100 and 𝑘final = 10.

Dataset, Index and RerankIndex compose. Any index implementing the contract, flat or one of KANNOLO’s, can serve as a RerankIndex’s first stage, and any V ECTORIUM dataset as its second. This yields many high-performance indexes, as Section 6 shows. Table 7 lists the meaningful compositions for each representation, with the published method each one corresponds to. The types allow many more. A first stage over f32 rescored with PQ gains nothing in a deployed index, but it can still be built. Such combinations can serve experiments and analysis, for example to isolate the effect of one stage. Figure 2 follows one composition from encoder to search. It is the dense two-stage index of Section 6: a 1-bit R AB ITQ first stage over a permuted H NSW graph, rescored on an f16 copy of the collection. The graph is built over the f16 vectors, and its dataset is then replaced by the codes, so the codes never take part in construction. The result is itself an Index, searched through the same call as its first stage. Each stage is named by a type alias. Replacing Codes with a PQ or R AB ITQ-ext dataset, and its configuration accordingly, yields another dense row of Table 7 with no other change.

5.2. Graph Index Framework K ANN OLO abstracts the storage layer of any proximity

graph [65, 67, 62, 61] behind a GraphTrait of six primitives: a node’s neighbor list, the node count, the edge count, the degree bound, memory accounting, and a translation from a local graph node ID to the corresponding identifier in the underlying dataset. The identifier translation decouples the graph’s internal node ordering from the external ID space of the application. The graph can therefore be permuted without invalidating any external reference, e.g. to improve cache locality or edge compressibility, or to map the local IDs of H NSW ’s upper levels onto the ground level. When no mapping is S. Martinico et al.: Preprint submitted to Elsevier

configured, the translation is a no-op. A mapping holds one entry per node, bit-packed at ⌈log2 𝑁⌉ bits. The two search kernels, a greedy walk and a beam search, are implemented once on GraphTrait. The greedy walk is the single-best-neighbor descent used for 1-NN search and for traversing the upper levels of H NSW . The beam search is used by every graph index and by the ground level of H NSW . It maintains a max-heap of the current top-𝑘 candidates and a min-heap of the traversal frontier, and processes nodes in best-first order until the termination criterion fires. A graph is built in a growable fixed-degree layout, where each node owns a fixed number of contiguous slots. While the graph is still growing, a neighbor list is located, read and extended with no offsets indirection. This layout needs no locking (Section 5.4.1). Once built, the graph is frozen into one of three search-time backends. The plain backend is a CSR-like layout: neighbor IDs plus an offsets array delimiting each node’s list. The fixed-degree backend pads every list to a common stride and drops the offsets array, which gives direct random access to any node’s list. The third is the compressed layout of Section 5.3.4. Its delta coding needs each list sorted by ID, and the conversion sorts the lists. The plain layout is the default for search, and the rest of this section writes the graph as Graph.

5.3. Graph Optimizations The beam-search kernel is implemented on GraphTrait, so the optimizations below apply to every graph backend and to any graph index added later. We measure all of them on the H NSW index. They are presented in the order in which they act during a query.

5.3.1. Visited Set Design and Optimization A graph index does not partition the search space, so the search has to track which nodes it has already evaluated. It keeps them in a “visited set”, which prevents redundant distance computations and infinite loops. The set is queried and updated in the innermost loop of the search, so its cost is on the critical path of every expansion. We evaluate two structurally different alternatives: • BitVec, a dense bitmap over the whole graph. Access is 𝑂(1) by direct indexing, but each query pays an 𝑂(𝑛) initialization. • FxHashSet, a fast non-cryptographic hash table sized from the expected number of visited nodes. Access is 𝑂(1) amortized and initialization is bounded by the table’s capacity, independent of 𝑛. A third common alternative is a global generation array. It stores, for each node, the generation of the last query that visited it. A new query increments the generation and does not clear the array, so there is no per-query initialization. We exclude it because it requires 𝑂(𝑛) permanent memory, 200 MB at 𝑛 = 100M with 16-bit generations, and introduces synchronization between threads querying the same index. Page 16 of 39

Modular and Efficient Approximate Nearest Neighbors Search with Vectorium and kANNolo Table 7 Index compositions expressible with Vectorium and kANNolo. Known as names the published method with the same decomposition, where one exists; Evidence points to where the composition is measured. Scalar stands for f32, f16, bf16 and 16/8-bit fixed point. RaBitQ stands for both the 1-bit RaBitQ and RaBitQ-ext. IVF clusters store either the vectors or, on dense data, their residuals from the centroid. A rescoring encoder is more precise than the first-stage one. † The Hnsw and IVF indexes of the RaBitQ authors’ library. Both prune candidates with a cheap estimate and its error bound before computing the full one; kANNolo does not (Sections 4.4.4 and 5.5.2). First-stage index

First-stage encoder

Rescoring

Known as

Evidence

Dense Flat Flat Flat Flat Hnsw Hnsw Hnsw Hnsw Hnsw IVF, flat centroids IVF, Hnsw centroids IVF, flat or Hnsw centroids IVF, IVF centroids IVF, any centroids

scalar PQ RaBitQ PQ, RaBitQ scalar PQ RaBitQ PQ RaBitQ scalar, PQ scalar, PQ RaBitQ scalar, PQ scalar, PQ, RaBitQ

— — — scalar, PQ, RaBitQ — — — scalar, PQ, RaBitQ scalar, PQ, RaBitQ — — — — scalar, PQ, RaBitQ

— PQ [37] — — Hnsw [64] HNSWPQ — HNSWPQ + refine Hnsw-RaBitQ † IVFFlat, IVFPQ [37] IVF-Hnsw IVF-RaBitQ † [35, 36] — —

Table 3 Tables 4, 6 Tables 5, 6 — Sections 5.4.2, 6 Figure C.2 — Section 6.3.3 Sections 6, 6.3.3 Figure 5 Figure 5 Table 9 — —

scalar, blocked, DotPacking8 scalar, blocked, DotPacking8 scalar, DotPacking8 scalar, DotPacking8 top-𝑎% of 𝓁1 mass, scalar or blocked

—

—

Table 3

—

—

Section 6.3.2

— scalar scalar or blocked residuals

— — —

Section 5.5 — Section 6.3.2

Sparse Flat Hnsw IVF, any centroids Flat, Hnsw, IVF Flat, Hnsw, IVF Multi-vector Flat, Hnsw Flat, Hnsw

scalar, PQ, two-level PQ scalar, PQ, two-level PQ

— scalar, PQ, two-level PQ

— —

— —

Hybrid Flat, Hnsw, IVF Flat, Hnsw, IVF Flat, Hnsw, IVF

any sparse any sparse any dense

any dense any multi-vector any sparse, any multi-vector

— — —

— [1] —

Cascade RerankIndex

any RerankIndex

any

—

—

Sizing and choice. Sizing the hash table requires an es-

timate of the number of nodes a query visits, 𝑛visited . Rehashing is expensive, while an oversized table wastes memory and cache. Only the ground level of H NSW keeps a visited set, since the greedy walk on the upper levels cannot revisit a node. On the ground level, a query expands about 𝑒𝑓search nodes and scores their unvisited neighbors. Distance-adaptive termination (Section 5.3.3) multiplies the expansions by a factor that grows faster than linearly in 𝜆. The table is therefore sized as ( )8 capacity = 𝑑̄ ⋅ 𝑒𝑓search ⋅ 1 + 1𝜆>0 (𝜆 + 0.05) + 200, where 𝑑̄ = ⌊|𝐸|∕𝑛⌋ is the average ground-level degree and 1𝜆>0 is 1 when early termination is on. Appendix C.4 derives this form and checks it per query. Without early termination, at most 7.2% of the queries rehash, and almost none rehash twice. At 𝑛 = 100M, clearing a BitVec costs 407–455 𝜇s per query, more than the whole query budget. BitVec is 2.2×

S. Martinico et al.: Preprint submitted to Elsevier

to 4.4× slower than FxHashSet in every configuration at this size. At 𝑛 ≤ 1M it is ahead by at most 32 𝜇s (Appendix C.4). K ANN OLO therefore uses FxHashSet unconditionally.

5.3.2. Neighbor Processing The innermost loop of the graph beam search reads a candidate’s neighbor list, fetches the bodies of the neighbors not yet visited, and scores them against the query. It executes hundreds to thousands of times per query, so its efficiency governs the throughput of the whole index. Two optimizations act on it, each on a different bottleneck: • Batched distance computation. The expansion buffers up to six unvisited neighbors and scores them in a single pass over six independent accumulators, through the batched kernel of Section 4.3. • Prefetching. Before scoring, the loop issues explicit PREFETCHT2 hints, targeting the L3 cache, for the batch it is about to score. The prefetch window is capped Page 17 of 39

Modular and Efficient Approximate Nearest Neighbors Search with Vectorium and kANNolo

at the first 256 bytes (4 cache lines) of each vector. Any remaining cache lines are left to the hardware prefetcher. Both are built from the Dataset access primitives of Section 4.5 and neither adds a primitive of its own. We measure them on four collections whose vectors span 2 to 48 cache lines (Appendix C.5). Enabling both is the fastest or tied-fastest configuration on all four. Prefetching alone gives up to +35% QPS where a whole vector fits in the prefetch window. Batching gives +34% on the largest vectors, where it cuts retired instructions by about a quarter at constant IPC. Hardware counters show that the 4-cacheline cap is the point where two costs balance: the first-access stall that the hints remove, and the issue slots that the hints occupy.

5.3.3. Early Termination Strategies K ANN OLO implements a relaxation-based termination criterion (Section 2.3.1) that generalizes Adaptive Beam Search [71]. Adaptive Beam Search keeps only the 𝑘 best candidates found so far. It stops when the best unexpanded candidate is farther from the query than 1 + 𝛾 times the 𝑘-th of them. KANNOLO keeps the beam of 𝑒𝑓search candidates and relaxes the comparison against its worst element. A candidate is expanded while its distance is within a factor 1+ 𝜆 of the 𝑒𝑓search -th best, through the relaxation predicate of Section 4.3. The authors of Adaptive Beam Search leave this hybrid of beam width and distance threshold open. It contains both rules: 𝜆 = 0 is the standard H NSW criterion, and 𝑒𝑓search = 𝑘 is Adaptive Beam Search, with 1 + 𝜆 = (1 + 𝛾)2 on squared distances. Tuning both parameters can therefore only match or improve on either rule. The approximation guarantee of Adaptive Beam Search holds for 𝑒𝑓search = 𝑘 and Euclidean distance on navigable graphs. H NSW graphs are not provably navigable, so neither form carries the guarantee here. We compare this rule against the saturation-based family, represented by Patience-in-Proximity (PiP) [72]. For this comparison we implemented PiP in the search loop of K ANN OLO ; it is not part of the released library. The search stops after Δ consecutive stable expansions, each of which admits at most ⌊(1 − 𝑝) 𝑘⌋ new candidates into the result set. We evaluate the two strategies on SIFT1M and COCONDENSER, which differ in representation and in metric (Table 1). We measure QPS and distance computations per query at Accuracy@10 targets from 90% to 99%. For each target we report the fastest configuration that reaches it and, separately, the one that reaches it with the fewest distance computations. Standard H NSW sweeps 𝑒𝑓search from 10 to 400. The relaxation grid crosses 𝑒𝑓search from 𝑘 = 10 to 20 with 𝜆 from 0.005 to 0.3, so Adaptive Beam Search is among its configurations. The PiP grid crosses 𝑒𝑓search from 10 to 400 with Δ from 3 to 100 and 𝑝 ∈ {0.95, 0.85, 0.75}. All strategies run on the same indexes. Figure 3 reports QPS and the reduction in distance computations against standard H NSW termination. Distanceadaptive termination is ahead on both representations and at every accuracy target. On SIFT1M it reduces distance S. Martinico et al.: Preprint submitted to Elsevier

computations by 5% at 90% Accuracy@10 and by 18% at 99%. On the sparse COCONDENSER the reduction goes from 5% at 90% to 34% at 99%, where QPS rises by 50%. The advantage grows with the accuracy target, and retrieval systems operate at high targets. The winning configurations keep 𝑒𝑓search between 10 and 17, close to 𝑘, and the beam adds little. Restricting the grid to 𝑒𝑓search = 𝑘, that is, to Adaptive Beam Search, costs at most 5.4% of the throughput on SIFT1M and 7.4% on COCONDENSER at any target. The single-parameter rule is therefore a sound default when one parameter is preferred. PiP is close to neutral. It saves at most 7.5% of the distance computations on SIFT1M and 4.5% on COCONDENSER, and under 5% at most targets. Its QPS stays within 8% of standard H NSW everywhere. At one target PiP computes more distances, 1.6% more at 93% on COCONDENSER. This reflects the grid resolution. An unbounded patience reduces PiP to the standard criterion at the same 𝑒𝑓search , and the best configurations use a patience between 5 and 30, well inside the grid. The difference follows from what each criterion can change. A fixed beam stops once 𝑒𝑓search candidates closer than the frontier have been found, however many near-ties remain around the query. A distance threshold adapts to the query [71]. It spends distance computations on hard queries with many near-ties and stops early on easy ones. A small beam with a relaxation thus reaches a given accuracy with fewer computations than the larger beam that the standard criterion needs. PiP can only stop a search earlier than the standard criterion would at the same 𝑒𝑓search ; it can never extend it. The result set stops changing at about the time the standard criterion fires. PiP therefore saves only the tail of each search and cannot lower the 𝑒𝑓search that a target requires. At high accuracy, the last true neighbors are found in that tail. The patience has to grow to keep them, and PiP converges to the standard criterion.

5.3.4. Graph Compression The adjacency lists are the other half of the memory cost of a graph index (Section 2.4): with the plain backend, every edge costs 4 bytes and every node an 8-byte offset. The compressed backend implements the scheme of [90], which was originally developed in KANNOLO. After the build, the EGB permutation [88] of the ground level relabels every level and the dataset, so graph-adjacent nodes receive close identifiers. Each sorted list is then stored as small gaps with STREAM VBYTE [89], and the per-node offsets with Elias–Fano [105, 106]. The permutation can also be applied to an uncompressed graph. It improves the locality of the neighbor accesses, so a compressed index can be as fast as an uncompressed, unpermuted one, or faster. Only the final top-𝑘 results are mapped back to the original identifiers. A STREAM VBYTE block holds at most 256 values, so 𝑀 is capped at 128.

5.4. H NSW K ANN OLO ’s H NSW index [64] is generic over the dataset and the graph storage. Every level is built in the

Page 18 of 39

Modular and Efficient Approximate Nearest Neighbors Search with Vectorium and kANNolo CoCondenser

15000

3000

10000

2000

5000

1000

0

0

90

95

Accuracy@10 (%)

100

Sift1M Distances saved (%)

QPS

Sift1M

CoCondenser 30

15

20

10

10

5

0

0 90

95

100

90

Accuracy@10 (%) Standard Hnsw

95

100

Accuracy@10 (%)

Distance-Adaptive

90

95

100

Accuracy@10 (%)

Patience-in-Proximity

Figure 3: Early termination strategies on Sift1M and CoCondenser, against standard Hnsw termination and Patience-inProximity.

mutable fixed-degree representation and frozen into the target graph type once complete. The layered structure is stored as an array of graphs, from the topmost level down to the ground level with all 𝑁 nodes. A query starts from the single node of the top level, descends the upper levels with the greedy walk, and runs the beam search on the ground level. The upper levels work in local IDs, since the nodes of each are a prefix of the level below. A mapping array translates the last of them to the ground level, whose order may be permuted (Section 5.3.4). K ANN OLO carries H NSW ’s three standard parameters of Section 2.3.1, with defaults 𝑀 = 16, 𝑒𝑓construction = 150 and 𝑒𝑓search = 100. The value of 𝑒𝑓search is clamped to max(𝑒𝑓search , 𝑘) so that at least 𝑘 candidates are evaluated. A second search parameter selects the termination criterion: the standard H NSW rule, which is the default, or the distanceadaptive rule of Section 5.3.3, with its own 𝜆.

5.4.1. Index Construction Standard implementations parallelize H NSW construction by inserting nodes in parallel, assigning the next node to a free thread. They use locks on neighbor lists to make concurrent access and updates safe. The locks are a bottleneck on large collections, and they also make the resulting graph depend on thread timing. KANNOLO follows the lock-free batched scheme of PARLAYANN [6], applied to each level of the hierarchy in turn. It departs from PARLAYANN in two places: the rule that sets the batch size, discussed at the end of this section, and the way the merge is made race-free. Each batch proceeds in two phases: 1. Phase 1 (parallel search). Every node in the batch searches the current graph for its candidate neighbors. It then prunes them to the degree bound with the pruning heuristic of the original algorithm. The graph is not modified during this phase, so it needs no locks and scales across all cores. 2. Phase 2 (merge and commit). The batch’s edges are written. Forward edges go to slots that belong to the new nodes. Reverse edges, i.e. the links from existing nodes to the new ones, are first merged per target, as described below. S. Martinico et al.: Preprint submitted to Elsevier

Reverse-edge insertion. Within a batch, several new

nodes often pick the same existing node 𝑣 as a neighbor, so several reverse edges arrive at 𝑣 together. The capacity of 𝑣’s neighbor list is fixed, so the edges cannot simply be appended. They have to be pruned against the current content of the list, and against each other. KANNOLO therefore prunes per target: all the edges arriving at a node are considered together. The batch’s reverse edges are collected as (target, new node) pairs and sorted. Each target is then pruned once, with all of its incoming edges in the candidate set. The sort has two effects. First, it makes the merges independent, so they run in parallel over distinct targets. Second, it puts every arriving edge into the single pruning decision that concerns it, so the only edges dropped are those the heuristic rejects. Distances are computed only where the degree bound forces a choice. With this procedure, adding reverse edges cannot cause a data race. Every write is deferred to a commit that runs afterwards over disjoint slots. The reverse targets of a batch are always already-inserted points, while its forward edges go to slots owned by the new points. The Rust borrow checker therefore verifies the absence of data races, and the build-time adjacency needs neither atomics nor interior mutability.

Reproducible construction. Determinism is a design goal

of KANNOLO, because prototyping and testing depend on it. Two things must not depend on thread timing. The first is the set of nodes each insertion sees, which the batch boundary already fixes. The second is the order of a target’s candidate set during the merge. Sorting the reverse edges by target and then by arriving node makes this order canonical. With levels and insertion order drawn from a fixed seed, two builds of the same dataset with the same parameters produce byteidentical indexes on the same target.

Adaptive batch sizing. The batch size balances two costs.

Small batches leave threads idle, because the time a vector takes to find its neighbors varies, much more on sparse collections. Large batches cost accuracy, because nodes in the same batch search a frozen graph and cannot discover each other. The damage depends on the ratio ∕𝑁 between Page 19 of 39

Modular and Efficient Approximate Nearest Neighbors Search with Vectorium and kANNolo

The batch targets a fixed fraction 𝜌 of the graph, which bounds the part a batch cannot see. It never falls below a floor 𝐵floor , stated in tasks per worker thread, so the first rounds already use every core. It never exceeds 𝜌max 𝑁, with 𝜌max = 0.05, so the floor cannot make the batch too large while the graph is small. We set 𝜌 on SIFT1M, whose run-to-run variation is narrow enough to resolve differences of this size. At the worst point of a sweep over 𝑒𝑓search , Accuracy@10 is unaffected at 𝜌 = 0.005, and falls by at most 0.09 points at 𝜌 = 0.01 and 0.15 at 𝜌 = 0.02. Build time falls as 𝜌 grows, with diminishing returns. On the 8.8M collections the differences stay within run-to-run variation. We therefore take 𝜌 = 0.005, four times smaller than PARLAYANN’s 𝜃.

5.4.2. Validation against hnswlib We validate our H NSW implementation against hnswlib [2], the reference implementation maintained by the algorithm’s authors. We use the simplest configuration that both systems support: dense f32, uncompressed edges, no early termination, no reranking. Both sides use 𝑀 = 32 and 𝑒𝑓construction = 200. hnswlib’s inner-product space returns 1 − ⟨𝑞, 𝑣⟩, which is rank-equivalent to the dot product K ANN OLO scores with. H NSW ’s guarantees are average-case over a random insertion order. Both the expected logarithmic search and the navigability argument assume that the inserted points arrive in random order. hnswlib randomizes each point’s level but inserts in the order the caller supplies, so it satisfies only half of that precondition. This makes no difference on a collection stored in random order, but it matters on one that is not, and text collections typically are not. On D RAGON, which is stored by document, consecutive rows have mean cosine 0.478 against 0.184 for random pairs. On SIFT1M, consecutive rows are no closer than random ones. We therefore report hnswlib twice: in file order as shipped, and under a random permutation of the insertion sequence. KANNOLO shuffles before assigning levels and is unaffected. On D RAGON, permuting the insertion order alone adds +1.98 points of Accuracy@10 on average over the 𝑒𝑓search ladder, from +3.71 at 𝑒𝑓search = 25 to +0.83 at 520 (Figure 4, dashed). On SIFT1M, which has no order to exploit, it adds +0.12. Measured in distance computations, the file order needs ≈1.4× as many as the permuted order for equal accuracy. Against hnswlib inserted in random order, the two implementations are indistinguishable on graph quality. Accuracy@10 differs by +0.051 points on average on SIFT1M and +0.074 on D RAGON. Both differences are inside hnswlib’s own spread over three independent builds. At matched 𝑒𝑓search , the two perform the same number S. Martinico et al.: Preprint submitted to Elsevier

Table 8 Validation of kANNolo’s Hnsw against hnswlib [2]; conventions of Section 3.5. original inserts in file order, fixed in permuted order. Footprint is peak resident memory. Collection

System

Sift1M (𝑑=128, 𝓁2 ) Dragon (𝑑=768, IP)

Build (s)

Bytes/vec

hnswlib (original) hnswlib (fixed)

kANNolo

25.6 25.6 28.0

653.8 895.1 894.9

kANNolo hnswlib (original) hnswlib (fixed)

1 710 1 498 1 726

3 214.6 3 444.7 3 444.8

Sift1M

Dragon 3000

10000

QPS

the batch size and the number of nodes already inserted. PARLAYANN doubles the batch until it reaches a single cap 𝜃, set empirically to 0.02𝑛. We keep the doubling and replace the cap with three rules: ( ) ( ) (𝑁) = min max 𝜌𝑁, 𝐵floor , 𝜌max 𝑁 .

2000 5000 1000 0

90

95

100

0

90

Accuracy@10 (%) kANNolo

hnswlib (original)

95

100

Accuracy@10 (%) hnswlib (fixed)

Figure 4: Accuracy@10 against throughput for kANNolo and

hnswlib, sweeping 𝑒𝑓search at 𝑘 = 10. original (dashed) inserts

in file order, fixed in permuted order.

of distance computations to within 4%. On SIFT1M, K ANN OLO and hnswlib in file order build in the same time, 25.56 and 25.63 seconds, and hnswlib in permuted order takes 27.95 (Table 8). On D RAGON its build time matches that of hnswlib once the insertion orders match. On D RAGON, the non-permuted hnswlib builds faster, in 1 498 seconds, but it produces a lower-quality graph. The two systems do not build the same graph, so this time is not directly comparable to KANNOLO’s. Throughput differs: the two are at parity at 𝑑 = 128 (0.99×), and KANNOLO is 1.15× faster at 𝑑 = 768. This difference tracks the advantage of the batched kernel over sequential scoring and is not a property of the graph. K ANN OLO stores more edges in less memory, 289.7M against 281.1M on D RAGON. hnswlib reserves all 2𝑀 ground-level slots for every node, half of them unused at the measured 48% occupancy, and keeps four auxiliary structures per element (≈86 bytes per vector). KANNOLO’s adjacency is variable-length, and its auxiliary arrays are sized to the upper levels. On resident memory, KANNOLO is ahead by 1.37× at 𝑑 = 128 and 1.07× at 𝑑 = 768.

5.5. Inverted Vector File Index An IVF index (Section 2.3.1) is at its core a pair: a centroid index that routes a query to a handful of clusters, and a quantizer that encodes the vectors inside them. Around that pair the design space is large. The centroid index can be flat or graph-based. The cluster vectors can be stored Page 20 of 39

Modular and Efficient Approximate Nearest Neighbors Search with Vectorium and kANNolo

at full precision or compressed. The cluster contents can be either the original vectors or centroid-relative residuals, under either metric. FAISS [5] exposes this space through a large family of classes (Section 2.5), one per combination. K ANN OLO represents every combination as one generic structure, IVF<DQ, D, CI>. The structure factors the space into three orthogonal axes, so a single search method and a single build method cover all of it, for every collection k-means can cluster. Like H NSW , the index names no representation. The first parameter, DQ, is the centroid and query dataset. Its encoder fixes the metric the index reports. The second, D, is the cluster dataset, which holds every vector in the collection. The vectors are reordered so that each cluster occupies a contiguous range, delimited by a cluster-offsets array. This CSR layout turns a cluster scan into sequential memory access. The third, CI, is the centroid index. It can be a flat index, an H NSW or even an IVF index. A second IVF gives a hierarchical, multi-level coarse quantizer through type composition alone. Both reported metrics decompose in the same form, whatever the centroid index. Let 𝑐 be the centroid of the cluster being scanned and 𝑟̂ the residual stored for a candidate, so that the candidate reconstructs as 𝑣̂ = 𝑐 + 𝑟̂. The inner-product decomposition is straightforward: ⟨𝑞, 𝑣⟩ ̂ = ⟨𝑞, 𝑐⟩ + ⟨𝑞, 𝑟̂⟩. The squared-Euclidean one is slightly more involved: ‖𝑞 − 𝑣‖ ̂ 2 = ‖(𝑞 − 𝑐) − 𝑟̂‖2 = ‖𝑞 − 𝑐‖2 + ‖̂𝑟‖2 − 2⟨𝑞 − 𝑐, 𝑟̂⟩ = ‖𝑞 − 𝑐‖2 + ‖̂𝑟‖2 + 2⟨𝑐, 𝑟̂⟩ − 2⟨𝑞, 𝑟̂⟩ . ⏟⏞⏟⏞⏟ ⏟⏞⏞⏞⏞⏞⏞⏟⏞⏞⏞⏞⏞⏞⏟ ⏟⏟⏟ centroid

offline correction

(1)

candidate

The two decompositions have the same structure if inner product is given a zero correction term. The centroid contribution is the query’s score against the centroid under the reported metric, and the centroid search has already produced it. The offline correction depends on the stored vector alone and is precomputed once per vector at build time. Under either metric, the only term that depends on both the query and the candidate is the dot product ⟨𝑞, 𝑟̂⟩. The cluster dataset therefore never evaluates the reported metric, and KANNOLO stores it in dot-product form throughout. Plain storage (no residuals) uses the same formula with 𝑐 = 0 and 𝑟̂ = 𝑣. ̂ The leading term is then the query-constant ‖𝑞‖2 , and KANNOLO does not compute it. It shifts every candidate equally and cannot change a ranking. Adding it once per query recovers the exact distance where needed. The metric the index reports fixes, at compile time, the sign of the dot product (−2 under squared-Euclidean distance and +1 under inner product) and whether a correction is present. The cluster encoder itself computes the offline correction, as ‖̂𝑟 + 𝑐‖2 − ‖𝑐‖2 . For scalar and PQ storage this equals the correction of Equation (1), with both terms acting on the decoded vector. For R AB ITQ and R AB ITQ-ext it is the encoder’s own estimate. On dense data, clusters can thus be stored as f32, any scalar type, PQ, R AB ITQ or R AB ITQ-ext codes, plain or as residuals, under either metric. S. Martinico et al.: Preprint submitted to Elsevier

Residuals pay off on dense data because the centroid already carries part of each vector. The residual can therefore be compressed more strongly at a smaller cost in accuracy. On sparse data the effect is reversed. A residual 𝑣 − 𝑐 has the nonzeros of both 𝑣 and 𝑐, so it takes more memory than the vector it replaces. Sparse clusters therefore store the vectors themselves.

The single-table residual decomposition. Under resid-

ual encoding, the candidate contribution ⟨𝑞, 𝑟̂⟩ does not depend on the cluster. For PQ storage, one 𝑚 × 256 table per query therefore serves every probed cluster. For any stored vector, the offline correction is a single number. KANNOLO and FAISS differ in when they materialize it. FAISS adopts the same decomposition. It enables the use_ precomputed_table mode whenever the table fits under its default 2 GB cap, as in our benchmarks. In this mode it builds the query table once and precomputes the correction during training, before any vector has been assigned to a cluster. Its inverted lists have no room for a per-vector scalar, so the correction must cover every (cluster, code) pair, as an nlist × 𝑚 × 256 array. Reading that array costs an 𝑂(𝑚 ⋅ 256) combine on every probed cluster, linear in nprobe. K ANN OLO assigns every vector’s cluster and code in one batch, before the index is assembled, so the correction does not have to stay generic. At build time it becomes a single stored f32 scalar per dataset vector. The encoder computes it from that vector’s own code, so it stays consistent with the ⟨𝑞, 𝑟̂⟩ that the search evaluates. Mixing the exact ‖𝑟‖2 with the approximate ⟨𝑞, 𝑟̂⟩ would introduce a 2⟨𝑞, 𝑟−𝑟̂⟩ ranking error. The query-side table is built once per query and reused across every probed cluster. Since the correction is already stored per vector, nothing is left to combine per cluster. At nprobe = 512 with 𝑚 = 64, FAISS’s combine alone performs ≈8.4M operations, against ≈2M multiply–adds for scoring the ≈31K candidates retrieved, which explains the growing advantage in Section 5.5.1.

Construction. Building an IVF index runs V ECTORIUM ’s

k-means, optionally trained on a random subsample of the collection. It then assigns every vector to its nearest centroid and reorders the dataset into the CSR layout above. The centroid encoder fixes the clustering metric. On sparse data, an optional refinement prunes the centroids during training. Centroid assignment, both during the k-means iterations and in the final pass, can switch from an exhaustive scan to an approximate search through H NSW . This pays off when the cluster count makes scanning tens of thousands of centroids per point the dominant cost. Which of the two is cheaper depends on the cluster count, the dimensionality and the collection size together, so the choice is a build parameter. The convenience entry point uses the H NSW path by default at 32 768 clusters and above. On sparse data, IVF is far slower than KANNOLO’s H NSW : on the first 1M passages of COCONDENSER, the best sparse IVF we built reaches 92.1%

Page 21 of 39

Modular and Efficient Approximate Nearest Neighbors Search with Vectorium and kANNolo Sift1M IVF

Sift1M IVFPQ64

Dragon IVF-Hnsw

Dragon IVFPQ128

600

4000

800

3000

QPS

3000 2000

400 200

1000

1000 0

600

400 2000

85

90

95

Accuracy@10 (%)

100

0

200

80

0

85

85

Accuracy@10 (%) kANNolo

90

Accuracy@10 (%) kANNolo (f16)

95

0

66

68

70

72

74

Accuracy@10 (%)

Faiss

Figure 5: Accuracy@10 against QPS for kANNolo and Faiss under matched parameters, single-query search, on the four configurations of Sift1M and Dragon.

Accuracy@10 at 75 QPS, against 93.4% at 2 053 QPS for H NSW .

5.5.1. Validation against FAISS We validated KANNOLO’s IVF against FAISS [5] to check that the generic, single-code-path design has no correctness or performance penalty. Conditions were matched: identical hardware, k-means and index parameters, accuracy targets, and single-query, single-threaded search timing on both systems. Four configurations span both metrics, both cluster-storage regimes and both kinds of centroid index. On SIFT1M (squared-Euclidean, flat centroid index) we test IVFFlat, with f32 and f16 cluster storage, and IVFPQ with residuals. On D RAGON (inner product, H NSW centroid index at 65 536 clusters) we test IVF and IVFPQ with residuals. Figure 5 reports accuracy against QPS for all four. Accuracy matches FAISS within 0.3 percentage points at every operating point across all four configurations. Under the same conditions, KANNOLO’s single-query latency is at or below FAISS’s. On SIFT1M, f32 cluster storage is 1.24–1.37× faster and f16 storage 1.6–1.9× faster. The margin grows with nprobe as the scan becomes memorybandwidth-bound. In the residual-PQ configuration the gap widens further, from 2.0× to 4.1× as nprobe grows. Here the single-table decomposition above removes FAISS’s percluster combine. On the inner-product D RAGON configurations, FAISS’s own residual formulation is already singletable, and the gap is a flatter 1.3–1.5×. Our decomposition offers no structural advantage on that metric, as the algebra above predicts. 5.5.2. Validation against the R AB ITQ reference The authors of R AB ITQ provide a reference implementation of IVF over R AB ITQ codes: the IVF index of their library. We compare it with KANNOLO’s IVF on one collection per metric, at 4 bits per component. We use SIFT1M under squared-Euclidean distance with 4 096 clusters, and SNOW FLAKE under inner product with 16 384 clusters. Both sides store centroid-relative codes. The reference takes its partition from outside the index. We therefore give it a KANNOLO k-means partition with the settings of KANNOLO’s

S. Martinico et al.: Preprint submitted to Elsevier

own build: 10 iterations, and on SNOW FLAKE spherical kmeans trained on a 4.2M-vector sample. Search runs one query at a time on one core, with 𝑘 = 10, following the conventions of Section 3.5. The reference library routes with exhaustively scanned centroids below 20 000 clusters, and with an H NSW over the centroids above. Both cluster counts here therefore use its flat routing. We compare it first with K ANN OLO ’s flat centroid index, and also report K ANN OLO with an H NSW centroid index. Table 9 reports the results. Accuracy matches. At every nprobe, the two flat-centroid indexes differ by at most 0.35 points on SIFT1M and 0.18 on SNOW FLAKE, and by at most 0.14 points on SIFT1M at 8 bits. Throughput is split. The split depends on where each design spends the time of a query: routing it to the clusters, or scanning their candidates. Before the reference reads a candidate’s full code, it scores the candidate with a 1-bit fast-scan estimate. It discards the candidate if the error bound of this estimate rules it out of the top 𝑘. Only the remaining candidates need a multibit estimate. KANNOLO has no such first pass and scores every candidate with the full code. Section 4.4.4 notes the same difference for the flat scan. The cost per candidate of K ANN OLO is therefore higher, and the reference is faster as soon as the candidates dominate the query time. On SIFT1M this happens from nprobe = 16, with a lead of 1.9× at 90% Accuracy@10. On SNOW FLAKE it happens from nprobe = 128, with a lead of up to 1.9× at the highest nprobe. The reference scans its centroids in f32, while K ANN OLO stores them in f16, which halves the bytes read. On SIFT1M the centroids are small, and this shows only at nprobe = 8, where KANNOLO is 1.12× faster. On SNOW FLAKE, 16 384 centroids of dimension 1 024 give the reference a minimum cost of about 4.4 ms per query, and K ANN OLO is 1.6× faster at nprobe = 16. Replacing the flat centroid index with an H NSW changes one type parameter and leaves the rest of the index unchanged. It removes most of that minimum cost. KANNOLO is then 5× faster than the reference at nprobe = 16 and 1.6× at 64, at a loss of 1.1 and 0.75 points of accuracy. At nprobe = 128 the scan dominates again. KANNOLO’s index files are 12% smaller on SIFT1M and 4% on SNOW FLAKE. Page 22 of 39

Modular and Efficient Approximate Nearest Neighbors Search with Vectorium and kANNolo Table 9 IVF over 4-bit RaBitQ-ext residual codes: kANNolo against the IVF of the RaBitQ authors’ library (ref.), singlequery search, Accuracy@10 (%) and QPS. flat and Hnsw name kANNolo’s centroid index; the reference scans its centroids exhaustively at these cluster counts. Accuracy@10 (%) nprobe

ref.

flat Hnsw

QPS ref.

flat Hnsw

Sift1M (𝓁2 , 4 096 clusters) 8 69.8 70.2 — 10 008 11 236 16 80.3 80.6 — 8 732 7 874 32 86.8 87.0 — 7 119 4 926 64 90.0 90.0 — 5 343 2 841 128 91.1 91.2 — 3 708 1 565 SnowFlake (IP, 16 384 clusters) 16 76.7 76.7 75.6 222 32 81.5 81.4 80.1 214 64 85.4 85.3 84.7 201 128 88.6 88.5 88.1 179 256 91.2 91.1 90.9 148

362 288 207 135 77

— — — — — 1 104 633 326 167 89

6. Experiments This section compares KANNOLO end to end against state-of-the-art systems on dense and sparse collections.

6.1. Experimental Setting For the dense domain, the comparison runs inside V IBE [107], a benchmark for vector indexes. V IBE fixes the collections, the ground truth, the parameter grids of every system it ships, and the timing loop. We run V IBE at commit 1b60239 with four patches, each applied to all systems in the same way. The first patch concerns build parallelism. V IBE pins construction to one core, while our builds use 64 threads. Parallel builds make this volume of benchmarks feasible. They also measure construction under realistic conditions: indexes are built in parallel in practice, and making construction scale with the number of available threads is part of designing an index. The other three patches concern measurement conventions and are described in Appendix A. For the sparse domain there is no equivalent of V IBE, so we use our own harness with the same conventions: construction on 64 threads under a per-collection budget, and search on a single core. Two sparse competitors, L SP and SINDI, can only be searched through their commandline tools. We therefore time every sparse system in the same way, as one search process that runs over the whole query set (Appendix A.2). Thus, dense and sparse systems cannot be compared with each other.

6.1.1. Collections Table 10 lists the four dense and three sparse collections. None of them was used in the design studies of Sections 4 and 5. Three of them use the M S M ARCO-V1 corpus of the development collections, but with a different encoder. S. Martinico et al.: Preprint submitted to Elsevier

On the dense side, AGN EWS [108] is a collection of news articles, embedded with mxbai-embed-large [109]. M SM ARCO-Q WEN is the passage corpus of M S M ARCOV 1 [98], drawn from Bing search results and embedded with Qwen3-Embedding-0.6B [110]. D PR is the Wikipedia passage collection of Dense Passage Retrieval [19], embedded with jina-embeddings-v5-text-nano [111]. On these three, V IBE holds out 1 000 random documents of the corpus as queries. LAION contains the CLIP [112] embeddings of one million images from LAION-400M [113]. It is searched with 1 000 embeddings of image captions, so its queries are outof-distribution. The set covers four encoder families and 0.77M to 21M points. M SM ARCO-Q WEN and D PR are the two largest collections that V IBE ships. They are closer to the scale of modern applications. They are also a challenge for some competitors: two of the four published degree settings of SYMPHONYQG cannot be built on this machine at 21M points, and a third exceeds the build budget (Section 6.3.1). On the sparse side, the vectors are learned term weights over a vocabulary, and queries are real user queries. M SM ARCO-SPLADE V3 encodes the M S M ARCO-V1 passages with SPLADE-v3 [30]. M SM ARCO-E-SPLADE encodes them with E-SPLADE [114], an efficiency-oriented SPLADE model whose query encoder is trained to produce queries with few nonzero entries.3 N Q is the BEIR [115] version of Natural Questions [116]: Wikipedia passages searched with Google queries, encoded with SPLADE COCONDENSER [28, 99]. M SM ARCO-E-SPLADE is the hardest case for a graph index because its queries are much shorter than its documents, making the search similar to an out-of-distribution task. Every collection is searched with all of its queries, against a ground truth we verified with an exhaustive scan.

6.1.2. Competitors The dense competitors cover the main designs of modern ANN indexes. Each design is represented by the strongest implementation we know of. hnswlib [64, 2] is the reference implementation of H NSW , over uncompressed f32 vectors. G LASS [8] has the same design as our KANNOLO-R ERANK: a graph searched over compressed vectors, with a rescoring stage on top. It is also the broadest library here, with four graph algorithms and fourteen quantizers, any of which can serve either stage. SVS-L EANV EC [74] combines VAMANA [67] with dimensionality reduction. It uses LVQ [34] codes in the first stage and a wider secondary dataset for rescoring. H NSW -R AB ITQ is the reference implementation of R AB ITQ [35, 36], which KANNOLO also uses in its first stage. Unlike our implementation, it stores additional metadata with each code to drive early exit. SYMPHONYQG [77] fuses the R AB ITQ codes into the adjacency lists of the graph. L ORANN [117] builds an IVF [4] index with a reducedrank regression per cluster. It is the only system that can be trained on the query distribution, from a learn set, and only L AION ships one. FAISS is not among the competitors. Its 3 The V-large checkpoints,

https://huggingface.co/naver/efficient-splade-V-large-doc and https://huggingface.co/naver/efficient-splade-V-large-query.

Page 23 of 39

Modular and Efficient Approximate Nearest Neighbors Search with Vectorium and kANNolo Table 10 Evaluation collections. nnz is the average number of nonzero entries. Size is the collection in raw format: f32 vectors for dense data, a 32-bit component and an f32 value per nonzero entry for sparse data. 𝐵 is the per-collection construction budget. † ships a learn set drawn from the query distribution. Collection

AgNews Ms Marco-v1 Dpr Laion †

Encoder

mxbai Qwen3 Jina CLIP

Ms Marco-v1 Splade-v3 Nq CoCondenser Ms Marco-v1 E-Splade

Type

𝑑

Metric

𝑛 Queries

nnz

Size (GiB) 𝐵 (s)

doc

query

dense dense dense dense

1 024 1 024 768 512

𝓁2 IP IP IP

769 382 8 840 823 20 969 760 1 000 448

1 000 1 000 1 000 1 000

— — — —

— — — —

2.9 33.7 60.0 1.9

300 2 400 5 400 300

sparse sparse sparse

— — —

IP IP IP

8 841 823 2 680 893 8 841 823

6 980 3 452 6 980

169.8 153.6 182.8

24.5 50.7 5.9

11.2 3.1 12.0

5 400 1 800 6 000

IVF indexes are compared with ours in Section 5.5.1, and its H NSW index was already slower than the first version of K ANN OLO [14]. We benchmark two configurations of KANNOLO. K ANN OLO -H NSW is our H NSW index built on f16 representations. KANNOLO-R ERANK is a two-stage index: a graph over quantized vectors retrieves a set of candidates, which are then reranked with a higherprecision representation. We use KANNOLO 0.15.1 with V ECTORIUM 0.11.1, hnswlib 0.8.0, G LASS 2.1.0, SYMPHONY QG 0.1.0, SVS-L EAN V EC 0.4.0, H NSW -R A B IT Q 0.1.1 and LORANN 0.4.6. G LASS, SYMPHONYQG and H NSW R AB ITQ are built from source at the commits pinned by V IBE. On sparse data, we compare against four inverted indexes: SEISMIC [10] and SEISMICWAVE [84] (both at 0.6.0), SINDI [81] (commit 4d22e55) and L SP [79, 80] (fa0ede5). We also compare against one graph index, PYANN [82] (cbb5c81). Our group develops SEISMIC and SEISMICWAVE on top of V ECTORIUM , so both also exercise the performance of V ECTORIUM under different indexing strategies. All the systems built on V ECTORIUM store sparse collections with the blocked encoder (Section 4.4.5). Our two sparse configurations mirror the dense ones. K ANN OLO -H NSW is an H NSW graph over the collection. The collection is stored with f16 values and u16 components. K ANN OLO -R ERANK splits each vector into two parts: the first part holds its largest entries, up to 𝑎% of its 𝓁1 norm, and the residual holds the remaining ones. The graph is built over the first part. The residual rescores the candidates through the residual mode of R ERANKINDEX (Section 5.1), and the two partial scores sum to the exact inner product.

6.1.3. Measures Accuracy@𝑘, throughput and build time are as defined in Section 3.2, with 𝑘 = 10 throughout. A query with no results has zero accuracy. The following conventions differ from those of the design studies.

Footprint. Footprint is the increase in the unique set size

(USS) of the search process across index construction: the private memory the process holds after the build minus the memory it held before. This measure charges an index S. Martinico et al.: Preprint submitted to Elsevier

for what it allocates. For every system here but one, this is the system’s own representation of the collection. The exception is L ORANN, as it does not copy the collection from V IBE, and it re-ranks its candidates against the array of the harness. The difference therefore does not include the memory its queries read. We add the raw f32 collection to its footprint, since any deployment would have to hold it. On sparse data, indexes are built and searched in separate processes, so for sparse data we report the resident memory of the search process. Both measures differ from those of the design studies (Section 3.2), so footprints are compared only within a table or figure.

Construction time. On dense data, construction time is the

time V IBE measures to build an index from the in-memory collection. On sparse data, it excludes the conversion of the collection into the input format of each library. In both cases it includes any processing of the data that an index relies on, such as encoding, pruning or document reordering.

Throughput per unit of memory. This measure is QPS per GiB of index. A system can gain throughput by using more memory, and Section 6.3 shows that the memory use of the systems here differs by up to an order of magnitude.

Repetition. Every search configuration is timed over two

passes of the query set, and we report the higher QPS of the two. Where a warm-up pass was run, it was discarded. The median spread between the two passes is 0.11–0.19% per collection, except on M SM ARCO-E-SPLADE (0.97%). Our graph builds are deterministic (Section 5.4.1), and rebuilds reproduce throughput within 1.0%. We therefore treat smaller differences as parity.

6.2. Methodology For each system and accuracy target, we report one configuration, chosen from a parameter grid within a construction budget.

6.2.1. Parameter Grids For dense benchmarks, grids come from V IBE verbatim. Where the shipped grid omits the documented default of a library, we add that default, which can only improve the Page 24 of 39

Modular and Efficient Approximate Nearest Neighbors Search with Vectorium and kANNolo

results of the competitors. Grids do not change across collections. The grid of KANNOLO-H NSW is 𝑀 ∈ {8, 16, 32}, 𝑒𝑓construction ∈ {200, 500}, and graph representation in {standard, permuted}. The grid of KANNOLO-R ERANK has 𝑀 ∈ {16, 32, 64}, the same values of 𝑒𝑓construction and the permuted representation (Section 6.3.3), and sweeps the firststage encoder over 1-bit R AB ITQ and R AB ITQ-ext at 2 and 4 bits. This gives 18 builds, against 36 for H NSW -R AB ITQ, the largest competitor grid. It is searched over the number of candidates, the beam width, 𝛼 and 𝜆, and with one-bit R AB ITQ also over the query width. The search ladders of V IBE are designed for 𝑘 = 100, and we run them at 𝑘 = 10. They do not limit the competitors. No selected configuration of a competitor is at an end of its search ladder, with two exceptions, both of L ORANN: at 90% on AGN EWS it uses its smallest reranking depth, and at 99% on L AION its largest number of probed clusters. On sparse data there is no shared benchmark, and only SEISMIC comes with configurations released by its authors for each collection, which we use unchanged. We define the grids of the other competitors. To avoid weakening a system because we know it less well than its authors, we start from the settings they recommend and extend an axis on a collection wherever a selected configuration is at its end. Axes defined by the authors of a system are reported at their limits. The sparse grids therefore differ across collections. For SINDI we sweep its two pruning parameters and its search parameter around the configuration recommended in its paper. For L SP, which provides no grid, we sweep around the settings of its documentation, with 8-bit weights as its authors recommend. For SEISMIC WAVE we reuse the construction settings of S EISMIC and sweep the number of 𝜅-NN neighbors expanded at search time. PYANN ships no construction grid, and its published search points reach at most 50% accuracy on N Q. We build it with our H NSW construction parameters and search it over a wider beam, with which it reaches 99% on N Q, 96.52% on M SM ARCO-SPLADE V3 and 90.92% on M SM ARCO-E-SPLADE. KANNOLO-H NSW is built with the standard representation, with 𝑀 ∈ {16, 32, 64} and 𝑒𝑓construction ∈ {200, 500, 1 000, 2 000}. Not every combination is built on every collection. It is searched over the beam width and 𝜆. KANNOLO-R ERANK is built with 𝑀 ∈ {16, 32}, 𝑒𝑓construction ∈ {1 000, 2 000} and pruning fraction 𝑎 ∈ {60, 70, 80}. It is searched over the number of candidates, the beam width and 𝜆. The same extension rule applies to the search axes of our two configurations. The full grid of every system on every collection is listed in the supplementary material, and code, results and sparse collections are publicly available.4

6.2.2. Build Budget V IBE bounds every configuration with the same time limit of 24 hours, and builds every index on a single thread [107]. Our builds use 64 threads, and our collections 4 https://github.com/TusKANNy/papers_reproducibility/tree/main/ information-systems-vectorium-kannolo

S. Martinico et al.: Preprint submitted to Elsevier

range from 0.77M to 21M points. We therefore replace this limit with a construction budget 𝐵 per collection, scaled with its size (Table 10). Every system has the same budget, and a configuration whose build exceeds 𝐵 is excluded from the comparison. The budgets reflect the build time of graph indexes. On sparse data we derive 𝐵 from the dense budgets, since a graph builds in similar time at the same construction parameters, and scale it with the number of nonzero entries per vector, which sets the cost of a distance computation. Inverted indexes build much faster by design (Section 6.3.2). To measure the effect of the budget, Appendix A.1 re-selects every reported dense configuration at 𝐵∕2, 2𝐵 and 4𝐵. For each system and each accuracy target 𝑡 ∈ {90, 95, 99}, Tables B.1 and B.2 in Appendix B report the configuration maximizing QPS subject to Accuracy@𝑘 ≥ 𝑡 over that system’s whole grid. The configuration behind every cell is listed in the supplementary material.

6.3. Experimental Results We compare KANNOLO with the competitors on dense data, then on sparse data, and conclude with ablations.

6.3.1. Dense Comparisons Search. Figure 6 draws the frontiers of the graph systems on the four collections. It leaves out L ORANN, the only system that does not build a graph, and SYMPHONYQG, the only one whose index is several times larger than the collection. Table B.1 reports every system at the three accuracy targets. Which system leads depends on the size of the collection. Below 1.1M points the field is tight and we do not always lead. On AGN EWS, G LASS is faster than KANNOLOR ERANK at all three targets. At 99% it reaches 6 085 QPS against our 4 585 (1.33× more), with 2.0× the footprint. L AION is the harder of the two small collections because its queries are out-of-distribution. Only our two configurations and L ORANN reach 99%. The other five systems stop between 95.84 and 98.74% over their whole grids. We lead at 95 and 99%, with 1 398 and 220 QPS. G LASS leads at 90%, with 3 417 QPS against our 2 987. At 8.8M and 21.0M points, the systems separate clearly and we lead. Excluding SYMPHONYQG, KANNOLOR ERANK is the fastest system at every target on both collections. On M SM ARCO-Q WEN it reaches 11 114, 7 837 and 2 896 QPS, ahead of the next system by 1.39, 1.29 and 1.15×. On D PR it reaches 7 095, 5 039 and 1 776 QPS, ahead by 1.19, 1.25 and 1.27×. Against hnswlib at 99%, KANNOLO-H NSW reaches 2.84× its throughput on AGN EWS, 2.70× on M SM ARCOQ WEN and 2.17× on D PR. KANNOLO-R ERANK reaches 3.04, 3.26 and 2.47×. In every case our indexes are 1.9 to 2.2× smaller. On L AION, hnswlib never reaches 99% and tops out at 98.74%. The gap combines several factors. On f32 data, our H NSW is on par with hnswlib on SIFT1M and faster on D RAGON (Section 5.4.2). Here it also stores the vectors in f16, and it uses distance-adaptive termination and the graph permutation. Table 11 adds these ingredients one at a time on AGN EWS and M SM ARCO-Q WEN, timing one query at a time as in V IBE. On AGN EWS, the fastest measured Page 25 of 39

Modular and Efficient Approximate Nearest Neighbors Search with Vectorium and kANNolo AgNews

MsMarco-Qwen

20000

10000

QPS

15000

Laion

6000

3000

4000

2000

2000

1000

7500

10000

5000

5000

2500

0

Dpr

90

95

0

100

Accuracy@10 (%) kANNolo-Hnsw

90

95

0

100

Accuracy@10 (%) kANNolo-Rerank

90

95

100

0

Accuracy@10 (%) hnswlib

Glass

Hnsw-RaBitQ

90

95

100

Accuracy@10 (%) Svs-LeanVec

Figure 6: Accuracy–throughput frontiers of the graph systems on the four dense collections. Each point is a system’s fastest configuration at the accuracy on its axis.

Table 11 Throughput at Accuracy@10 = 99% as the ingredients of kANNolo-Hnsw are added one at a time, every index built with 𝑀 = 16 and 𝑒𝑓construction = 200. Each value is interpolated between the two configurations of its row that bracket 99% on the accuracy–throughput frontier. Step is the ratio to the row above, and for the random-order row the ratio to file order. AgNews

MsMarco-Qwen

QPS

Step

QPS

Step

1 501 2 050 2 803 + early termination 4 201 + permutation 4 449

1.37× 1.37× 1.50× 1.06×

716 1 016 1 401 2 253 2 529

1.42× 1.38× 1.61× 1.12×

1 470 0.98×

634

0.89×

hnswlib, file order kANNolo, f32 + f16

hnswlib, random order

configuration of the last row that reaches 99% runs at 4 405 QPS, within 3% of the 4 277 QPS that V IBE measures for the same configuration. Inserting the points in random order makes no measurable difference for hnswlib on these collections. At equal 𝑒𝑓search , the two builds have the same throughput within 1%. Their accuracy differs by at most 0.41 points on AGN EWS and 1.06 on M SM ARCO-Q WEN, with no consistent sign, so the 0.89× on M SM ARCO-Q WEN is variation between builds. In total, the ingredients give 2.96 and 3.53× at these build parameters, against 2.84 and 2.70× with each system at its fastest configuration. We also compare against G LASS at the same f16 rescoring precision, which isolates the index from the choice of rescoring representation. In this setting KANNOLOR ERANK is 1.39, 1.27 and 1.27× faster on D PR, with 1.17 to 1.36× less memory. On M SM ARCO-Q WEN it is 1.33× faster at 99%, with 1.27× less memory. The 95% point of G LASS on D PR uses an f32 refine stage (73.44 GiB), and it is within the pass-to-pass spread of its f16 point.

Memory. Figure 7 places every system by footprint and

throughput at 95%, where all eight reach the target on every collection. Among the competitors, only H NSW -R AB ITQ builds a smaller index than KANNOLO-R ERANK with its f16 S. Martinico et al.: Preprint submitted to Elsevier

rescoring stage, on every collection and at every target it reaches. At 99%, KANNOLO-R ERANK holds 1.75 GiB on AGN EWS and 36.96 GiB on D PR, 0.62× the size of that collection. On D PR, hnswlib holds 74.64 GiB, G LASS 50.40, SVS-L EANV EC 43.23 and H NSW -R AB ITQ 22.58. SYMPHONYQG needs much more memory than the other systems. In throughput per GiB, both our configurations are ahead of it at every target on every collection. H NSW -R AB ITQ has the best trade-off among the competitors. It is ahead of our f16 configurations in throughput per GiB at 90 and 95% on AGN EWS, M SM ARCO-Q WEN and D PR, and at 99% on AGN EWS. It is not ahead on L AION, where it stops at 95.84%. In K ANN OLO the rescoring dataset is chosen independently of the index, so changing it requires no new index. If the f16 rescoring dataset is replaced with R A B ITQ-ext at 8 bits (Section 6.3.3), our configurations lead every system in throughput per GiB in 9 of the 12 collection–target pairs. R AB ITQ-ext rescoring leads in eight of them. In the ninth, L AION at 99%, R AB ITQ-ext rescoring does not reach the target and KANNOLO-H NSW leads. H NSW -R AB ITQ still leads on AGN EWS at 90 and 95% and on D PR at 90%. At 99% these configurations are 25–38% slower than with f16.

Construction. We do not lead on construction time (Ta-

ble B.1). At 99%, SVS-L EANV EC builds fastest wherever it reaches the target, and hnswlib slowest on the two large collections. Our builds are in the middle: across collections and targets they rank between second and sixth of the eight systems, most often fourth or fifth. At 99% on D PR, ours take 4 814 and 4 891 s, against 2 982 s for SVS-L EANV EC and 5 238 s for hnswlib. SVS-L EANV EC builds its VAMANA graph on vectors reduced to 384 dimensions and compressed with LVQ to one byte per dimension. The collections have 768 or 1 024 dimensions at four bytes each.

Summary. On the two small collections no system domi-

nates. G LASS is the fastest system after SYMPHONYQG on AGN EWS at every target, and the fastest on LAION at 90%, by 1.14 to 1.33× over us. We lead on L AION at 95 and 99%. On the two large collections, KANNOLO-R ERANK is faster than every system except SYMPHONYQG at every target, by 1.15 Page 26 of 39

Modular and Efficient Approximate Nearest Neighbors Search with Vectorium and kANNolo AgNews

MsMarco-Qwen

Dpr

QPS at 95%

20000 100/GiB

10000

20000 10000

1000/GiB

5000

0.5

1

2

5

10

Memory (GiB)

20

50

100

100/GiB

500 1/GiB

20

Memory (GiB)

kANNolo-Hnsw kANNolo-Rerank

hnswlib Glass

1000

10/GiB

200

500

2000

50

100

200

Memory (GiB)

SymphonyQG LoRANN

1000/GiB

2000

2000 1000 500

10/GiB

1000

100/GiB

20000 10000 5000

2000

5000

Laion

1000/GiB

1000/GiB

10000/GiB

50000

Hnsw-RaBitQ Svs-LeanVec

1

2

5

Memory (GiB) kANNolo-Rerank (ext-8)

Figure 7: Footprint (increase in USS, Section 6.1.3) against throughput at Accuracy@𝑘 ≥ 95% on the four dense collections, one point per system. Dotted lines are constant QPS per GiB of index. ext-8 marks kANNolo-Rerank with RaBitQ-ext at 8 bits as the rescoring dataset (Section 6.3.3). MsMarco-Splade

Nq

QPS

10000

MsMarco-E-Splade

6000

kANNolo-Hnsw kANNolo-Rerank Seismic SeismicWave PyANN Sindi Lsp

10000

7500 4000 5000

5000 2000

2500 0

90

95

100

Accuracy@10 (%)

0

90

95

Accuracy@10 (%)

100

0

90

95

100

Accuracy@10 (%)

Figure 8: Accuracy–throughput frontiers on the three sparse collections. Each point is a system’s fastest configuration at the accuracy on its axis.

to 1.39×. Its index is 0.57 to 0.62× the size of the collection, and its build time ranks third to sixth of the eight systems. SYMPHONYQG is the fastest system at every target on AGN EWS, M SM ARCO-Q WEN and D PR. This speed comes at a cost in memory, and the extra memory is proportionally larger than the extra throughput. On M SM ARCO-Q WEN it has 3.82, 3.62 and 2.39× the throughput of KANNOLOH NSW at the three targets, for 4.00× its footprint. On D PR it has 3.36 down to 1.81× the throughput of KANNOLOR ERANK, for 4.06× its footprint. At 99% it is therefore behind both our configurations in throughput per GiB: 80.1 against 134.1 and 151.9 on M SM ARCO-Q WEN, and 21.4 against 46.9 and 48.0 on D PR. Each node stores the R AB ITQ code of every neighbor next to its own vector. Each code is therefore replicated about as many times as the average degree. On D PR the index takes 149.99 GiB, 2.50× the size of the collection. At degree 64 the index grows to 240 GiB and the build to 6 485 s, over 𝐵. The two largest settings would need about 490 and 600 GiB, close to or above the memory of the machine, and were not run. The smallest setting was measured with its memory interleaved across both nodes, because it could not be placed on one. This makes its D PR row a lower bound. On L AION it is behind G LASS at 90% and behind us at 95%, and stops at 97.60%.

6.3.2. Sparse Comparisons S. Martinico et al.: Preprint submitted to Elsevier

Search. Figure 8 draws the frontiers on the three col-

lections, and Table B.2 reports every system at the three accuracy targets. The ranges below are taken over the whole frontier from 90 to 99%, so they include targets between those of the table. On M SM ARCO-SPLADE V3 and N Q the ranking is the same at every target. The two SEISMIC variants lead. SEISMICWAVE is 2.32–5.04× faster than KANNOLOR ERANK on M SM ARCO-SPLADE V3 and 1.57–2.50× on N Q. SEISMIC is 1.93–3.71× and 1.28–1.89× faster. The gap grows with the target. KANNOLO-R ERANK is faster than every other system. It is 1.26–1.87× faster than L SP wherever L SP reaches the target (it never reaches 99%). It is 1.46– 4.82× faster than SINDI, and 5.2–44.5× faster than PYANN, the other graph index. KANNOLO-H NSW is compared with K ANN OLO -R ERANK in Section 6.3.3. On M SM ARCO-E-SPLADE the inverted indexes pull further ahead. SEISMIC and SEISMICWAVE are within 9% of each other up to 95%, and SEISMICWAVE leads from 96 to 99%. They are 4.13–17.6× and 3.79–19.2× faster than KANNOLO-R ERANK, and again the gap grows with the target. L SP is third up to 98%, 1.25–2.95× faster than KANNOLOR ERANK, and never reaches 99%. KANNOLO-R ERANK is 2.90× faster than SINDI at 90% and stays ahead of it up to 98%. At 99%, SINDI is faster, with 463 against 255 QPS. K ANN OLO -H NSW does not reach 98%. P Y ANN reaches 90% at 68 QPS, 48.4× slower than KANNOLO-R ERANK. Page 27 of 39

Modular and Efficient Approximate Nearest Neighbors Search with Vectorium and kANNolo MsMarco-Splade

QPS at 95%

Nq

1000/GiB

10000 5000

MsMarco-E-Splade 1000/GiB

1000/GiB 100/GiB

5000

5000

2000 1000 500

kANNolo-Hnsw kANNolo-Rerank Seismic SeismicWave PyANN Sindi Lsp

10000 100/GiB

2000 100/GiB 10/GiB

1000

1000

500

200 100 5

10

20

RSS (GiB)

2000

500 2

5

10

RSS (GiB)

10/GiB

5

10

20

50

RSS (GiB)

Figure 9: Resident memory (RSS) of the search process against throughput at Accuracy@𝑘 ≥ 95% on the three sparse collections, one point per system. Dotted lines are constant QPS per GiB.

Only the largest build of PYANN, of degree 64, fails on the two M S M ARCO-V1 collections, where it returns no correct result: its graph loader computes the size of the graph in bytes with a signed 32-bit integer, which overflows at 8.8M documents. The PYANN results above come from its other builds. L SP is not robust here. Every configuration that keeps at most half of the query terms aborts (241 of 504), because the shortest queries are pruned to no term at all. Of the other configurations, 8 leave two queries unanswered. The configurations in Table B.2 are not affected.

Memory. Our indexes use the least memory on every

collection and at every target (Figure 9). On M SM ARCOSPLADE V3 they take 6.8–7.6 GiB, against 7.9–10.6 GiB for SEISMIC, 9.2–11.8 for SEISMICWAVE, 13.4 for SINDI, 23.0 for PYANN and 30.1 for L SP. On N Q they take 2.0–2.1 GiB, against 3.4–8.7 GiB for every other system. This is less than half the memory of SEISMIC, and enough to give K ANN OLO -R ERANK the highest throughput per GiB at all three targets. On M SM ARCO-E-SPLADE they take 7.0–8.0 GiB, against 8.4–11.5 GiB for SEISMIC, 10.1–11.6 for SEISMICWAVE, 13.2–14.2 for SINDI, 25.9 for PYANN and 31.3 for L SP.

Construction. Graph indexes take far longer to build than

inverted indexes, which build in seconds to minutes. The fastest builds are those of SINDI (26–104 s) and SEISMIC (63–156 s). SEISMICWAVE and L SP, which add a 𝜅-NN graph and a document reordering, take 292–1 844 s and 259– 1 118 s. Our builds take 481 to 5 428 s. The two-stage index builds its graph over the pruned vectors, and is therefore cheaper to build than the single-stage one. On M SM ARCOSPLADE V3, at the same 𝑀 = 16 and 𝑒𝑓construction = 1 000, it takes 2 268 s against 4 835 s. Within 𝐵 = 5 400 s it can therefore afford 𝑀 = 32 or 𝑒𝑓construction = 2 000, which push the single-stage build beyond 𝐵. On M SM ARCOSPLADE V3, the selected builds of PYANN are cheaper than those of our single-stage index, 2 054 against 4 835 s, because they use smaller construction parameters. They also search 7.6–20.3× more slowly.

Summary. On sparse data the inverted indexes are faster.

SEISMIC and SEISMICWAVE beat KANNOLO-R ERANK at S. Martinico et al.: Preprint submitted to Elsevier

every target on every collection they run on, and on M SM ARCO-E-SPLADE L SP does too, up to 98%. KANNOLOR ERANK is the fastest graph index we measured. Together with KANNOLO-H NSW , it is also the most compact index on every collection, but it takes much longer to build.

6.3.3. Ablations Two-stage over single-stage search. At 99%, K ANN OLO -R ERANK is 1.20× faster than K ANN OLO H NSW on M SM ARCO-Q WEN, 1.14× on D PR and 1.07× on AGN EWS and LAION. It gains more at the two lower targets: +41% and +30% on AGN EWS, +37% and +34% on M SM ARCO-Q WEN, +24% and +24% on D PR, +16% and +11% on L AION. On sparse data the gap is larger. KANNOLO-R ERANK is 1.63 to 3.97× faster than KANNOLO-H NSW on M SM ARCOSPLADE V3 and 1.42 to 1.61× on N Q, and at equal construction parameters it is also faster to build, because its graph is built over the first part only. Its selected builds can take longer (Table B.2), because within 𝐵 it can afford larger construction parameters. Memory is roughly the same, since the two parts together hold all the nonzero entries. On M SM ARCO-E-SPLADE it is 2.41 and 3.69× faster at 90 and 95%. On this collection KANNOLO-H NSW reaches neither 98 nor 99% within 𝐵.

Pruning fraction. The fastest configurations use 𝑎 = 60

or 70 on M SM ARCO-SPLADE V3, 𝑎 = 70 or 80 on N Q and 𝑎 = 60 on M SM ARCO-E-SPLADE. On M SM ARCO-ESPLADE, at 𝑀 = 16 and 𝑒𝑓construction = 1 000, 𝑎 = 60 is the fastest at 90, 95 and 98%, by 1–10% over 𝑎 = 70 and 20–36% over 𝑎 = 80. Pruning further hurts: on M SM ARCOSPLADE V3, at the same parameters, 𝑎 = 60 is 30–71% faster than 𝑎 = 40. With 𝑎 = 100 the residual is empty and the index reduces to KANNOLO-H NSW .

Graph representation. We compare

permuted with standard on AG N EWS at 𝑀 = 16 and 𝑒𝑓construction = 200, over 117 matched pairs of search parameters. permuted gains

+11.81% median throughput, and never less than +9.04%. It costs +20.9% build time and +4.2% USS. The permutation is a separate pass after graph construction, and it scales worse than the graph: at 21.0M points it takes 1 391–2 386 Page 28 of 39

Modular and Efficient Approximate Nearest Neighbors Search with Vectorium and kANNolo

s for 𝑀 = 8 to 32, which roughly doubles the build time. permuted is still faster at every accuracy target at 21.0M, and our D PR rows in Table B.1 use it. We leave the compressed layout of Section 5.3.4 out of the grid. At 𝑑 ∈ {768, 1 024} the adjacency is 4.5% of the index, so compressing it by 51.6% saves 2.3% overall. On sparse data we use only standard, because build times there are already high and the permutation would increase them further.

Early termination. For each collection and target we compare the fastest configuration with the fastest one with 𝜆 disabled. For KANNOLO-H NSW on the dense collections, 𝜆 gives 4–9% more throughput at 90% and 25–28% at 95%. On M SM ARCO-SPLADE V3 and N Q it gives 20–67% at every target, and on M SM ARCO-E-SPLADE 93 and 161% at 90 and 95%. At the other dense targets, our grid reaches the target only with 𝜆 on, because its beam widths are sized for search with early termination. On LAION this also holds beyond our grid. With 𝜆 off and 𝑒𝑓search swept up to 1 000, KANNOLO-H NSW plateaus at 98.36% and hnswlib at 98.74%. With early termination, KANNOLO-H NSW reaches 99.11% at 𝑒𝑓search = 12. KANNOLO-R ERANK depends less on 𝜆 on dense data: 𝜆 improves its fastest configuration by more than 1% in 4 of the 10 pairs measured both ways, by 3.5–9%, and on L AION it reaches 95 and 99% only with 𝜆 on. It depends more on 𝜆 on sparse data. There 𝜆 improves its fastest configuration in 5 of the 6 pairs measured both ways, by 7–34%, and KANNOLO-R ERANK reaches 99% only with 𝜆 on, on all three collections. Tail latency. Early termination widens the distribution of

per-query latency, as the spread of visited nodes in Section 5.3.1 predicts. We read the per-query latencies that V IBE records for the best pass of each reported configuration, over the 1 000 queries. A 99th percentile therefore rests on the ten slowest queries. On AGN EWS at 95%, 𝜆 lowers the median latency of KANNOLO-H NSW from 141 to 100 𝜇s but raises its 99th percentile from 226 to 308 𝜇s. At 99%, which our grid reaches only with 𝜆 on, the 99th percentile of KANNOLO-H NSW is 4.2–5.4× its median on AGN EWS, M SM ARCO-Q WEN and D PR, against 1.7–2.0× for hnswlib. It still stays below that of hnswlib: 957 against 1 211 𝜇s, 1 402 against 2 144 𝜇s and 2 180 against 2 310 𝜇s. K ANN OLO -R ERANK keeps the ratio between 1.5 and 2.4 on the same three collections. On L AION at 99%, a target no other graph system reaches, the ratio is 5.3 for KANNOLO-R ERANK and 8.0 for KANNOLO-H NSW . Sparse data shows the same pattern. Here we time each query inside the search tool, in a re-run of the reported configuration. The re-run reproduces the accuracy of Table B.2 exactly and its throughput within 3.2%. At 99%, the 99th percentile of KANNOLO-H NSW is 8.7× its median on M SM ARCOSPLADE V3 and 5.6× on N Q. KANNOLO-R ERANK keeps the ratio at 1.7 and 2.4×, below SEISMIC (3.9 and 3.8×) and SEISMICWAVE (3.2 and 3.1×). On M SM ARCO-E-SPLADE, where KANNOLO-H NSW does not reach 99%, the ratio is 4.4× for KANNOLO-R ERANK, against 5.0× for SEISMIC and S. Martinico et al.: Preprint submitted to Elsevier

4.1× for SEISMICWAVE. PYANN reaches 99% only on N Q. There it keeps the ratio at 1.5×, but its 99th percentile is 47 ms, 18× that of KANNOLO-R ERANK. L SP and SINDI can only be timed per pass and are not included.

Rescoring precision. The rescoring dataset dominates the

footprint of KANNOLO-R ERANK, so we test whether it needs to be near-exact. We replace f16 with R AB ITQ-ext at 8 bits per component on all four collections, keeping the selected builds and the search grid, so the candidate lists are unchanged. At the same search parameters the quantized stage is as fast as f16, within 2%. It loses, however, an almost constant amount of accuracy: a median of 0.25–0.49 points on AGN EWS, M SM ARCO-Q WEN and D PR, and 1.1 points on L AION. At 90% it keeps 0.95–0.99× the throughput of f16 on the first three collections, with 1.69–1.79× less memory. At 99% it costs 25, 25 and 38% of the throughput. Throughput per GiB is 1.42× better on AGN EWS (3 721 against 2 617), 1.21× better on M SM ARCO-Q WEN (184.6 against 151.9) and 1.04× better on D PR (49.9 against 48.0). On L AION, whose queries are out-of-distribution, it costs 24% of the throughput already at 95%, for a footprint 1.51× smaller, and stops at 97.73%, against 99.45% with f16. We report f16 because it keeps 1.3–1.6× the throughput at 99% and is the only one of the two that reaches 99% on L AION. Rescoring with encoded residuals, through the residual mode of R ERANKINDEX, gives mixed results: 58.9 QPS/GiB at 99% on D PR, but 99% is never reached on M SM ARCO-Q WEN.

The first-stage encoder does not transfer. The fastest

first-stage encoder differs across collections. One-bit R AB ITQ is selected on AGN EWS and M SM ARCO-Q WEN at almost every target, R AB ITQ-ext at two bits on D PR and at four bits on L AION. The share of the exact top-𝑘 that the top 100 candidates of each encoder cover explains this choice. One-bit R AB ITQ covers at least 99.7% on AGN EWS and M SM ARCO-Q WEN. On D PR it saturates at 99.11%, against 99.93% for R AB ITQ-ext at two bits. On LAION it collapses to 35–65% at every query width. The failure on D PR is not due to modality. Like M SM ARCO-Q WEN, D PR is text retrieval under inner product, but the encoder that is sufficient on M SM ARCO-Q WEN is not sufficient on D PR. We do not sweep PQ. At the code width of one-bit R AB ITQ (128 subspaces at 𝑑 = 1 024), it was 1.4–2.3× slower than R AB ITQ on AGN EWS. At 64 subspaces it covers at most 96.59% on every collection.

7. Conclusions and Future Work We presented V ECTORIUM and KANNOLO, two Rust libraries that factorize approximate nearest neighbor search into composable abstractions. Users compose them into indexes for dense, sparse and multi-vector representations, covering many trade-offs between construction time, memory, accuracy and speed. Because the compiler specializes every combination, the factoring adds no run-time overhead and still reaches state-of-the-art performance. Our Product Quantizer matches the accuracy of FAISS at competitive Page 29 of 39

Modular and Efficient Approximate Nearest Neighbors Search with Vectorium and kANNolo

speed, our R AB ITQ and H NSW match or outperform their reference implementations, and our IVF outperforms FAISS. End to end (Section 6), our two-stage H NSW is faster than every dense system except SYMPHONYQG on the two large collections, with an index of 0.57–0.62× the size of the collection. SYMPHONYQG needs about four times our memory. On sparse data, the systems built on V ECTORIUM lead. SEISMIC and SEISMICWAVE are the fastest systems on every collection, our indexes are the most compact on every collection, and our two-stage index is the fastest graph index we measured. Our indexes are not the fastest to build. On dense data they rank in the middle of the field, and on sparse data they take far longer to build than inverted indexes. In KANNOLO, the same components build both the dense and the sparse two-stage indexes, with no index code specific to either representation. This is useful in practice. The first-stage encoder that is sufficient differs from one collection to the next, and the rescoring dataset trades accuracy for memory. In both cases, a change requires only a new configuration and no new index. Future work goes in three directions. First, H NSW is so far the only algorithm built on the shared graph layer. We plan to add further graph indexes on it and to support dynamic collections, with insertions and deletions after construction. Second, we will keep both libraries current with new quantization and indexing techniques. We will also evaluate the capabilities that are already implemented but not measured here. Third, and most important, we plan to refine and extend the abstractions themselves, so that they cover designs that currently fall outside them. Specialized indexes such as SYMPHONYQG, SEISMIC and TACHIOM cannot yet be expressed with our framework’s abstractions, even though SEISMIC and TACHIOM are already built on V ECTORIUM . Bringing them inside the framework would let us study and develop these techniques further, and combine their components in ways not explored so far.

A. Harness and Measurement Details Three of the four V IBE patches of Section 6.1 are measurement conventions, applied to every system alike.

Page placement. Our machine has two NUMA nodes, and

the 64 build threads run on both. If an index were allocated on one node only, the threads of the other node would read remote memory at every access. Graph builds consist mostly of such reads and would take roughly 50% longer, while the build of L ORANN would be almost unchanged. This would distort the comparison of build costs across index families. Builds therefore spread their memory over both nodes, and the finished index is then migrated to the node of the core that runs the search. We discard any run whose migration cannot be verified from the process’s page map, and every reported configuration is at least 99% local to its node. With migrated pages, graph indexes lose 1.5% (ours) to 2.5% (hnswlib) of the throughput they reach with pages allocated locally from the start, so our reported ratio against hnswlib is inflated by S. Martinico et al.: Preprint submitted to Elsevier

about 1%. Only the index of SYMPHONYQG on D PR cannot be migrated (Section 6.3.1).

Memory measurement and allocator settings. Build

and search are recorded in two separate passes, so that the memory we report belongs to the index and not to the benchmark harness. Two settings of the C memory allocator, applied to every system alike, make this measurement reliable. The first caps the number of allocation arenas (MALLOC_ARENA_MAX); without it, a multithreaded build leaves memory reserved but unused. The second (MALLOC_MMAP_THRESHOLD_) serves allocations of up to 256 MiB from the allocator’s own pools. Without it, the 64 build threads queue on the lock that guards the process’s memory map, which slowed the construction of G LASS on D PR by a factor of 14. Rebuilt with both settings, the cheapest configuration of every system measured before the change moved by −4.6 to +1.6% in build time, so those measurements stand.

A.1. Budget Sensitivity We also measure configurations over the budget, and this makes the analysis possible. On dense data, each library’s grid is built cheapest first, in ascending order of its cost-monotone parameters. A build that exceeds 𝐵 is still completed and measured, up to a time limit of about 2𝐵, but it is excluded from the selection. All configurations at least as large as it in every one of those parameters are then skipped without being built, since they cannot be cheaper. A configuration that is larger in some parameters and smaller in others is still built. Because of this order, the set of skipped configurations does not depend on chance. We re-select every reported configuration on D PR at 𝐵∕2, 2𝐵 and 4𝐵 from the measured configurations, taking the maximum QPS at each accuracy target. On AGN EWS our whole grid already builds within 𝐵. On dense data, skipped configurations cannot be recovered, so the analysis bounds the sensitivity of the selection but not that of the full grid. Above 𝐵 the ranking changes only once: on D PR at 99%, H NSW -R AB ITQ overtakes SVS-L EANV EC from 2𝐵 on (1 109 against 993 QPS). The largest gains leave every position unchanged. SYMPHONYQG, already the fastest on D PR, goes from 3 213 to 6 030 QPS at 99%, and at 4𝐵 our configurations gain at most 0.1%. At 𝐵∕2, KANNOLO-H NSW loses 13–17%, again without changing its position, while K ANN OLO -R ERANK , G LASS, SYMPHONY QG, L O RANN and H NSW -R AB ITQ have no configuration that reaches 90% on D PR.

A.2. Sparse Timing Every sparse system is timed in the same way: one search process runs over the whole query set on a single core, with its memory local to that core’s node. The commandline tools of L SP and SINDI pass over the query set several times in one invocation. We take the end of each pass from their output, discard the first pass and report the best of the next two. We time every other system with the same rule, through a command-line driver of our own. Systems Page 30 of 39

Modular and Efficient Approximate Nearest Neighbors Search with Vectorium and kANNolo Table B.1 Throughput, footprint and construction time of each system’s fastest configuration at each accuracy target on the dense collections, two collections per cell; footprint is the increase in unique set size (USS), for LoRANN plus the raw collection (Section 6.1.3); – where a system does not reach the target, ⋅ where it was not run. On Laion the LoRANN row is its learned variant. + ext-8 is kANNolo-Rerank rescoring with RaBitQ-ext at 8 bits instead of f16 (Section 6.3.3). † : control configuration, not drawn in Figure 6: on Laion LoRANN without the learned variant. QPS System

90%

Memory (GiB)

95%

99%

AgNews / Laion kANNolo-Hnsw 11 665 /2 578 8 930 /1 254 4 277 kANNolo-Rerank 16 464 /2 987 11 595 /1 398 4 585 + ext-8 15 949 /2 835 11 450 /1 069 3 442 hnswlib 5 307 /1 504 3 624 /682 1 507 Glass 20 547 /3 417 14 374 /1 032 6 085 SymphonyQG 38 898 /3 116 31 356 /1 223 13 685 LoRANN 6 218 /1 696 2 989 /823 912 Hnsw-RaBitQ 13 090 /1 547 9 802 /550 3 188 Svs-LeanVec 10 357 /1 984 8 924 /617 3 608 LoRANN† ⋅ /1 593 ⋅ /761 ⋅ MsMarco-Qwen / Dpr kANNolo-Hnsw 8 092 /5 737 5 845 /4 065 kANNolo-Rerank 11 114 /7 095 7 837 /5 039 + ext-8 10 510 /7 046 6 953 /4 659 hnswlib 4 015 /3 355 2 686 /2 149 Glass 8 020 /5 098 6 058 /3 977 SymphonyQG 30 937/23 858 21 148/15 061 LoRANN 1 841 /550 656 /227 Hnsw-RaBitQ 7 568 /5 953 5 168 /4 047 Svs-LeanVec 5 395 /3 361 4 400 /1 981

90%

/206 /220 /– /– /– /– /211 /– /– /208

1.55 1.65 0.92 3.16 2.45 6.25 3.28 0.59 1.83 ⋅

Build (s)

95%

/1.15 /1.33 /0.87 /2.26 /1.64 /4.30 /2.05 /0.87 /1.59 /2.07

1.55 1.65 0.92 3.25 3.52 6.25 3.20 0.68 1.87 ⋅

/1.15 /1.33 /0.88 /2.31 /2.74 /4.30 /2.05 /0.88 /1.59 /2.05

99% 1.67 1.75 0.92 3.39 3.48 9.55 3.23 0.96 1.83 ⋅

90%

/1.15 /1.34 /– /– /– /– /2.07 /– /– /2.05

60 79 80 58 39 136 45 214 35 ⋅

95%

/127 /129 /130 /142 /57 /93 /153 /103 /88 /198

60 79 80 124 80 136 85 215 40 ⋅

99%

/127 /129 /146 /176 /63 /93 /153 /239 /88 /140

78 97 80 157 63 197 148 217 35 ⋅

/127 /144 /– /– /– /– /192 /– /– /140

2 405/1 562 17.71 /31.66 17.71 /32.46 17.92 /33.33 1 738/4 058 1 738/3 983 1 335/4 814 2 896/1 776 19.06 /36.15 19.06 /36.15 19.06 /36.96 1 622/3 938 1 622/3 938 1 622/4 891 2 165/1 106 10.70 /21.42 10.70 /21.42 11.73 /22.16 1 629/3 970 1 629/3 970 1 631/5 034 889 /719 36.76 /66.55 36.76 /68.40 38.86 /74.64 1 292/2 226 1 292/3 663 2 107/5 238 2 527/1 400 24.17 /42.40 24.17 /73.44 42.18 /50.40 964/3 949 964/4 311 1 384/4 358 5 742/3 213 71.67/149.99 71.67/149.99 71.67/149.99 1 944/4 331 1 944/4 331 1 944/4 331 135 /– 34.68 /61.73 34.68 /61.73 34.68 /– 1 990/4 996 1 990/4 996 1 990 /– 1 365 /986 10.61 /16.95 11.66 /22.58 10.61 /22.58 1 413/4 841 2 030/4 868 1 413/4 868 1 601 /993 22.45 /40.72 21.39 /40.72 22.45 /43.23 1 102/1 377 454/1 377 1 102/2 982

Table B.2 Throughput, resident memory of the search process and construction time of each system’s fastest configuration at each accuracy target on the sparse collections, two collections per cell in the first block; – where a system does not reach the target. QPS System

90%

95%

RSS (GiB) 99%

90%

95%

Build (s) 99%

90%

95%

99%

MsMarco-SPLADEv3 / Nq kANNolo-Hnsw 2 867/2 804 1 570/1 839 206 /563 6.75/1.96 6.75/2.05 6.75/2.05 4 835 /481 4 835/1 640 4 835/1 640 kANNolo-Rerank 4 682/4 506 3 183/2 811 818 /861 7.15/2.04 7.21/2.04 7.60/2.05 2 268 /994 4 160 /994 3 218/1 574 Seismic 9 877/5 755 7 091/3 804 3 032/1 631 7.94/4.39 8.66/4.39 10.59/6.40 98 /63 110 /63 139 /92 SeismicWave 11 407/7 462 8 113/5 345 4 121/2 151 9.18/6.37 9.68/4.74 11.82/6.37 858 /364 1 157 /292 1 742 /364 PyANN 357 /866 112 /395 – /33 23.00/6.42 23.00/6.42 –/6.42 2 054/1 248 2 054/1 248 –/1 248 Sindi 1 036/1 858 660/1 091 449 /588 13.35/3.63 13.35/3.98 13.35/3.98 92 /26 92 /28 92 /28 Lsp 3 402/2 761 2 040/1 857 – /– 30.08/3.42 30.08/8.74 – /– 1 045 /259 1 045 /270 – /– MsMarco-E-Splade kANNolo-Hnsw kANNolo-Rerank Seismic SeismicWave PyANN Sindi Lsp

1 363 3 286 13 585 12 457 68 1 132 4 118

504 1 862 10 302 9 885 – 836 3 201

– 255 4 501 4 908 – 463 –

that can be called in-process (ours, SEISMIC, SEISMICWAVE and PYANN) are tuned in-process, and only their fastest configurations at each target are timed from the command line. The systems built on V ECTORIUM use the blocked encoder. PYANN is called one query at a time. Its batch interface prunes each query against the largest weight in the whole batch, which changes its results. Calling it per query costs it 5.9% of its throughput on M SM ARCO-SPLADE V3.

B. Full Results Tables B.1 and B.2 report, for every system and collection of Section 6, the fastest configuration at each accuracy target, and Table B.3 the tail latency of the same configurations (Section 6.3.3). S. Martinico et al.: Preprint submitted to Elsevier

7.04 7.92 8.37 10.11 25.86 13.23 31.28

7.04 7.92 10.39 10.11 – 14.23 31.28

– 7.95 11.49 11.63 – 14.23 –

4 150 5 428 103 1 562 1 474 99 1 118

4 150 5 428 135 1 562 – 104 1 118

– 5 428 156 1 844 – 104 –

C. Design Studies This appendix reports the measurements behind design decisions stated in Sections 4 and 5: the distance kernels, the Product-Quantization table build, dataset addressing, the visited set and the neighbor-processing loop of the graph search.

C.1. Distance Kernels We compile each arithmetic kernel in two variants that differ only in their reduction operators. We run both on a handful of cache-resident vectors, which isolates pure compute. This measurement and the batched one below sweep synthetic vectors from 64 to 2 048 dimensions, and

Page 31 of 39

Modular and Efficient Approximate Nearest Neighbors Search with Vectorium and kANNolo Table B.3 Tail latency of the configurations of Tables B.1 and B.2 at 95 and 99% accuracy: 99th-percentile latency in µs, and its ratio to the median in parentheses. Dense: per-query times recorded by Vibe for the best pass, 1 000 queries. Sparse: each query timed inside the search tool over the timed passes of the reported configuration (in-process for PyANN); Lsp and Sindi can only be timed per pass and are not included. – where a system does not reach the target. AgNews

kANNolo-Hnsw kANNolo-Rerank + ext-8 hnswlib

Glass SymphonyQG LoRANN Hnsw-RaBitQ Svs-LeanVec

MsMarco-Qwen

Dpr

95%

99%

95%

99%

95%

308 (3.1) 121 (1.4) 119 (1.4) 446 (1.7) 102 (1.5) 47 (1.5) 362 (1.1) 170 (1.7) 156 (1.4)

957 (5.4) 488 (2.4) 450 (1.6) 1 211 (1.9) 228 (1.4) 121 (1.7) 1 193 (1.1) 442 (1.4) 346 (1.3)

392 (2.5) 211 (1.7) 226 (1.6) 641 (1.8) 237 (1.5) 77 (1.7) 1 805 (1.2) 349 (1.9) 336 (1.5)

1 402 (4.3) 519 (1.6) 1 439 (3.7) 2 144 (2.0) 756 (1.9) 290 (1.8) 8 265 (1.1) 1 195 (1.7) 894 (1.5)

486 (2.1) 301 (1.5) 326 (1.5) 691 (1.5) 333 (1.3) 103 (1.6) 5 132 (1.2) 426 (1.8) 710 (1.4)

MsMarco-SPLADEv3 95%

99%

kANNolo-Hnsw 2 082 (3.6) 28 446 (8.7) kANNolo-Rerank 574 (1.9) 2 047 (1.7) Seismic 425 (3.6) 1 073 (3.9) SeismicWave 327 (3.0) 683 (3.2) PyANN 14 864 (1.7) –

Nq 95% 1 106 (2.1) 623 (1.8) 760 (3.3) 492 (3.0) 4 040 (1.6)

we confirm them on real SIFT, COCONDENSER and higherdimensional embeddings (micro-benchmark conventions of Section 3.5). Since the loop shape is the same, the reassociation barrier alone decides whether the loop vectorizes; Section 4.3 reports the resulting speedups. The strict kernel retires several times more instructions per distance and stalls on its dependency chain: its IPC falls from 1.6 to below 0.9 as the reduction lengthens. On the exhaustive top-10 results, whose rankings are identical, the mean relative score deviation is ≈3 × 10−7 on D RAGON and exactly 0 on the integervalued SIFT vectors, whose sums stay exactly representable.

evaluation. The batch-of-six kernel of Section 4.3 cuts the instructions per f32 distance from 392 to 207 at 𝑑=1 024, while IPC stays flat or drops slightly. Its gain therefore comes from query reuse and not from instruction-level parallelism. Table C.1 reports the speedup over six single-candidate calls in three access regimes. In the cache-resident regime, a few vectors are rescored without leaving L1/L2. In the sequential regime, a scan runs over a buffer far larger than the last-level cache, as in exhaustive search. In the random regime, access is scattered over the same buffer, as among graph neighbors. With cache-resident vectors, the f32 gain grows with dimension, from 1.22× at 𝑑=64 to 1.73× at 𝑑=1 024, because there is more compute per query load to amortize. f16 gains less (1.14–1.45× up to 𝑑=2 048). Batching frees the single memory operand of the fused multiply–add for the vector, but an f16 value cannot use it without a separate widening instruction. The sequential scan is bandwidth-bound and the prefetcher already hides latency, so executing fewer instructions gains little. Under random access every vector is a cold miss, and the batched kernel wins by issuing six independent loads that overlap. The densified-query sparse kernel (COCONDENSER, nominal dimension ≈30k) has no query reuse: at each step, its six documents index six Batched

S. Martinico et al.: Preprint submitted to Elsevier

Laion 99%

95%

99%

2 180 (4.2) 2 342 (3.4) 26 235 (8.0) 850 (1.6) 1 445 (2.1) 18 187 (5.3) 1 777 (2.1) 1 894 (2.1) – 2 310 (1.7) 1 996 (1.3) – 1 024 (1.4) 1 341 (1.4) – 454 (1.4) 1 096 (1.3) – – 1 385 (1.1) 5 247 (1.1) 1 621 (1.6) 3 283 (1.9) – 1 367 (1.4) 2 020 (1.2) –

MsMarco-E-Splade 99%

95%

99%

7 778 (5.6) 11 287 (8.4) – 2 635 (2.4) 1 302 (2.6) 14 992 (4.4) 1 944 (3.8) 301 (3.6) 891 (5.0) 1 259 (3.1) 287 (3.2) 698 (4.1) 47 220 (1.5) – –

Table C.1 Speedup of the batch-of-six dot-product kernel over six singlecandidate calls, by access regime; micro-benchmark conventions of Section 3.5. Sparse: densified-query kernel on CoCondenser. f32

𝑑

f16

128 768 1 024 128 768 1 024 Sparse

cache 1.32 1.55 sequential 1.03 1.16 random 1.24 1.54

1.73 1.14 1.19 1.25 0.99 0.97 1.43 1.41 1.59

1.15 1.16 1.58

1.02 0.95 1.22

scattered positions of the query buffer. Batching therefore retires more instructions, and its only benefit comes from overlapping loads, so it helps under scattered access alone. This is the setting where the system uses it: the neighbors of a graph node, whose learned-sparse embeddings also share much of their support.

Other instruction sets. The same source vectorizes for

other instruction sets without changes. We rebuilt the binaries for AVX2 alone, with no AVX-512 instruction. The algebraic f32 kernels are then 13× faster than the strict ones at 𝑑=128 and 18–20× at higher dimensions, as in the native build. The batch-of-six kernel gains 1.20–1.84× on cache-resident vectors, against 1.12–1.76× for the native build in the same run. The PQ table build of Section 4.4.3 is 4.5–5.2× faster than the natural AoS build, against 4.4–5.2× for the native build. Absolute times differ by at most 7%, because on this CPU the compiler emits 256-bit vectors for these kernels even in the native build. On a Raspberry Pi 5 (Armv8, 128bit NEON), the algebraic f32 kernels are 4.4–6.3× faster than the strict ones, the batch-of-six kernel gains 1.65–1.70×, and the PQ table build is 2.6–4.2× faster than the natural AoS

Page 32 of 39

Modular and Efficient Approximate Nearest Neighbors Search with Vectorium and kANNolo Table C.2 Product-Quantization distance-table build, µs per table; microbenchmark conventions of Section 3.5, minimum over three builds. fm: fast-math; 𝑇 : SoA tile width, none for the untiled SoA build. AoS Shape

Metric plain

Sift1M 𝑑=128, 𝑚=16

𝓁2 IP

12.4 11.0

Dragon 𝑑=768, 𝑚=96

𝓁2 IP

SnowFlake 𝓁2 𝑑=1 024, 𝑚=128 IP

SoA fm none 𝑇 =8 𝑇 =32 9.8 9.8

2.6 2.5

3.8 3.4

2.4 2.2

74.2 58.9 66.1 58.8

14.5 22.9 14.5 20.4

14.2 13.3

99.3 78.6 88.7 78.5

26.7 33.5 26.4 30.3

20.4 19.8

build. The f16 kernels widen each value with Rust’s builtin f16 conversion, which the compiler lowers to NEON’s widening instruction. Their algebraic version is 3.5–3.8× faster than the strict one.

C.2. Product-Quantization Table Build We time the table build alone. We sweep 𝑚 over the full divisor range at 𝑑 ∈ {128, 768, 1 024} for both metrics. The 256𝑑 multiply–adds do not depend on 𝑚, so the sweep changes only the reduction length dsub and the loop nest around it. The natural AoS build has a single IEEEstrict accumulator, which keeps the reduction serial. Besides this build, we compile an AoS variant with the algebraic operators of Section 4.3, which separates the effect of the layout from that of the reassociation. All variants produce the same table up to f32 round-off, ≈5 × 10−7 (microbenchmark conventions of Section 3.5). Loops this short are sensitive to code placement. The alignment of branch targets alone changes the time of a kernel by up to 1.65×, through instruction-fetch stalls and not through its arithmetic. Each figure is therefore the minimum over three builds that differ only in code alignment. Table C.2 reports the three deployed shapes, all at dsub = 8, the least favorable end of the sweep (Section 4.4.3 gives their speedups). Across the 36 shapes, the gain over natural AoS ranges from 4.5× to 12.5×. Against the fast-math variant it never falls below 2.5× at any dsub ∈ {4, 8, 16}. Transposing alone, with neither build tiled, gives 3.4–5.1× on the deployed shapes, 3.4–13.8× across the sweep and 2.9–4.1× against the fast-math variant. Tiling with 𝑇 =32 adds a further 1.02–1.33×, while 𝑇 =8 is slower than no tiling. Fast-math itself speeds up AoS by 1.26× on 𝓁2 but only by 1.12× on IP. At dsub = 8 both AoS builds are bound by instruction issue, and the two fast-math kernels run in nearly the same time. The gain therefore measures the extra work of the strict kernel. The strict 𝓁2 kernel also pays a subtraction per element, which the vectorized form amortizes across lanes. AoS wins only at dsub ≳ 96, where a long reduction finally amortizes the horizontal sum. This regime means 𝑚 ≤ 8, a code of at most eight bytes for a 768-dimensional vector, and is of no practical interest. S. Martinico et al.: Preprint submitted to Elsevier

C.3. Identifier and Range Addressing We isolate the effect on the full COCONDENSER collection stored at f16. We perform five million accesses in random order, which is the scattered graph-neighbor pattern, and in storage order as a control. Each access is scored by the same QueryEvaluator. The three strategies fetch by identifier, by a range precomputed once, and by range with a prefetch issued 32 vectors ahead of use. We follow the microbenchmark conventions of Section 3.5. Under random access, an access costs 331 ns by identifier, 226 ns by range and 160 ns by range with prefetch. Range addressing removes the ≈1.8 last-level-cache misses per access caused by the offsets table (≈70 MB), which does not fit in cache. The prefetch removes ≈4 more misses. It executes slightly more instructions in fewer cycles, and IPC rises from 0.18 to 0.43. The gain therefore comes from hiding latency, and the amount of work does not decrease. In storage order the three strategies tie at 142 ns, so the benefit is specific to scattered access. Section 5.3.2 measures the same prefetch inside the graph beam search.

C.4. Visited Set Sizing the visited set. Only the ground level of H NSW

keeps a visited set. On the upper levels, the greedy walk moves strictly closer to the query and cannot revisit a node. On the ground level without early termination (𝜆 = 0), the search expands about 𝑒𝑓search nodes and scores the unvisited neighbors of each. We therefore expect 𝑛visited to follow 𝑒𝑓search and the degree of the graph, reduced by the overlap between neighbor lists. We instrumented KANNOLO to record, per query, 𝑛visited , the number of expanded nodes and the length of the neighbor lists they scan. We vary 𝑒𝑓search ∈ {10, 20, 40, 80, 160, 320} and 𝜆 ∈ {0, 0.05, 0.1, 0.15, 0.2, 0.3} on a 10M sample of SIFT1B, on D RAGON and on COCONDENSER, all indexed with 𝑀 = 32 and 𝑒𝑓construction = 200. At 𝜆 = 0 a query expands between 1.0 and 1.5 times 𝑒𝑓search nodes, at most 7 beyond it. Between 62% and 90% of the neighbors it scans are new, and this share decreases as 𝑒𝑓search grows. 𝑛visited therefore follows 𝑒𝑓search times the average ground-level degree. The estimate reads this degree from the index, because the pruning heuristic leaves most nodes below the 2𝑀 bound. Adaptive termination (Section 5.3.3) keeps expanding while a candidate is within a factor 1+𝜆 of the 𝑒𝑓search -th best distance. Its extra work therefore also scales with 𝑒𝑓search . At a fixed 𝜆, the number of expansions is multiplied by a factor that barely depends on 𝑒𝑓search : between 3.8 and 4.7 at 𝜆 = 0.15 on SIFT at 10M. The factor grows faster than linearly in 𝜆. At 𝑒𝑓search = 40 it is 1.6, 2.7, 4.6, 7.7 and 19.7 for 𝜆 from 0.05 to 0.3. Its growth rate depends on the collection: at 𝑒𝑓search = 20 and 𝜆 = 0.3 the factor is 19 on SIFT at 10M, 30 on COCONDENSER and 108 on D RAGON at 8.8M. The capacity estimate of Section 5.3.1 follows from these observations. The term 𝑑̄ ⋅ 𝑒𝑓search counts the neighbors of the expanded nodes, ignoring their overlap. The constant Page 33 of 39

Modular and Efficient Approximate Nearest Neighbors Search with Vectorium and kANNolo Table C.3 Sizing quality of the capacity estimate, 𝑒𝑓search from 10 to 320; the 𝜆 > 0 rows pool 𝜆 ∈ {0.05, 0.1, 0.15}. Rehash: share of queries rehashing at least once and at least twice. Overshoot: over the queries that do not rehash. Rehash (%)

Overshoot (× 𝑛visited )

any

≥𝟐

med.

p99

max

= 0 5.1 > 0 10.9 Dragon = 0 5.5 > 0 21.0 CoCondenser = 0 7.2 > 0 7.4

0.0 0.7 0.0 4.0 0.2 0.6

1.7 2.1 1.8 1.9 2.1 2.5

4.8 10.3 4.3 6.1 5.5 8.4

11.6 41.0 8.4 16.3 9.4 20.1

𝝀 Sift

200 covers the few expansions beyond 𝑒𝑓search . The 𝜆 term models the growth of the expansions. Suppose that the number of candidates within a threshold 𝑡 grows as 𝑡𝐷 , with 𝐷 the local intrinsic dimensionality of the data around the query [118, 119]. Then a factor 1+𝜆 on a squared-Euclidean score multiplies them by (1 + 𝜆)𝐷∕2 . Under inner product, scores are measured from zero and not from the query, so the effective exponent also grows with the ratio of the scores to their spread. Measured on the grid as ln(factor)∕ ln(1 + 𝜆), the exponent ranges from 8 to 11 on SIFT, 9 to 14 on COCONDENSER and 12 to 18 on D RAGON. We take the conservative value 8. It never over-predicts the growth, but it causes more rehashes where the exponent is larger. The offset 0.05, like the constant 200, gives small relaxations more room. We chose the constants by hand on these development collections and did not fit them. The estimate balances rehashing against oversized tables. Small sizing errors have little effect on search time, which is dominated by distance computations. We judge the estimate by two errors. A query that outgrows its table pays a rehash, which is a full copy of everything inserted so far. A table much larger than needed costs more to clear and spreads the inserts over more cache lines. For each query we record whether it rehashes. If it does not, we record its overshoot: the slots the table holds divided by the query’s 𝑛visited . The allocator rounds the capacity up to a power of two, so a table holds one to about two times the estimate. A good estimate keeps rehashes rare and the overshoot between 2× and 4×, with 10× as the exception. Figure C.1 compares the estimate with the measured mean, and Table C.3 reports both errors. Without early termination the estimate meets this target on all three collections: few queries rehash, almost none twice, and no table is oversized beyond 12×. Adaptive termination makes both errors more frequent, for two reasons. First, it widens the distribution of 𝑛visited within a configuration: its 99th percentile is a median 2.6 times the mean, against 1.9 at 𝜆 = 0. A table sized once per configuration therefore fits fewer queries well: more of them rehash, and more of them are oversized. Second, the growth with 𝜆 differs between collections. On SIFT, 1.0% of queries are oversized beyond 10×, against 0.03% on D RAGON and 0.35% on COCONDENSER. S. Martinico et al.: Preprint submitted to Elsevier

On D RAGON, a fifth of the queries rehash. At 𝜆 = 0.2 and 0.3 the difference grows: 4.2 and 7.6% of SIFT queries are oversized beyond 10×, and 48 and 92% of D RAGON queries rehash. A constant that served one collection better would serve the other worse. Even a perfect estimate of the mean would leave the easiest queries oversized and the hardest ones rehashing, since within one configuration the 99th percentile of 𝑛visited is 1.6 to 3.8 times its mean. BitVec against FxHashSet. With the capacity settled, we

can compare the two structures end to end. At 𝑛 = 108 and 𝑒𝑓search = 20, BitVec has the cheaper insert and hit (4.2 against 4.9 ns and 1.3 against 2.3 ns), because every hashset access computes a hash. It has the more expensive miss (2.8 against 2.1 ns). It also needs an 𝑂(𝑛) initialization on every query. Figure C.1 (right) isolates this initialization. The initialization of FxHashSet is bounded by the capacity: 36 to 62 ns for tables of 896 to 3 584 slots, at every collection size. BitVec clears 𝑛∕8 bytes per query. This takes 1.7 𝜇s at 𝑛 = 1M, 23 to 24 𝜇s at 10M, and 407 to 455 𝜇s at 100M, where it is a 12.5 MB memset. The end-to-end comparison crosses over at the same point as these two curves. At 𝑛 ≤ 1M, BitVec is ahead by at most 10 𝜇s for 𝑒𝑓search ≤ 80 and by up to 32 𝜇s at 𝑒𝑓search = 160. At 𝑛 = 8.8M, FxHashSet is faster in eight of nine configurations, by 2 to 69%; BitVec wins the ninth by 13%. At 𝑛 = 100M, BitVec is slower in every configuration by 414 to 424 𝜇s per query (Section 5.3.1).

C.5. Neighbor Processing Four collections span the memory footprints the loop can encounter, from a vector that fits inside the prefetch window to one far larger: SNOW FLAKE compressed with PQ128 (2 cache lines per vector), SIFT1M uncompressed at f16 (4), the sparse COCONDENSER (≈7.5 on average), and D RAGON uncompressed at f32 (48). Figure C.2 reports QPS with each feature enabled alone and with both together. Which single feature wins depends on the collection’s physical size and vector type. The four collections separate the effects of the two mechanisms. On the compact SNOW FLAKE codes, prefetching alone gives +35% QPS. The PQ codes span only 128 bytes, so one hint covers the whole vector. Batching adds nothing, since the PQ encoder has no single-pass batched kernel. The vectors of SIFT1M also fit inside the window, and as plain dense vectors they also use the batched kernel; the two contributions add up to 29%. On COCONDENSER both features help, prefetching more than batching. A batch of sparse vectors has six different sets of active components, so the gain of batching is limited to memory-level parallelism (Appendix C.1). Combining the two still costs nothing compared with the better of them. On D RAGON, whose vectors are much larger than the window, batching gives the largest gain at +34% QPS, and prefetching contributes +16%. Hardware performance counters (IPC, last-level-cache load misses and retired instructions) show the mechanism behind each gain. As for the isolated kernel of Appendix C.1, batching improves instruction efficiency. On D RAGON the Page 34 of 39

mean nvisited

λ=0

efsearch = 20

efsearch = 80

104

10

efsearch = 320

5

105

104 104

103

104

103 10

40

160

0

efsearch Sift

0.1

0.2

0.3

0

0.1

relaxation λ Dragon

0.2

0.3

0

relaxation λ

CoCondenser

measured

0.1

0.2

0.3

construction time (µs)

Modular and Efficient Approximate Nearest Neighbors Search with Vectorium and kANNolo

relaxation λ

capacity estimate

100 10 1 0.1 104 105 106 107 108

collection size n

FxHashSet

O(n)

BitVec

Figure C.1: Mean 𝑛visited (solid) against the capacity estimate (dashed) on Sift at 10M and Dragon and CoCondenser at 8.8M: against 𝑒𝑓search at 𝜆 = 0, then against 𝜆 at three fixed 𝑒𝑓search . Right: per-query construction cost of the two visited sets against 𝑛. SnowFlake (PQ128)

Sift1M (f16)

6000

CoCondenser

Dragon (f32) 4000

3000

15000

QPS

3000 4000

10000

2000

2000

5000

1000

0

70

75

Accuracy@10 (%)

80

0

2000

90

95

100

0

1000 90

Accuracy@10 (%) Baseline

Batch

95

100

Accuracy@10 (%) Prefetch

0

90

95

100

Accuracy@10 (%)

Batch + Prefetch

Figure C.2: Batching and prefetching against QPS on four collections spanning the range of per-vector memory footprints.

IPC stays at ≈3.3 whether batching is on or off, while retired instructions fall by ≈26%. Prefetching helps only where the CPU stalls on memory latency. On COCONDENSER the irregular sparse accesses hold the baseline IPC down to 2.1, and the 4-cache-line prefetch cuts LLC misses by 45% and raises QPS by 37%. On D RAGON we can test the cap itself: prefetching all 48 cache lines of a vector cuts LLC misses by 46% but gains only 0.5%. The out-of-order engine already hid the misses that this prefetch removes. The 1 536 PREFETCHT2 instructions it issues per expansion also lower the IPC enough to cancel the rest of the gain.

Declaration of generative AI and AI-assisted technologies in the manuscript preparation process

CRediT authorship contribution statement

Acknowledgements

Silvio Martinico: Conceptualization, Methodology, Software, Validation, Formal analysis, Investigation, Data curation, Visualization, Writing – original draft, Writing – review & editing. Franco Maria Nardini: Conceptualization, Methodology, Software, Validation, Formal analysis, Investigation, Data curation, Visualization, Writing – review & editing. Cosimo Rulli: Conceptualization, Methodology, Software, Validation, Formal analysis, Investigation, Data curation, Visualization, Writing – review & editing. Rossano Venturini: Conceptualization, Methodology, Software, Validation, Formal analysis, Investigation, Data curation, Visualization, Writing – review & editing.

We thank Leonardo Delfino and Domenico Erriquez for their work on the preliminary version of KANNOLO [14]. The authors were partially supported by the European Union’s Horizon Europe Research and Innovation program, under the RenAI project (GA No. 101291378) and by the ENRESFOOD project, funded under the FutureFoodS partnership. The project is co-funded by the Austrian Science Fund (FWF; 10.55776/KIN4661525), the Dutch Ministry of Agriculture, Fisheries, Food Security and Nature (BO-43-219-026), the Italian Ministry of University and Research (MUR), the German Federal Ministry of Research, Technology and Space (031B1717), and the European Union’s Horizon Europe research and innovation programme via FutureFoodS, the European Partnership for a Sustainable Future of Food Systems (GA No. 101136361).

S. Martinico et al.: Preprint submitted to Elsevier

During the preparation of this work the authors used Claude (Anthropic) in order to assist with parts of the software implementation, the setup of the experiments and the collection of their results, the preparation of the figures, and the preparation of the manuscript. After using this tool, the authors reviewed and edited the content as needed and take full responsibility for the content of the published article.

Page 35 of 39

Modular and Efficient Approximate Nearest Neighbors Search with Vectorium and kANNolo

The funders had no role in the study design, in the collection, analysis and interpretation of the data, in the writing of the article, or in the decision to submit it for publication.

References [1] Silvio Martinico, Franco Maria Nardini, Cosimo Rulli, and Rossano Venturini. Multivector reranking in the era of strong first-stage retrievers. In Advances in Information Retrieval: 48th European Conference on Information Retrieval (ECIR), pages 49–65, Berlin, Heidelberg, 2026. Springer-Verlag. ISBN 978-3-032-21323-5. doi: 10.1007/978-3-032-21324-2_4. [2] Yury Malkov et al. Hnswlib – fast approximate nearest neighbor search. https://github.com/nmslib/hnswlib, 2023. Version 0.8.0, accessed 2026-10-22. [3] Harsha Vardhan Simhadri, Ravishankar Krishnaswamy, Gopal Srinivasa, Suhas Jayaram Subramanya, Andrija Antonijevic, Dax Pryce, David Kaczynski, Shane Williams, Siddarth Gollapudi, Varun Sivashankar, Neel Karia, Aditi Singh, Shikhar Jaiswal, Neelam Mahapatro, Philip Adams, Bryan Tower, and Yash Patel. DiskANN: Graph-structured indices for scalable, fast, fresh and filtered approximate nearest neighbor search. https://github.com/ microsoft/DiskANN, 2023. Version 0.6.1, accessed 2026-10-22. [4] Jeff Johnson, Matthijs Douze, and Hervé Jégou. Billion-scale similarity search with GPUs. IEEE Transactions on Big Data, 7 (3):535–547, 2019. [5] Matthijs Douze, Alexandr Guzhva, Chengqi Deng, Jeff Johnson, Gergely Szilvasy, Pierre-Emmanuel Mazaré, Maria Lomeli, Lucas Hosseini, and Hervé Jégou. The FAISS library. IEEE Transactions on Big Data, 2025. [6] Magdalen Dobson Manohar, Zheqi Shen, Guy E. Blelloch, Laxman Dhulipala, Yan Gu, Harsha Vardhan Simhadri, and Yihan Sun. ParlayANN: Scalable and deterministic parallel graph-based approximate nearest neighbor search algorithms. In Proceedings of the 29th ACM SIGPLAN Annual Symposium on Principles and Practice of Parallel Programming (PPoPP), pages 270–285. ACM, 2024. doi: 10.1145/3627535.3638475. [7] Zheqi Shen, Zijin Wan, Jingbo Su, Yan Gu, and Yihan Sun. ANNLib: A development framework for efficient approximate nearest neighbor search. arXiv preprint arXiv:2607.17582, 2026. doi: 10.48550/ arXiv.2607.17582. [8] Zihao Wang. Graph library for approximate similarity search. https: //github.com/zilliztech/pyglass, 2025. Version 2.1.0, commit d2296ec, accessed 2026-10-22. [9] Leonid Boytsov and Bilegsaikhan Naidan. Engineering efficient and effective non-metric space library. In Similarity Search and Applications - 6th International Conference, SISAP 2013, volume 8199 of Lecture Notes in Computer Science, pages 280–293. Springer, 2013. doi: 10.1007/978-3-642-41062-8\_28. [10] Sebastian Bruch, Franco Maria Nardini, Cosimo Rulli, and Rossano Venturini. Efficient inverted indexes for approximate retrieval over learned sparse representations. In Proceedings of the 47th International ACM SIGIR Conference on Research and Development in Information Retrieval, 2024. [11] Sebastian Bruch, Franco Maria Nardini, Cosimo Rulli, and Rossano Venturini. Efficient sketching and nearest neighbor search algorithms for sparse vector sets. arXiv preprint arXiv:2509.24815, 2025. doi: 10.48550/arXiv.2509.24815. [12] Silvio Martinico, Franco Maria Nardini, Cosimo Rulli, and Rossano Venturini. Efficient multivector retrieval with token-aware clustering and hierarchical indexing. In Proceedings of the 49th International ACM SIGIR Conference on Research and Development in Information Retrieval, SIGIR ’26, pages 3982–3987, New York, NY, USA, 2026. Association for Computing Machinery. ISBN 9798400725999. doi: 10.1145/3805712.3809927. [13] Franco Maria Nardini, Luca Rizzo, Cosimo Rulli, and Rossano Venturini. Compact and efficient indexes for learned sparse retrieval.

S. Martinico et al.: Preprint submitted to Elsevier

arXiv preprint arXiv:2610.12300. Accepted at ICDE 2027, 2026. doi: 10.48550/arXiv.2610.12300. [14] Leonardo Delfino, Domenico Erriquez, Silvio Martinico, Franco Maria Nardini, Cosimo Rulli, and Rossano Venturini. kANNolo: Sweet and smooth approximate 𝑘-nearest neighbors search. In Advances in Information Retrieval: 47th European Conference on Information Retrieval (ECIR), pages 400–406, Berlin, Heidelberg, 2025. Springer-Verlag. ISBN 978-3-031-887161. doi: 10.1007/978-3-031-88717-8_29. [15] Stephen Robertson and Hugo Zaragoza. The probabilistic relevance framework: BM25 and beyond. Found. Trends Inf. Retr., 3(4):333– 389, 2009. ISSN 1554-0669. doi: 10.1561/1500000019. [16] D.G. Lowe. Object recognition from local scale-invariant features. In Proceedings of the Seventh IEEE International Conference on Computer Vision, volume 2, pages 1150–1157, 1999. doi: 10.1109/ ICCV.1999.790410. [17] Ashish Vaswani, Noam Shazeer, Niki Parmar, Jakob Uszkoreit, Llion Jones, Aidan N. Gomez, Łukasz Kaiser, and Illia Polosukhin. Attention is all you need. In Proceedings of the 31st International Conference on Neural Information Processing Systems (NeurIPS), pages 6000–6010, 2017. [18] Jacob Devlin, Ming-Wei Chang, Kenton Lee, and Kristina Toutanova. BERT: Pre-training of deep bidirectional transformers for language understanding. In Proceedings of the 2019 Conference of the North American Chapter of the Association for Computational Linguistics: Human Language Technologies, pages 4171– 4186, 2019. doi: 10.18653/v1/N19-1423. [19] Vladimir Karpukhin, Barlas Oguz, Sewon Min, Patrick Lewis, Ledell Wu, Sergey Edunov, Danqi Chen, and Wen-tau Yih. Dense passage retrieval for open-domain question answering. In Proceedings of the 2020 Conference on Empirical Methods in Natural Language Processing (EMNLP), pages 6769–6781, 2020. doi: 10.18653/v1/2020.emnlp-main.550. [20] Sheng-Chieh Lin, Akari Asai, Minghan Li, Barlas Oguz, Jimmy Lin, Yashar Mehdad, Wen-tau Yih, and Xilun Chen. How to train your dragon: Diverse augmentation towards generalizable dense retrieval. In Findings of the Association for Computational Linguistics: EMNLP 2023, pages 6385–6400. Association for Computational Linguistics, 2023. doi: 10.18653/v1/2023.findings-emnlp. 423. [21] Puxuan Yu, Luke Merrick, Gaurav Nuti, and Daniel Campos. ArcticEmbed 2.0: Multilingual retrieval without compromise. arXiv preprint arXiv:2412.04506, 2024. doi: 10.48550/arXiv.2412.04506. [22] Yang Bai, Xiaoguang Li, Gang Wang, Chaoliang Zhang, Lifeng Shang, Jun Xu, Zhaowei Wang, Fangshan Wang, and Qun Liu. SparTerm: Learning term-based sparse representation for fast text retrieval. arXiv preprint arXiv:2010.00768, 2020. doi: 10.48550/ arXiv.2010.00768. [23] Zhuyun Dai and Jamie Callan. Context-aware sentence/passage term importance estimation for first stage retrieval. arXiv preprint arXiv:1910.10687, 2019. doi: 10.48550/arXiv.1910.10687. [24] Zhuyun Dai and Jamie Callan. Context-aware term weighting for first stage passage retrieval. In Proceedings of the 43rd International ACM SIGIR Conference on Research and Development in Information Retrieval, pages 1533–1536, 2020. [25] Jimmy Lin and Xueguang Ma. A few brief notes on DeepImpact, COIL, and a conceptual framework for information retrieval techniques. arXiv preprint arXiv:2106.14807, 2021. doi: 10.48550/ arXiv.2106.14807. [26] Xueguang Ma, Ronak Pradeep, Rodrigo Nogueira, and Jimmy Lin. Document expansion baselines and learned sparse lexical representations for MS MARCO v1 and v2. In Proceedings of the 45th International ACM SIGIR Conference on Research and Development in Information Retrieval, pages 3187–3197, 2022. [27] Rodrigo Nogueira and Jimmy Lin. From doc2query to docTTTTTquery. Online preprint, https://cs.uwaterloo.ca/~jimmylin/ publications/Nogueira_Lin_2019_docTTTTTquery-v2.pdf, 2019.

Page 36 of 39

Modular and Efficient Approximate Nearest Neighbors Search with Vectorium and kANNolo [28] Thibault Formal, Benjamin Piwowarski, and Stéphane Clinchant. SPLADE: Sparse lexical and expansion model for first stage ranking. In Proceedings of the 44th International ACM SIGIR Conference on Research and Development in Information Retrieval, SIGIR ’21, pages 2288–2292, New York, NY, USA, 2021. Association for Computing Machinery. ISBN 9781450380379. doi: 10.1145/ 3404835.3463098. [29] Thibault Formal, Carlos Lassance, Benjamin Piwowarski, and Stéphane Clinchant. SPLADE v2: Sparse lexical and expansion model for information retrieval. arXiv preprint arXiv:2109.10086, 2021. doi: 10.48550/arXiv.2109.10086. [30] Carlos Lassance, Hervé Déjean, Thibault Formal, and Stéphane Clinchant. SPLADE-v3: New baselines for SPLADE. arXiv preprint arXiv:2403.06789, 2024. doi: 10.48550/arXiv.2403.06789. [31] Thibault Formal, Carlos Lassance, Benjamin Piwowarski, and Stéphane Clinchant. Towards effective and efficient sparse neural information retrieval. ACM Transactions on Information Systems, 42(5):1–46, 2024. [32] Omar Khattab and Matei Zaharia. ColBERT: Efficient and effective passage search via contextualized late interaction over BERT. In Proceedings of the 43rd International ACM SIGIR Conference on Research and Development in Information Retrieval, pages 39–48, 2020. ISBN 9781450380164. doi: 10.1145/3397271.3401075. [33] Keshav Santhanam, Omar Khattab, Jon Saad-Falcon, Christopher Potts, and Matei Zaharia. ColBERTv2: Effective and efficient retrieval via lightweight late interaction. In Proceedings of the 2022 Conference of the North American Chapter of the Association for Computational Linguistics: Human Language Technologies, pages 3715–3734, 2022. doi: 10.18653/v1/2022.naacl-main.272. [34] Cecilia Aguerrebere, Ishwar Singh Bhati, Mark Hildebrand, Mariano Tepper, and Theodore Willke. Similarity search in the blink of an eye with compressed indices. Proc. VLDB Endow., 16(11): 3433–3446, July 2023. ISSN 2150-8097. doi: 10.14778/3611479. 3611537. [35] Jianyang Gao and Cheng Long. RaBitQ: Quantizing highdimensional vectors with a theoretical error bound for approximate nearest neighbor search. Proc. ACM Manag. Data, 2(3), May 2024. doi: 10.1145/3654970. [36] Jianyang Gao, Yutong Gou, Yuexuan Xu, Yongyi Yang, Cheng Long, and Raymond Chi-Wing Wong. Practical and asymptotically optimal quantization of high-dimensional vectors in Euclidean space for approximate nearest neighbor search. Proc. ACM Manag. Data, 3(3):202:1–202:26, 2025. doi: 10.1145/3725413. [37] Hervé Jégou, Matthijs Douze, and Cordelia Schmid. Product quantization for nearest neighbor search. IEEE Transactions on Pattern Analysis and Machine Intelligence, 33(1):117–128, 2011. doi: 10.1109/TPAMI.2010.57. [38] Tiezheng Ge, Kaiming He, Qifa Ke, and Jian Sun. Optimized product quantization. IEEE Trans. Pattern Anal. Mach. Intell., 36 (4):744–755, 2014. doi: 10.1109/TPAMI.2013.240. [39] Yongjian Chen, Tao Guan, and Cheng Wang. Approximate nearest neighbor search by residual vector quantization. Sensors, 10(12): 11259–11273, 2010. doi: 10.3390/s101211259. [40] Iris Huijben, Matthijs Douze, Matthew Muckley, Ruud van Sloun, and Jakob Verbeek. Residual quantization with implicit neural codebooks. In ICML, 2024. [41] Artem Babenko and Victor Lempitsky. Additive quantization for extreme vector compression. In CVPR, 2014. [42] Ruiqi Guo, Philip Sun, Erik Lindgren, Quan Geng, David Simcha, Felix Chern, and Sanjiv Kumar. Accelerating large-scale inference with anisotropic vector quantization. In Proceedings of the 37th International Conference on Machine Learning (ICML), volume 119 of Proceedings of Machine Learning Research, pages 3887–3896. PMLR, 2020. URL https://proceedings.mlr.press/v119/guo20h. html. [43] Thomas Cover and Peter Hart. Nearest neighbor pattern classification. IEEE Transactions on Information Theory, 13(1):21–27, 1967.

S. Martinico et al.: Preprint submitted to Elsevier

[44] Piotr Indyk and Rajeev Motwani. Approximate nearest neighbors: towards removing the curse of dimensionality. In Proceedings of the Thirtieth Annual ACM Symposium on Theory of Computing, pages 604–613, 1998. [45] Mayur Datar, Nicole Immorlica, Piotr Indyk, and Vahab S Mirrokni. Locality-sensitive hashing scheme based on p-stable distributions. In Proceedings of the Twentieth Annual Symposium on Computational Geometry, pages 253–262, 2004. [46] Hakan Ferhatosmanoglu, Ertem Tuncel, Divyakant Agrawal, and Amr El Abbadi. Approximate nearest neighbor searching in multimedia databases. In Proceedings 17th International Conference on Data Engineering, pages 503–511. IEEE, 2001. [47] Martin Aumüller, Erik Bernhardsson, and Alexander Faithfull. ANN-Benchmarks: A benchmarking tool for approximate nearest neighbor algorithms. Information Systems, 87:101374, 2020. URL https://github.com/erikbern/ann-benchmarks. [48] Jon Louis Bentley. Multidimensional binary search trees used for associative searching. Communications of the ACM, 18(9):509–517, 1975. [49] Peter N. Yianilos. Data structures and algorithms for nearest neighbor search in general metric spaces. In Proceedings of the Fourth Annual ACM-SIAM Symposium on Discrete Algorithms (SODA), pages 311–321. Society for Industrial and Applied Mathematics, 1993. [50] Stephen M. Omohundro. Five balltree construction algorithms. Technical Report TR-89-063, International Computer Science Institute, 1989. [51] Paolo Ciaccia, Marco Patella, and Pavel Zezula. M-tree: An efficient access method for similarity search in metric spaces. In Proceedings of the 23rd International Conference on Very Large Data Bases (VLDB), pages 426–435. Morgan Kaufmann, 1997. [52] Marius Muja and David G. Lowe. Scalable nearest neighbor algorithms for high dimensional data. IEEE Transactions on Pattern Analysis and Machine Intelligence, 36(11):2227–2240, 2014. doi: 10.1109/TPAMI.2014.2321376. [53] Andrei Z. Broder. On the resemblance and containment of documents. In Proceedings of the Compression and Complexity of Sequences, pages 21–29, 1997. doi: 10.1109/SEQUEN.1997.666900. [54] Moses S. Charikar. Similarity estimation techniques from rounding algorithms. In Proceedings of the Thirty-Fourth Annual ACM Symposium on Theory of Computing (STOC), pages 380–388, 2002. doi: 10.1145/509907.509965. [55] Franz Aurenhammer. Voronoi diagrams—a survey of a fundamental geometric data structure. ACM Comput. Surv., 23(3):345–405, 1991. doi: 10.1145/116873.116880. [56] Gerard Salton. Automatic Information Organization and Retrieval. McGraw-Hill, 1968. [57] Philip Sun, David Simcha, Dave Dopson, Ruiqi Guo, and Sanjiv Kumar. SOAR: Improved indexing for approximate nearest neighbor search. In Advances in Neural Information Processing Systems, volume 36, pages 3189–3204. Curran Associates, Inc., 2023. doi: 10.52202/075280-0141. [58] Qi Chen, Haidong Wang, Mingqin Li, Gang Ren, Scarlett Li, Jeffery Zhu, Jason Li, Chuanjie Liu, Lintao Zhang, and Jingdong Wang. SPTAG: A Library for Fast Approximate Nearest Neighbor Search, 2018. URL https://github.com/Microsoft/SPTAG. [59] Ilias Azizi, Karima Echihabi, and Themis Palpanas. ELPIS: Graphbased similarity search for scalable data science. Proc. VLDB Endow., 16(6):1548–1559, 2023. doi: 10.14778/3583140.3583166. [60] Zeyu Wang, Peng Wang, Themis Palpanas, and Wei Wang. Graphand tree-based indexes for high-dimensional vector similarity search: Analyses, comparisons, and future directions. IEEE Data Eng. Bull., 46:3–21, 2023. [61] Mengzhao Wang, Xiaoliang Xu, Qiang Yue, and Yuxiang Wang. A comprehensive survey and experimental comparison of graph-based approximate nearest neighbor search. Proceedings of the VLDB Endowment, 14(11):1964–1978, 2021. doi: 10.14778/3476249. 3476255.

Page 37 of 39

Modular and Efficient Approximate Nearest Neighbors Search with Vectorium and kANNolo [62] Ilias Azizi, Karima Echihabi, and Themis Palpanas. Graph-based vector search: An experimental evaluation of the state-of-the-art. Proc. ACM Manag. Data, 3(1), February 2025. doi: 10.1145/ 3709693. URL https://doi.org/10.1145/3709693. [63] Yury Malkov, Alexander Ponomarenko, Andrey Logvinov, and Vladimir Krylov. Approximate nearest neighbor algorithm based on navigable small world graphs. Information Systems, 45:61–68, 2014. doi: 10.1016/j.is.2013.10.006. [64] Yu A. Malkov and D. A. Yashunin. Efficient and robust approximate nearest neighbor search using hierarchical navigable small world graphs. IEEE Trans. Pattern Anal. Mach. Intell., 42(4):824–836, apr 2020. ISSN 0162-8828. doi: 10.1109/TPAMI.2018.2889473. [65] Cong Fu, Chao Xiang, Changxu Wang, and Deng Cai. Fast approximate nearest neighbor search with the navigating spreading-out graph. Proceedings of the VLDB Endowment, 12(5):461–474, 2019. doi: 10.14778/3303753.3303754. [66] Godfried T. Toussaint. The relative neighbourhood graph of a finite planar set. Pattern Recognition, 12(4):261–268, 1980. doi: 10.1016/0031-3203(80)90066-7. [67] Suhas Jayaram Subramanya, Fnu Devvrit, Harsha Vardhan Simhadri, Ravishankar Krishnawamy, and Rohan Kadekodi. DiskANN: Fast accurate billion-point nearest neighbor search on a single node. Advances in Neural Information Processing Systems, 32, 2019. [68] Wei Dong, Moses Charikar, and Kai Li. Efficient 𝑘-nearest neighbor graph construction for generic similarity measures. In Proceedings of the 20th International Conference on World Wide Web (WWW), pages 577–586, 2011. doi: 10.1145/1963405.1963487. [69] Masajiro Iwasaki. Pruned bi-directed 𝑘-nearest neighbor graph for proximity search. In Similarity Search and Applications — 9th International Conference (SISAP), Lecture Notes in Computer Science, pages 20–33. Springer, 2016. doi: 10.1007/978-3-319-46759-7\_2. [70] Cong Fu, Changxu Wang, and Deng Cai. High dimensional similarity search with satellite system graph: Efficiency, scalability, and unindexed query compatibility. IEEE Transactions on Pattern Analysis and Machine Intelligence, 2021. doi: 10.1109/TPAMI. 2021.3067706. [71] Yousef Al-Jazzazi, Haya Diwan, Jinrui Gou, Cameron Musco, Christopher Musco, and Torsten Suel. Distance adaptive beam search for provably accurate graph-based nearest neighbor search. In Advances in Neural Information Processing Systems 38, pages 110958–110987, 2025. doi: 10.52202/085713-3346. [72] Tommaso Teofili and Jimmy Lin. Patience in proximity: A simple early termination strategy for HNSW graph traversal in approximate 𝑘-nearest neighbor search. In Advances in Information Retrieval – 47th European Conference on Information Retrieval (ECIR), pages 401–407. Springer, 2025. doi: 10.1007/978-3-031-88714-7_39. [73] Conglong Li, Minjia Zhang, David G. Andersen, and Yuxiong He. Improving approximate nearest neighbor search through learned adaptive early termination. In Proceedings of the 2020 ACM SIGMOD International Conference on Management of Data, pages 2539–2554, 2020. doi: 10.1145/3318464.3380600. [74] Mariano Tepper, Ishwar Singh Bhati, Cecilia Aguerrebere, Mark Hildebrand, and Theodore Willke. LeanVec: Searching vectors faster by making them fit. Transactions on Machine Learning Research, 2024. ISSN 2835-8856. [75] Masajiro Iwasaki and Daisuke Miyazaki. Optimization of indexing based on k-nearest neighbor graph for proximity search in highdimensional data. arXiv preprint arXiv:1810.07355, 2018. doi: 10.48550/arXiv.1810.07355. [76] Kohei Sugawara, Hayato Kobayashi, and Masajiro Iwasaki. On approximately searching for similar word embeddings. In Katrin Erk and Noah A. Smith, editors, Proceedings of the 54th Annual Meeting of the Association for Computational Linguistics (Volume 1: Long Papers), pages 2265–2275, Berlin, Germany, August 2016. Association for Computational Linguistics. doi: 10.18653/v1/P16-1214. URL https://aclanthology.org/P16-1214/.

S. Martinico et al.: Preprint submitted to Elsevier

[77] Yutong Gou, Jianyang Gao, Yuexuan Xu, and Cheng Long. SymphonyQG: Towards symphonious integration of quantization and graph for approximate nearest neighbor search. Proc. ACM Manag. Data, 3(1), February 2025. doi: 10.1145/3709730. URL https: //doi.org/10.1145/3709730. [78] Antonio Mallia, Torsten Suel, and Nicola Tonellotto. Faster learned sparse retrieval with block-max pruning. In Proceedings of the 47th International ACM SIGIR Conference on Research and Development in Information Retrieval, pages 2411–2415, 2024. [79] Parker Carlson, Wentai Xie, Shanxiu He, and Tao Yang. Dynamic superblock pruning for fast learned sparse retrieval. In Proceedings of the 48th International ACM SIGIR Conference on Research and Development in Information Retrieval, pages 3004–3009, 2025. doi: 10.1145/3726302.3730183. [80] Parker Carlson, Wentai Xie, Rohil Shah, and Tao Yang. Efficient learned sparse retrieval with lightweight superblock pruning. In The 49th International ACM SIGIR Conference on Research and Development in Information Retrieval (SIGIR ’26). ACM, 2026. [81] Ruoxuan Li, Xiaoyao Zhong, Jiabao Jin, Peng Cheng, Wangze Ni, Zhitao Shen, Wei Jia, Xiangyu Wang, Heng Tao Shen, and Jingkuan Song. SINDI: An efficient index for sparse vector approximate maximum inner product search. In 2026 IEEE 42nd International Conference on Data Engineering (ICDE), pages 2793–2806. IEEE, 2026. doi: 10.1109/ICDE65706.2026.00208. [82] Zihao Wang. PyAnns: A python library for approximate nearest neighbor search. https://github.com/hhy3/pyanns, 2025. Commit cbb5c81 (2025-06-17), doi: 10.5281/zenodo.15331856, accessed 2026-10-22. [83] Harsha Vardhan Simhadri, Martin Aumüller, Matthijs Douze, Dmitry Baranchuk, Amir Ingber, Edo Liberty, George Williams, Ben Landrum, Magdalen Manohar, Mazin Karjikar, Laxman Dhulipala, Meng Chen, Yue Chen, Rui Ma, Kai Zhang, Yuzheng Cai, Jiayang Shi, Weiguo Zheng, Yizhuo Chen, Jie Yin, and Ben Huang. Results of the Big ANN: NeurIPS’23 competition. In Advances in Neural Information Processing Systems 38, Datasets and Benchmarks Track. Curran Associates, Inc., 2025. doi: 10.52202/085713-2303. [84] Sebastian Bruch, Franco Maria Nardini, Cosimo Rulli, and Rossano Venturini. Pairing clustered inverted indexes with 𝜅-nn graphs for fast approximate retrieval over learned sparse representations. In Proceedings of the 33rd ACM International Conference on Information and Knowledge Management (CIKM), pages 3642–3646, New York, NY, USA, 2024. Association for Computing Machinery. ISBN 9798400704369. doi: 10.1145/3627673.3679977. [85] Daniel Severo, Giuseppe Ottaviano, Matthew Muckley, Karen Ullrich, and Matthijs Douze. Lossless compression of vector IDs for approximate nearest neighbor search. arXiv preprint arXiv:2501.10479, 2025. doi: 10.48550/arXiv.2501.10479. [86] Justin Zobel and Alistair Moffat. Inverted files for text search engines. ACM Comput. Surv., 38(2):6–es, July 2006. ISSN 03600300. doi: 10.1145/1132956.1132959. [87] Laxman Dhulipala, Igor Kabiljo, Brian Karrer, Giuseppe Ottaviano, Sergey Pupyrev, and Alon Shalita. Compressing graphs and indexes with recursive graph bisection. In Proceedings of the 22nd ACM SIGKDD International Conference on Knowledge Discovery and Data Mining, pages 1535–1544. ACM, 2016. doi: 10.1145/2939672. 2939862. [88] J. Mackenzie, M. Petri, and A. Moffat. Faster index reordering with bipartite graph partitioning. In Proc. SIGIR, pages 1910–1914, 2021. [89] Daniel Lemire, Nathan Kurz, and Christoph Rupp. Stream VByte: Faster byte-oriented integer compression. Information Processing Letters, 130:1–6, 2018. [90] Francesco Benocci, Leonard Dahlmann, Francesco Gargiulo, Tobias Heuer, Silvio Martinico, Franco Maria Nardini, Cosimo Rulli, and Rossano Venturini. Lossless compression of proximity graph indexes for efficient approximate nearest neighbors search. In Proceedings of the 35th ACM International Conference on Information and Knowledge Management (CIKM ’26). Association for Computing Machinery, 2026. ISBN 979-8-4007-2539-5. doi: 10.1145/3799682.

Page 38 of 39

Modular and Efficient Approximate Nearest Neighbors Search with Vectorium and kANNolo 3840116. [91] Benjamin Coleman, Santiago Segarra, Alex Smola, and Anshumali Shrivastava. Graph reordering for cache-efficient near neighbor search. In Proceedings of the 36th International Conference on Neural Information Processing Systems, NIPS ’22, Red Hook, NY, USA, 2022. Curran Associates Inc. ISBN 9781713871088. [92] Kakao. TOROS N2: Lightweight approximate nearest neighbor library. https://github.com/kakao/n2, 2020. Version 0.1.7, released 2020-10-16, accessed 2026-11-01. [93] Jianguo Wang, Xiaomeng Yi, Rentong Guo, Hai Jin, Peng Xu, Shengjun Li, Xiangyu Wang, Xiangzhou Guo, Chengming Li, Xiaohai Xu, Kun Yu, Yuxing Yuan, Yinghao Zou, Jiquan Long, Yudong Cai, Zhenxiang Li, Zhifeng Zhang, Yihua Mo, Jun Gu, Ruiyi Jiang, Yi Wei, and Charles Xie. Milvus: A purpose-built vector data management system. In Proceedings of the 2021 International Conference on Management of Data (SIGMOD), 2021. doi: 10. 1145/3448016.3457550. [94] Qdrant. Qdrant: Vector search engine. https://github.com/qdrant/ qdrant, 2026. Version 1.20.0, released 2026-10-22, accessed 202611-01. [95] Vespa.ai. Vespa: The AI search platform. https://github.com/ vespa-engine/vespa, 2026. Version 8.778.45, released 2026-10-27, accessed 2026-11-01. [96] Elastic. Elasticsearch. https://github.com/elastic/elasticsearch, 2026. Version 9.5.5, released 2026-10-06, accessed 2026-11-01. [97] Weaviate. Weaviate. https://github.com/weaviate/weaviate, 2026. Version 1.40.3, released 2026-10-28, accessed 2026-11-01. [98] Tri Nguyen, Mir Rosenberg, Xia Song, Jianfeng Gao, Saurabh Tiwary, Rangan Majumder, and Li Deng. MS MARCO: A human generated machine reading comprehension dataset. In Proceedings of the Workshop on Cognitive Computation: Integrating neural and symbolic approaches, volume 1773 of CEUR Workshop Proceedings. CEUR-WS.org, 2016. URL https://ceur-ws.org/Vol-1773/ CoCoNIPS_2016_paper9.pdf. [99] Thibault Formal, Carlos Lassance, Benjamin Piwowarski, and Stéphane Clinchant. From distillation to hard negative sampling: Making sparse neural IR models more effective. In Proceedings of the 45th International ACM SIGIR Conference on Research and Development in Information Retrieval, pages 2353–2359, 2022. doi: 10.1145/3477495.3531857. [100] Mark Kac. Random walk and the theory of Brownian motion. The American Mathematical Monthly, 54(7):369–391, 1947. [101] Keshav Santhanam, Omar Khattab, Christopher Potts, and Matei Zaharia. PLAID: An efficient engine for late interaction retrieval. In Proceedings of the 31st ACM International Conference on Information & Knowledge Management (CIKM), pages 1747–1756, 2022. doi: 10.1145/3511808.3557325. [102] Franco Maria Nardini, Cosimo Rulli, and Rossano Venturini. Efficient multi-vector dense retrieval with bit vectors. In Proceedings of the 46th European Conference on Information Retrieval (ECIR), pages 3–17, 2024. doi: 10.1007/978-3-031-56060-6_1. [103] Jan Luca Scheerer, Matei Zaharia, Christopher Potts, Gustavo Alonso, and Omar Khattab. WARP: An efficient engine for multivector retrieval. In Proceedings of the 48th International ACM SIGIR Conference on Research and Development in Information Retrieval, pages 2504–2512, 2025. doi: 10.1145/3726302.3729904. [104] Zheng Bian, Man Lung Yiu, and Bo Tang. IGP: Efficient multivector retrieval via proximity graph index. In Proceedings of the 48th International ACM SIGIR Conference on Research and Development in Information Retrieval, pages 2524–2533, 2025. doi: 10.1145/3726302.3730004. [105] Peter Elias. Efficient storage and retrieval by content and address of static files. Journal of the ACM (JACM), 21(2):246–260, 1974. [106] Giulio Ermanno Pibiri and Rossano Venturini. Techniques for inverted index compression. ACM Computing Surveys (CSUR), 53 (6):125:1–125:36, 2021. doi: 10.1145/3415148. [107] Elias Jääsaari, Ville Hyvönen, Matteo Ceccarello, Teemu Roos, and Martin Aumüller. VIBE: Vector index benchmark for embeddings.

S. Martinico et al.: Preprint submitted to Elsevier

Journal of Data-centric Machine Learning Research, 2026. URL https://openreview.net/forum?id=6Sx6ra1QLz. [108] Xiang Zhang, Junbo Zhao, and Yann LeCun. Character-level convolutional networks for text classification. In Advances in Neural Information Processing Systems, volume 28, pages 649–657, 2015. [109] Sean Lee, Aamir Shakir, Darius Koenig, and Julius Lipp. Open source strikes bread – new fluffy embedding model. Mixedbread blog, https://www.mixedbread.com/blog/mxbai-embed-large-v1, 2024. Accessed 2026-10-22. [110] Yanzhao Zhang, Mingxin Li, Dingkun Long, Xin Zhang, Huan Lin, Baosong Yang, Pengjun Xie, An Yang, Dayiheng Liu, Junyang Lin, Fei Huang, and Jingren Zhou. Qwen3 Embedding: Advancing text embedding and reranking through foundation models. arXiv preprint arXiv:2506.05176, 2025. doi: 10.48550/arXiv.2506.05176. [111] Mohammad Kalim Akram, Saba Sturua, Nastia Havriushenko, Quentin Herreros, Michael Günther, Maximilian Werk, and Han Xiao. jina-embeddings-v5-text: Compact and robust text embedding models using task-targeted distillation. In Proceedings of the 49th International ACM SIGIR Conference on Research and Development in Information Retrieval, pages 4454–4458. ACM, 2026. doi: 10.1145/3805712.3808455. [112] Alec Radford, Jong Wook Kim, Chris Hallacy, Aditya Ramesh, Gabriel Goh, Sandhini Agarwal, Girish Sastry, Amanda Askell, Pamela Mishkin, Jack Clark, Gretchen Krueger, and Ilya Sutskever. Learning transferable visual models from natural language supervision. In Proceedings of the 38th International Conference on Machine Learning (ICML), volume 139 of Proceedings of Machine Learning Research, pages 8748–8763, 2021. [113] Christoph Schuhmann, Richard Vencu, Romain Beaumont, Robert Kaczmarczyk, Clayton Mullis, Aarush Katta, Theo Coombes, Jenia Jitsev, and Aran Komatsuzaki. LAION-400M: Open dataset of CLIP-filtered 400 million image-text pairs. arXiv preprint arXiv:2111.02114, 2021. doi: 10.48550/arXiv.2111.02114. [114] Carlos Lassance and Stéphane Clinchant. An efficiency study for SPLADE models. In Proceedings of the 45th International ACM SIGIR Conference on Research and Development in Information Retrieval, pages 2220–2226, 2022. doi: 10.1145/3477495.3531833. [115] Nandan Thakur, Nils Reimers, Andreas Rücklé, Abhishek Srivastava, and Iryna Gurevych. BEIR: A heterogeneous benchmark for zero-shot evaluation of information retrieval models. In Thirty-fifth Conference on Neural Information Processing Systems Datasets and Benchmarks Track (Round 2), 2021. URL https://openreview.net/ forum?id=wCu6T5xFjeJ. [116] Tom Kwiatkowski, Jennimaria Palomaki, Olivia Redfield, Michael Collins, Ankur Parikh, Chris Alberti, Danielle Epstein, Illia Polosukhin, Jacob Devlin, Kenton Lee, Kristina Toutanova, Llion Jones, Matthew Kelcey, Ming-Wei Chang, Andrew M. Dai, Jakob Uszkoreit, Quoc Le, and Slav Petrov. Natural questions: A benchmark for question answering research. Transactions of the Association for Computational Linguistics, 7:452–466, 2019. doi: 10.1162/tacl_a_ 00276. [117] Elias Jääsaari, Ville Hyvönen, and Teemu Roos. LoRANN: Lowrank matrix factorization for approximate nearest neighbor search. In Proceedings of the 38th International Conference on Neural Information Processing Systems, NIPS ’24, Red Hook, NY, USA, 2024. Curran Associates Inc. ISBN 9798331314385. [118] Elizaveta Levina and Peter Bickel. Maximum likelihood estimation of intrinsic dimension. In Advances in Neural Information Processing Systems, volume 17, 2004. [119] Martin Aumüller and Matteo Ceccarello. The role of local dimensionality measures in benchmarking nearest neighbor search. Information Systems, 101:101807, 2021. ISSN 0306-4379. doi: 10.1016/j.is.2021.101807. URL https://www.sciencedirect.com/ science/article/pii/S0306437921000569.

Page 39 of 39

Related documents

• Unified and Efficient Approach for Multi-Vector Similarity Search#2768
• Graph-free Retrieval with Graph-Modified Embeddings#1343921
• RSLM: Training-Free Vector Quantization for Approximate Nearest Neighbor Search#799588
• Latent space configuration for neural network finetuning, fast ranking, and embedding dimensionality reduction#1413900
• Retrieval-Centric Deep Learning in Growing Nonparametric Neural Networks#1199132