nx_sketch_geomean.nx source
↩ module page · 147 lines · 4862 B
1// sketch_geomean.nx -- streaming geometric mean via log-space accumulation.
2//
3// Geometric mean = (Π x_i)^(1/n) = 2^((Σ log_2 x_i) / n)
4//
5// Naive product overflows immediately (n=20, x=10 -> product ~ 10^20 >> i64).
6// We accumulate in log_2 domain instead: log space is additive and bounded.
7//
8// USE CASES (where geometric mean beats arithmetic mean):
9// - growth-rate averaging (CAGR / compound returns)
10// - ratios (speedup factors, multiplicative gain)
11// - normalized scores (e.g., averaging across scales)
12// - magnitude data (sizes, populations, intensities)
13//
14// INTEGER LOG-BASE-2 (no f64):
15// floor(log_2(x)) via bit_length(x) - 1
16// For finer fractional bits, mantissa interpolation could be added
17// (queued for v2).
18//
19// SUM_LOG_2 fixed-point: track in PPM scale. Each insert contributes
20// (bitlen - 1) * 1_000_000. Geometric mean recovered via 2^(sum/n_ppm).
21//
22// For 2^x where x is in PPM (fractional log_2): floor(x / 1_000_000) +
23// linear interpolation on the fractional part. Returns approximate
24// geometric mean.
25//
26// LOSSLESS-LANGUAGE DISCIPLINE: NX_ENV_REL_STDDEV with param_a
27// reflecting quantization (use 50_000_000 ppb = 5% as conservative
28// bound for floor-log estimator). ReferenceImpl tier.
29
30// nx_safety_envelope:
31// intended_use: AUTO_APPLIED -- primitive-specific tuning queued
32// sil_target: SIL1
33// evidence: [bulk_applied_2026-05-16, see-file-comment-for-detail]
34// verdict: NOT_YET_EVALUATED
35
36import "nx_syscalls.nx"
37import "nx_sketch_types.nx"
38
39struct Geomean {
40 sum_log_ppm: i64, // Σ log_2(x_i) * 1_000_000
41 n: i64,
42 n_positive: i64, // x_i must be > 0 for log; track invalid
43}
44
45// === construction =================================================
46
47func nx_gm_alloc() -> *Geomean {
48 let raw: *u8 = sys_mmap(32)
49 let g: *Geomean = raw as *Geomean
50 g.sum_log_ppm = 0
51 g.n = 0
52 g.n_positive = 0
53 return g
54}
55
56// === log_2 floor ==================================================
57
58func nx_gm_log2_floor(x: i64) -> i64 {
59 if x <= 1 { return 0 }
60 var n: i64 = 0
61 var t: i64 = x
62 while t > 1 {
63 t = t >> 1
64 n = n + 1
65 }
66 return n
67}
68
69// === add ==========================================================
70//
71// Returns 0 on success; -1 if value <= 0 (geometric mean undefined for
72// non-positive values).
73
74func nx_gm_add(g: *Geomean, value: i64) -> i64 {
75 if value <= 0 { return -1 }
76 g.sum_log_ppm = g.sum_log_ppm + nx_gm_log2_floor(value) * 1000000
77 g.n = g.n + 1
78 g.n_positive = g.n_positive + 1
79 return 0
80}
81
82// === 2^x reconstruction ===========================================
83//
84// Input: x_ppm (fractional log_2 in PPM, x = k + f where f in [0, 1)).
85// Output: 2^(k + f) = 2^k * 2^f.
86// 2^f for f in [0, 1) approximated via linear interpolation:
87// 2^f ≈ 1 + f (Bernoulli's inequality, conservative; underestimates by ~6%)
88// More accurate: 2^f ≈ 1 + 0.69 * f (using ln(2) factor)
89// We use the linear approximation for simplicity; v2 can use a 16-entry
90// lookup table.
91
92func nx_gm_pow2_floor(x_ppm: i64) -> i64 {
93 if x_ppm < 0 { return 0 }
94 let k: i64 = x_ppm / 1000000
95 let f_ppm: i64 = x_ppm - k * 1000000 // fractional, in PPM
96 // Cap k to prevent overflow (2^63 = 9.2e18).
97 if k > 62 { return 0x7FFFFFFFFFFFFFFF }
98 var base: i64 = 1
99 if k > 0 {
100 base = 1 << k
101 }
102 // 2^f ≈ (1_000_000 + f_ppm) / 1_000_000
103 return (base * (1000000 + f_ppm)) / 1000000
104}
105
106// === geometric mean ==============================================
107//
108// gmean = 2^(sum_log_ppm / n)
109
110func nx_gm_value(g: *Geomean) -> i64 {
111 if g.n_positive == 0 { return 0 }
112 let log_mean_ppm: i64 = g.sum_log_ppm / g.n_positive
113 return nx_gm_pow2_floor(log_mean_ppm)
114}
115
116// Equivalent log-domain query (useful for callers who want to combine).
117func nx_gm_log_mean_ppm(g: *Geomean) -> i64 {
118 if g.n_positive == 0 { return 0 }
119 return g.sum_log_ppm / g.n_positive
120}
121
122// === typed envelope ==============================================
123
124func nx_gm_query(g: *Geomean) -> *ApproxI64 {
125 let v: i64 = nx_gm_value(g)
126 return nx_approx_new(v, NX_ENV_REL_STDDEV, 50000000, // ~5% conservative
127 682700000,
128 NX_MATURITY_REFERENCE_IMPL,
129 NX_ADV_HONEST)
130}
131
132func nx_gm_count(g: *Geomean) -> i64 { return g.n }
133func nx_gm_count_positive(g: *Geomean) -> i64 { return g.n_positive }
134
135func nx_gm_memory_bytes(g: *Geomean) -> i64 {
136 return 32
137}
138
139// === merge =======================================================
140
141func nx_gm_merge(a: *Geomean, b: *Geomean) -> *Geomean {
142 let out: *Geomean = nx_gm_alloc()
143 out.sum_log_ppm = a.sum_log_ppm + b.sum_log_ppm
144 out.n = a.n + b.n
145 out.n_positive = a.n_positive + b.n_positive
146 return out
147}