Oregami
Repositories/oxedyne/fe2o3

oxedyne/fe2o3/fe2o3_data/src/hll/sketch.rs

6.7 KiB, 1 run

created by r1870400018:11204, which is this file's identity for as long as the history lasts, whatever it is later renamed to

download · who wrote it · its history

1//! The HyperLogLog cardinality sketch.
2//!
3//! A [`HyperLogLog`] sketch uses $m = 2^p$ single-byte registers to estimate
4//! the number of distinct 64-bit hashes it has observed. The sketch is fixed
5//! size regardless of the true cardinality and merges with any other sketch
6//! of the same precision by register-wise maximum.
7//!
8//! This crate does not bundle a hash function. Callers hash their inputs to
9//! `u64` using whichever algorithm suits them (SeaHash, SipHash, SHA3
10//! truncated, etc.) and call [`HyperLogLog::add_hash`]. Keeping the hash
11//! choice external preserves the primitive's independence from any particular
12//! authentication or cryptographic scheme.
13
14use oxedyne_fe2o3_core::prelude::*;
15
16
17/// The minimum precision parameter. With `p = 4` the sketch has 16 registers.
18pub const P_MIN: u8 = 4;
19
20/// The maximum precision parameter. With `p = 18` the sketch has 262144
21/// registers. Above this the leading-zero range left in the 64-bit hash
22/// collapses to fewer than 64 bits of entropy on the register side, which is
23/// both wasteful and risks undercounting.
24pub const P_MAX: u8 = 18;
25
26/// The precision used by the distributed Ozone layer. 16384 registers, each
27/// one byte -- a 16 KiB sketch, matching #raw("sec_ozone.typ") §"Network Size
28/// Estimation: HyperLogLog".
29pub const P_DEFAULT: u8 = 14;
30
31
32/// A HyperLogLog cardinality sketch.
33///
34/// Holds $m = 2^p$ single-byte registers. Each register stores
35/// $max(rho(h) : h mod m = j)$ where $rho(h)$ is the 1-based position of the
36/// first set bit in the hash suffix after the register-selecting prefix, and
37/// the max is taken over every hash observed so far that maps to register $j$.
38#[derive(Clone, Debug)]
39pub struct HyperLogLog {
40 /// The precision parameter.
41 p: u8,
42 /// The $m = 2^p$ register bytes.
43 registers: Vec<u8>,
44}
45
46impl HyperLogLog {
47 /// Builds an empty sketch with precision `p`.
48 ///
49 /// Validates `p ∈ [P_MIN, P_MAX]`. Allocates `2^p` zeroed bytes.
50 pub fn new(p: u8) -> Outcome<Self> {
51 if p < P_MIN || p > P_MAX {
52 return Err(err!(
53 "HyperLogLog precision p = {} out of range [{}, {}].",
54 p, P_MIN, P_MAX;
55 Invalid, Input));
56 }
57 let m = 1usize << p;
58 Ok(Self {
59 p,
60 registers: vec![0u8; m],
61 })
62 }
63
64 /// Constructs a sketch from raw register bytes.
65 ///
66 /// Validates `p ∈ [P_MIN, P_MAX]` and that `bytes.len() == 2^p`. The
67 /// bytes are copied into the sketch unchanged -- register values are not
68 /// clamped; a malformed peer could feed a register above `64 - p + 1` and
69 /// inflate the estimate. Callers exchanging sketches over the wire should
70 /// apply their own rate limiting and reputation accounting before merging.
71 pub fn from_bytes(p: u8, bytes: &[u8]) -> Outcome<Self> {
72 if p < P_MIN || p > P_MAX {
73 return Err(err!(
74 "HyperLogLog precision p = {} out of range [{}, {}].",
75 p, P_MIN, P_MAX;
76 Invalid, Input));
77 }
78 let m = 1usize << p;
79 if bytes.len() != m {
80 return Err(err!(
81 "HyperLogLog for p = {} needs {} register bytes, got {}.",
82 p, m, bytes.len();
83 Invalid, Input, Size));
84 }
85 Ok(Self {
86 p,
87 registers: bytes.to_vec(),
88 })
89 }
90
91 /// Returns the precision parameter.
92 pub fn precision(&self) -> u8 {
93 self.p
94 }
95
96 /// Returns the number of registers, `m = 2^p`.
97 pub fn m(&self) -> usize {
98 self.registers.len()
99 }
100
101 /// Borrows the register bytes for serialisation.
102 pub fn as_bytes(&self) -> &[u8] {
103 &self.registers
104 }
105
106 /// Incorporates a 64-bit hash into the sketch.
107 ///
108 /// The top `p` bits select the register; the remaining `64 - p` bits are
109 /// scanned for the position of the first set bit, which if greater than
110 /// the current register value replaces it.
111 pub fn add_hash(&mut self, hash: u64) {
112 let p = self.p as u32;
113 // Register index from the top p bits.
114 let idx = (hash >> (64 - p)) as usize;
115 // Suffix = low (64 - p) bits shifted into the top of a u64 so that
116 // the suffix's MSB sits at bit 63 and its LSB at bit p. The bottom p
117 // bits are zero by construction.
118 let suffix = hash << p;
119 let rho = if suffix == 0 {
120 // No set bit in the suffix; the conventional cap.
121 (64 - p + 1) as u8
122 } else {
123 // leading_zeros counts zero bits from bit 63 downward. Since the
124 // suffix has at least one set bit in positions p..64, the count
125 // sits in [0, 63 - p]; +1 for the 1-based rho convention leaves
126 // rho in [1, 64 - p].
127 suffix.leading_zeros() as u8 + 1
128 };
129 if rho > self.registers[idx] {
130 self.registers[idx] = rho;
131 }
132 }
133
134 /// Merges `other` into `self` by register-wise maximum.
135 ///
136 /// Returns an error if the two sketches have different precision.
137 pub fn merge(&mut self, other: &Self) -> Outcome<()> {
138 if self.p != other.p {
139 return Err(err!(
140 "Cannot merge HyperLogLog sketches with different precision \
141 (self p = {}, other p = {}).",
142 self.p, other.p;
143 Invalid, Input, Mismatch));
144 }
145 for (dst, src) in self.registers.iter_mut().zip(other.registers.iter()) {
146 if *src > *dst {
147 *dst = *src;
148 }
149 }
150 Ok(())
151 }
152
153 /// Returns the current cardinality estimate.
154 ///
155 /// The formula is:
156 ///
157 /// - Raw: $E = alpha_m dot m^2 / sum_j 2^(-M_j)$.
158 /// - Linear counting for small cardinalities: if $E <= 5/2 dot m$ and at
159 /// least one register is zero, return $m dot ln(m / z)$ where $z$ is
160 /// the number of zero registers. This is the well-known HLL small-range
161 /// correction.
162 ///
163 /// At `p = 14` the expected standard error is approximately
164 /// $1.04 / sqrt(m) ≈ 0.008$, i.e. around 0.8%. The spec's 2% target at
165 /// $10^6$ peers is comfortably within that.
166 pub fn estimate(&self) -> f64 {
167 let m = self.registers.len() as f64;
168 let alpha = alpha_m(self.registers.len());
169
170 let mut sum = 0.0f64;
171 let mut zeros = 0usize;
172 for &r in &self.registers {
173 if r == 0 {
174 zeros += 1;
175 }
176 // 2^(-r) computed as ldexp(1, -r). For r up to ~65 this is well
177 // within f64 precision.
178 sum += (-(r as f64)).exp2();
179 }
180 let raw = alpha * m * m / sum;
181
182 // Linear counting correction for small cardinalities.
183 let small_threshold = 2.5f64 * m;
184 if raw <= small_threshold && zeros > 0 {
185 return m * (m / zeros as f64).ln();
186 }
187 raw
188 }
189
190 /// Convenience wrapper around [`HyperLogLog::estimate`] that rounds to
191 /// the nearest non-negative integer.
192 pub fn estimate_rounded(&self) -> u64 {
193 let e = self.estimate();
194 if e < 0.0 {
195 0
196 } else {
197 e.round() as u64
198 }
199 }
200
201 /// Resets every register to zero without reallocating.
202 pub fn clear(&mut self) {
203 for r in self.registers.iter_mut() {
204 *r = 0;
205 }
206 }
207}
208
209/// The $alpha_m$ bias-correction constant from the original HyperLogLog
210/// paper. For $m in {16, 32, 64}$ the constants are tabulated; for larger
211/// $m$ the formula $alpha_m = 0.7213 / (1 + 1.079 / m)$ is used.
212fn alpha_m(m: usize) -> f64 {
213 match m {
214 16 => 0.673,
215 32 => 0.697,
216 64 => 0.709,
217 _ => 0.7213 / (1.0 + 1.079 / m as f64),
218 }
219}