echo.rs
8.2 kB · rust · 220 lines
1use crate::Fault;2use mrlycore::json;3use mrlynum::design;4use wasm_bindgen::prelude::*;56const SAMPLES: usize = 4096;7const ELEMENT_CAP: u128 = 1 << 18;8const VALUE_CAP: u128 = 1 << 27;9const ECHO_CAP: u128 = 1 << 24;10const FLOOR: usize = 101;11const BAND: (f64, f64) = (4.0, 60.0);12const THRESHOLD: f64 = 8.0;13const TOP: usize = 10;14const LOW: usize = 3;1516/// The design Mobius meter read at one base, digit set and depth: the meter drawn against log x, its density echo and residual, and the spectrum they carry.17#[wasm_bindgen(getter_with_clone)]18pub struct Echo {19 /// The log of x at every sample of the meter.20 pub logx: Vec<f32>,21 /// The meter M_F(x) over x to the alpha over two at every sample.22 pub meter: Vec<f32>,23 /// The density echo at every sample, empty when the depth is past the sieve cap.24 pub echo: Vec<f32>,25 /// The residual, the meter less its echo, empty when the depth is past the sieve cap.26 pub rest: Vec<f32>,27 /// The ordinate gamma of every spectral bin.28 pub gamma: Vec<f32>,29 /// The power over its local median floor at every spectral bin.30 pub score: Vec<f32>,31 /// The reading beside the curves, as JSON.32 pub read: String,33}3435fn digits_of(base: u32, mask: u32) -> Result<(u64, Vec<u64>), Fault> {36 if !(2..=10).contains(&base) {37 return Err(Fault::new(format!(38 "the base {base} is not between two and ten."39 )));40 }41 if mask >> base != 0 {42 return Err(Fault::new(format!(43 "the digit set names a digit at or above the base {base}."44 )));45 }46 let digits = design::digits_of(mask, u64::from(base));47 if digits.len() < 2 {48 return Err(Fault::new(49 "the design needs at least two digits.".to_string(),50 ));51 }52 if digits.iter().all(|&d| d == 0) {53 return Err(Fault::new(54 "the design needs a digit above zero.".to_string(),55 ));56 }57 Ok((u64::from(base), digits))58}5960fn depths(base: u64, digits: &[u64]) -> (usize, usize) {61 let mut deepest = 0;62 let mut sieved = 0;63 for depth in 1..64 {64 let span = match u128::from(base).checked_pow(depth as u32) {65 Some(span) => span,66 None => break,67 };68 if span > VALUE_CAP || design::size(digits, depth) > ELEMENT_CAP {69 break;70 }71 deepest = depth;72 if span <= ECHO_CAP {73 sieved = depth;74 }75 }76 (deepest, sieved)77}7879fn matched(peak: f64, list: &[f64], bin: f64) -> bool {80 design::nearest(peak, list) <= bin81}8283/// Returns the digit set the mask names, its exponent alpha and the depths the caps allow, as JSON.84///85/// The deepest depth is the last one whose element count stays inside `2^18` and whose span86/// `q^L` stays inside `2^27`; the sieved depth is the last one whose span stays inside `2^24`,87/// past which the density echo needs a Mobius sieve this page will not run.88#[wasm_bindgen]89pub fn echo_caps(base: u32, mask: u32) -> Result<String, Fault> {90 let (base, digits) = digits_of(base, mask)?;91 let (deepest, sieved) = depths(base, &digits);92 if deepest < LOW {93 return Err(Fault::new(format!(94 "base {base} with {} digits reaches only depth {deepest}.",95 digits.len()96 )));97 }98 Ok(json!({99 "base": base,100 "digits": digits,101 "k": digits.len(),102 "alpha": (digits.len() as f64).ln() / (base as f64).ln(),103 "deepest": deepest,104 "sieved": sieved,105 "least": LOW,106 "samples": SAMPLES,107 "caps": { "elements": ELEMENT_CAP as f64, "span": VALUE_CAP as f64, "echo": ECHO_CAP as f64 },108 })109 .to_string())110}111112/// Reads the design Mobius meter at the base, digit set and depth: the meter resampled uniformly in log x and scaled by x to the alpha over two, its density echo and residual when the span is inside the sieve cap, and the power spectrum of the meter or of the residual against its local median floor.113///114/// The elements of `S_F` are enumerated in order and their Mobius values taken by trial115/// division; the echo is `sum of mu(n) A_F(n)/n` over every whole number up to x, which needs116/// a Mobius sieve to the span and is refused past `2^24`. The spectrum is mean-removed,117/// Hann-windowed and read as `gamma = 2 pi j` over the log range, and a peak counts as landing118/// on a list when it sits within one bin of one of its entries.119#[wasm_bindgen]120pub fn echo_read(base: u32, mask: u32, depth: usize, subtract: bool) -> Result<Echo, Fault> {121 let (base, digits) = digits_of(base, mask)?;122 let (deepest, sieved) = depths(base, &digits);123 if !(LOW..=deepest).contains(&depth) {124 return Err(Fault::new(format!(125 "the depth must be between {LOW} and {deepest} at this base and digit set."126 )));127 }128 let values = design::elements(base, &digits, depth);129 let mu = design::mobius_of(&values);130 let running = design::meter(&mu);131 let k = digits.len() as f64;132 let alpha = k.ln() / (base as f64).ln();133 let exponent = alpha / 2.0;134 let logx = design::log_grid(&values, SAMPLES);135 let meter = design::resample(&values, &running, exponent, &logx);136 let sieve = depth <= sieved;137 let echo = if sieve {138 design::echo_series(&values, &logx, exponent)139 } else {140 Vec::new()141 };142 let rest: Vec<f64> = meter143 .iter()144 .zip(echo.iter())145 .map(|(whole, part)| whole - part)146 .collect();147 let taken = if subtract && sieve { &rest } else { &meter };148 let (gamma, power) = design::spectrum(&logx, taken);149 let score = design::score(&power, FLOOR);150 let found = design::peaks(&gamma, &score, BAND, THRESHOLD);151 let bin = gamma.get(1).copied().unwrap_or(0.0);152 let zeros: Vec<f64> = design::ZETA_ORDINATES153 .iter()154 .copied()155 .filter(|g| *g > BAND.0 && *g < BAND.1)156 .collect();157 let lattice: Vec<f64> = design::pole_lattice(base, BAND.1)158 .into_iter()159 .filter(|g| *g > BAND.0)160 .collect();161 let head: Vec<usize> = found.iter().copied().take(TOP).collect();162 let rows: Vec<mrlycore::Json> = head163 .iter()164 .map(|&at| {165 json!({166 "gamma": gamma[at],167 "score": score[at],168 "zeta": design::nearest(gamma[at], &zeros),169 "lattice": design::nearest(gamma[at], &lattice),170 })171 })172 .collect();173 let last = running.last().copied().unwrap_or(0);174 let peak = running.iter().map(|v| v.abs()).max().unwrap_or(0);175 let scale = design::upper_rms(&meter);176 let width = BAND.1 - BAND.0;177 let chance = |list: &[f64]| (2.0 * bin * list.len() as f64 / width).min(1.0);178 let read = json!({179 "base": base,180 "digits": digits,181 "k": digits.len(),182 "depth": depth,183 "deepest": deepest,184 "sieved": sieved,185 "sieve": sieve,186 "alpha": alpha,187 "count": values.len(),188 "last": last,189 "peak": peak,190 "theta": if values.len() > 1 { (peak.max(1) as f64).ln() / (values.len() as f64).ln() } else { 0.0 },191 "span": logx.last().copied().unwrap_or(0.0) - logx.first().copied().unwrap_or(0.0),192 "bin": bin,193 "samples": SAMPLES,194 "band": [BAND.0, BAND.1],195 "threshold": THRESHOLD,196 "zeros": zeros,197 "lattice": lattice,198 "peaks": rows,199 "found": found.len(),200 "hits": head.iter().filter(|&&at| matched(gamma[at], &zeros, bin)).count(),201 "lines": head.iter().filter(|&&at| matched(gamma[at], &lattice, bin)).count(),202 "share": if sieve && scale > 0.0 { json!(design::upper_rms(&echo) / scale) } else { json!(null) },203 "residual": if sieve && scale > 0.0 { json!(design::upper_rms(&rest) / scale) } else { json!(null) },204 "chance": chance(&zeros),205 "chanceLines": chance(&lattice),206 "rate": -alpha.min(1.0 - alpha) / 2.0,207 "caps": { "elements": ELEMENT_CAP as f64, "span": VALUE_CAP as f64, "echo": ECHO_CAP as f64 },208 })209 .to_string();210 let thin = |series: &[f64]| series.iter().map(|&v| v as f32).collect::<Vec<f32>>();211 Ok(Echo {212 logx: thin(&logx),213 meter: thin(&meter),214 echo: thin(&echo),215 rest: thin(&rest),216 gamma: thin(&gamma),217 score: thin(&score),218 read,219 })220}