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}