feat(core): portable seeded rotations and a frozen fast-mode scale table

CMakeLists.txt

@@ -84,6 +84,7 @@ endfunction()
 
 set(RABITQ_CORE_SOURCES
     src/core/api.cpp
+    src/core/fast_scale_table.cpp
     src/utils/cpu_features.cpp
     src/utils/thread_count.cpp
     src/simd/dispatch.cpp

docs/development.md

@@ -48,7 +48,7 @@ factors estimate L2 distance or inner product.
 | `include/rabitqlib/clustering/` | RaBitQKMeans flat assignment, QGKMeans graph assignment, and shared Lloyd training |
 | `include/rabitqlib/index/{ivf,hnsw,symqg}/` | Index construction, persistence, and search |
 | `include/rabitqlib/{core,fastscan,index}/detail/` | Module-owned implementation headers; no public stability promise |
-| `src/core/` | Compiled stable core v1 API implementation |
+| `src/core/` | Compiled stable core v1 API implementation and the frozen `faster_config` scale table |
 | `src/simd/` | Core kernels, dispatch, and the private shared resolver |
 | `src/fastscan/` | FastScan kernels, batch estimators, and paired backend selection |
 | `src/index/` | Compiled IVF/HNSW/SymphonyQG operations and private search kernels |
@@ -470,6 +470,8 @@ contract. For a local implementation fix, identify its affected paths from the m
 below and check those paths; broaden the audit if the shared contract changes.
 
 - `core/encoder.hpp` and `core/detail/rabitq_impl.hpp`;
+- the frozen `faster_config` scales in `core/detail/fast_scale_table.hpp` and
+  `src/core/fast_scale_table.cpp`; never regenerate them or change their values;
 - `core/data_layout.hpp` and `core/pack_excode.hpp`;
 - SIMD pack, inner-product, and FastScan kernels;
 - query-side factors in `core/query.hpp` and `fastscan/query.hpp`;

docs/docs/index/qg.md

@@ -131,7 +131,9 @@ this for its centroid graph; the cache is not persisted.
 
 ### Random seeds
 
-The graph seed controls FHT/Kac rotation; MatrixRotator samples its own matrix.
+The graph seed controls the rotation. An FHT/Kac seed gives the same rotation on
+every platform; a MatrixRotator seed reproduces the matrix only with the same
+standard library and build, so rely on the saved index rather than the seed.
 The builder seed controls random initialization and fallback edges. PiPNN uses
 fixed internal seeds.
 

docs/docs/rabitq/quantizer.md

@@ -45,9 +45,16 @@ intervals; a scalar event sweep remains the reference and fallback when the SIMD
 search cannot certify a winner. This path preserves the reference tie and rounding rules.
 
 The fast path reuses a calibrated scaling factor instead of searching per vector.
-`faster_config` samples uniformly distributed unit vectors, finds their scaling
-factors with the default path, and averages them. Calibrate once during setup
-and reuse the result for the same padded dimension and total bit width. This
+The factor is the average of the default path's scaling factors over sampled
+uniformly distributed unit vectors. `faster_config` reads it from a frozen,
+precomputed table for padded dimensions that are multiples of 32 up to 65536 and
+total bit widths 2 through 9, so the same arguments give the same factor on every
+platform, standard library and CPU, at no setup cost. For other dimensions,
+`faster_config` computes the factor from `std::normal_distribution` samples, which
+differ between standard libraries (MSVC and libc++ results differ from libstdc++ by
+several percent) and, in the last bits, between FMA and non-FMA builds and between
+math libraries; calibrate once and reuse the
+result for the same padded dimension and total bit width. This
 reduces encoding work at a possible accuracy cost; evaluate it on your data.
 
 ## Data format

docs/docs/rabitq/rotator.md

@@ -56,8 +56,11 @@ padded dimension explicitly if subsequent packed-code operations require a multi
 of 32. Stored space and per-vector work are proportional to `dim * padded_dim`
 (quadratic when both dimensions are equal).
 Constructing the random matrix also performs QR factorization. The `seed` argument
-of `choose_rotator` controls FHT/Kac signs only; MatrixRotator samples its own matrix.
-Save and restore the actual rotator state rather than relying on that seed to
+of `choose_rotator` seeds both rotators. FHT/Kac signs come from a fixed mapping of
+`std::mt19937` output, so a seed gives the same rotation on every platform and
+standard library. A seeded MatrixRotator is reproducible only with the same standard
+library and build, because `std::normal_distribution` and QR rounding differ between
+them. Save and restore the actual rotator state rather than relying on a seed to
 recreate a matrix rotation.
 
 ## FFHT + Kac’s walk

include/rabitqlib/core/detail/fast_scale_table.hpp

@@ -0,0 +1,23 @@
+#pragma once
+
+#include <cstddef>
+#include <optional>
+
+namespace rabitqlib::quant::rabitq_impl::ex_bits {
+
+// Frozen fast-mode scales, so faster_config gives the same t_const on every platform,
+// standard library, compiler and CPU. The table covers padded_dim in multiples of
+// kFastScaleDimStep up to kFastScaleMaxDim and ex_bits in [1, kFastScaleMaxExBits].
+// Values are const_scaling_factor(padded_dim, ex_bits) as computed once by
+// scripts/generate_fast_scale_table.cpp; see src/core/fast_scale_table.cpp. The
+// table is frozen data: never regenerate it or change its values.
+inline constexpr size_t kFastScaleDimStep = 32;
+inline constexpr size_t kFastScaleMaxDim = 65536;
+inline constexpr size_t kFastScaleMaxExBits = 8;
+
+// Returns the frozen scale, or no value when (padded_dim, ex_bits) is outside the table.
+std::optional<double> frozen_const_scaling_factor(
+    size_t padded_dim, size_t ex_bits
+) noexcept;
+
+}  // namespace rabitqlib::quant::rabitq_impl::ex_bits

include/rabitqlib/core/detail/rabitq_impl.hpp

@@ -547,7 +547,9 @@ inline T quantize_ex(const T* o_abs, TP* code, size_t dim, size_t ex_bits) {
     );
 }
 
-// For given dim and ex_bits, use random vectors to get the const rescale factor
+// For given dim and ex_bits, use random vectors to get the const rescale factor.
+// faster_config reads these values from the frozen table (fast_scale_table.hpp);
+// this computation generates that table and serves dimensions outside it.
 inline double const_scaling_factor(size_t dim, size_t ex_bits) {
     constexpr Eigen::Index kConstNum = 100;
 

include/rabitqlib/core/encoder.hpp

@@ -2,10 +2,12 @@
 
 #include <cstddef>
 #include <cstdint>
+#include <optional>
 #include <vector>
 
 #include "rabitqlib/core/data_layout.hpp"
 #include "rabitqlib/core/detail/eigen.hpp"
+#include "rabitqlib/core/detail/fast_scale_table.hpp"
 #include "rabitqlib/core/detail/rabitq_impl.hpp"
 #include "rabitqlib/core/metric.hpp"
 #include "rabitqlib/core/types.hpp"  // IWYU pragma: export
@@ -14,10 +16,19 @@ namespace rabitqlib::quant {
 
 // Prepare a constant rescaling factor for total_bits > 1. The one-bit path
 // needs no extra-code scaling and leaves the default configuration unchanged.
+// Padded dimensions in the frozen table give the same factor on every platform.
+// Others fall back to computing it from std::normal_distribution samples, which differ
+// between standard libraries (by several percent) and between FMA and non-FMA builds
+// or math libraries (in the last bits).
 inline RabitqConfig faster_config(size_t dim, size_t total_bits) {
     RabitqConfig config;
     if (total_bits > 1) {
-        config.t_const = rabitq_impl::ex_bits::const_scaling_factor(dim, total_bits - 1);
+        const size_t ex_bits = total_bits - 1;
+        const std::optional<double> frozen =
+            rabitq_impl::ex_bits::frozen_const_scaling_factor(dim, ex_bits);
+        config.t_const = frozen.has_value()
+                             ? *frozen
+                             : rabitq_impl::ex_bits::const_scaling_factor(dim, ex_bits);
     }
     return config;
 }

include/rabitqlib/core/rotator.hpp

@@ -89,9 +89,13 @@ class MatrixRotator : public Rotator<T> {
     RowMajorMatrix<T> rand_mat_;
 
    public:
-    explicit MatrixRotator(size_t dim, size_t padded_dim)
+    // The seed makes the matrix reproducible only with the same standard library
+    // (std::normal_distribution) and build (QR rounding); persist it with save().
+    explicit MatrixRotator(
+        size_t dim, size_t padded_dim, uint32_t seed = std::random_device{}()
+    )
         : MatrixRotator(dim, padded_dim, AllocateOnly{}) {
-        RowMajorMatrix<T> rand = random_gaussian_matrix<T>(padded_dim, padded_dim);
+        RowMajorMatrix<T> rand = random_gaussian_matrix<T>(padded_dim, padded_dim, seed);
         Eigen::HouseholderQR<RowMajorMatrix<T>> qr(rand);
         RowMajorMatrix<T> q_inv =
             qr.householderQ()
@@ -173,10 +177,13 @@ class FhtKacRotator : public Rotator<float> {
         size_t dim, size_t padded_dim, uint32_t seed = std::random_device{}()
     )
         : FhtKacRotator(dim, padded_dim, AllocateOnly{}) {
+        // Frozen, portable mapping: the top byte of each mt19937 output. Its output
+        // sequence is fixed by the standard, unlike std::uniform_int_distribution;
+        // GCC 11+ libstdc++'s uniform_int_distribution<int>(0, 255) yields these same
+        // bytes; older libstdc++ occasionally drew different ones.
         std::mt19937 gen(seed);
-        std::uniform_int_distribution<int> dist(0, 255);
         for (auto& i : flip_) {
-            i = static_cast<uint8_t>(dist(gen));
+            i = static_cast<uint8_t>(gen() >> 24);
         }
     }
 
@@ -251,8 +258,9 @@ inline std::unique_ptr<Rotator<float>> make_rotator_for_loading(
 }  // namespace rotator_impl
 
 // Return an owning pointer; wrap it in RAII and use size() for the output dimension.
-// FHT/Kac supports float and uses seed. Matrix rotation supports float/double;
-// it samples its own matrix and does not use this factory's seed argument.
+// FHT/Kac supports float; its seed gives the same rotation on every platform.
+// Matrix rotation supports float/double; its seed reproduces the matrix only with
+// the same standard library and build, so persist matrix rotators with save().
 template <typename T>
 Rotator<T>* choose_rotator(
     size_t dim,
@@ -271,7 +279,7 @@ Rotator<T>* choose_rotator(
     }
 
     if (type == RotatorType::MatrixRotator) {
-        return ::new rotator_impl::MatrixRotator<T>(dim, padded_dim);
+        return ::new rotator_impl::MatrixRotator<T>(dim, padded_dim, seed);
     }
 
     throw std::invalid_argument("Invalid rotator type in choose_rotator()");

scripts/generate_fast_scale_table.cpp

@@ -0,0 +1,114 @@
+// Historical record of how src/core/fast_scale_table.cpp, the frozen fast-mode
+// scales read by quant::faster_config, was generated. The table is frozen data and
+// is never regenerated: rerunning this after a change to the computation would
+// silently change every fast-mode code. A different table would be a new, separately
+// named table. The original run, from the repository root:
+//
+//   g++ -std=c++17 -O3 -ffp-contract=off -pthread -DEIGEN_MAX_ALIGN_BYTES=64 \
+//       -DEIGEN_MAX_STATIC_ALIGN_BYTES=64 -Iinclude \
+//       scripts/generate_fast_scale_table.cpp -o /tmp/generate_fast_scale_table
+//   /tmp/generate_fast_scale_table > src/core/fast_scale_table.cpp
+//   clang-format-23 -i src/core/fast_scale_table.cpp
+//
+// The frozen values are the ones portable x86-64 Linux builds computed before the
+// table existed: libstdc++'s std::normal_distribution, the scalar event sweep (the
+// double-precision search never uses the SIMD kernel) and no fused multiply-add.
+// Other standard libraries draw different samples and change the values by several
+// percent; FMA contraction or AVX vectorization (a different Eigen reduction order)
+// changes the last bits. The guards below refuse to build in those configurations.
+
+#include <algorithm>
+#include <atomic>
+#include <cstddef>
+#include <cstdio>
+#include <thread>
+#include <vector>
+
+#include "rabitqlib/core/detail/fast_scale_table.hpp"
+#include "rabitqlib/core/detail/rabitq_impl.hpp"
+
+#if !defined(__GLIBCXX__)
+#error "Generate the table with libstdc++; its std::normal_distribution defines the values"
+#endif
+#if !defined(__x86_64__) || defined(__AVX__) || defined(__FMA__)
+#error "Generate the table on baseline x86-64 (SSE2, no -m/-march flags), -ffp-contract=off"
+#endif
+
+namespace ex_bits = rabitqlib::quant::rabitq_impl::ex_bits;
+
+int main() {
+    constexpr size_t kRows = ex_bits::kFastScaleMaxDim / ex_bits::kFastScaleDimStep;
+    constexpr size_t kCols = ex_bits::kFastScaleMaxExBits;
+    std::vector<double> table(kRows * kCols);
+
+    // Largest dimensions first, so the slowest jobs do not finish last.
+    std::atomic<size_t> next{0};
+    auto worker = [&] {
+        for (size_t job = next++; job < table.size(); job = next++) {
+            const size_t index = table.size() - 1 - job;
+            const size_t padded_dim = (index / kCols + 1) * ex_bits::kFastScaleDimStep;
+            const size_t bits = index % kCols + 1;
+            table[index] = ex_bits::const_scaling_factor(padded_dim, bits);
+        }
+    };
+    std::vector<std::thread> threads(std::max(1U, std::thread::hardware_concurrency()));
+    for (auto& thread : threads) {
+        thread = std::thread(worker);
+    }
+    for (auto& thread : threads) {
+        thread.join();
+    }
+
+    std::printf(
+        "// Generated once by scripts/generate_fast_scale_table.cpp. Frozen data: never\n"
+        "// regenerate it or change its values; a different table would be a new, "
+        "separately\n"
+        "// named table.\n"
+        "// Each value is const_scaling_factor(padded_dim, ex_bits) from libstdc++ on\n"
+        "// baseline x86-64 (SSE2, no FMA), written as an exact hexadecimal "
+        "floating-point\n"
+        "// literal.\n"
+        "\n"
+        "#include \"rabitqlib/core/detail/fast_scale_table.hpp\"\n"
+        "\n"
+        "#include <array>\n"
+        "#include <cstddef>\n"
+        "#include <optional>\n"
+        "\n"
+        "namespace rabitqlib::quant::rabitq_impl::ex_bits {\n"
+        "namespace {\n"
+        "\n"
+        "constexpr size_t kFastScaleRows = kFastScaleMaxDim / kFastScaleDimStep;\n"
+        "\n"
+        "// Entry (row, col) holds padded_dim = (row + 1) * kFastScaleDimStep and\n"
+        "// ex_bits = col + 1, at index row * kFastScaleMaxExBits + col.\n"
+        "constexpr std::array<double, kFastScaleRows * kFastScaleMaxExBits> kFastScales = "
+        "{\n"
+    );
+    for (size_t row = 0; row < kRows; ++row) {
+        std::printf("    ");
+        for (size_t col = 0; col < kCols; ++col) {
+            std::printf("%a,%s", table[row * kCols + col], col + 1 < kCols ? " " : "\n");
+        }
+    }
+    std::printf(
+        "};\n"
+        "\n"
+        "}  // namespace\n"
+        "\n"
+        "std::optional<double> frozen_const_scaling_factor(\n"
+        "    size_t padded_dim, size_t ex_bits\n"
+        ") noexcept {\n"
+        "    if (padded_dim == 0 || padded_dim %% kFastScaleDimStep != 0 ||\n"
+        "        padded_dim > kFastScaleMaxDim || ex_bits == 0 ||\n"
+        "        ex_bits > kFastScaleMaxExBits) {\n"
+        "        return std::nullopt;\n"
+        "    }\n"
+        "    const size_t row = padded_dim / kFastScaleDimStep - 1;\n"
+        "    return kFastScales[(row * kFastScaleMaxExBits) + ex_bits - 1];\n"
+        "}\n"
+        "\n"
+        "}  // namespace rabitqlib::quant::rabitq_impl::ex_bits\n"
+    );
+    return 0;
+}