novelty.rs

1use crate::Fault;2use mrlyrs::num::factor::totients;3use mrlyrs::num::zeta::{4    novelty_main, novelty_wave, sharp_novelty, smoothed_novelty, Complex, Line,5};6use wasm_bindgen::prelude::*;78const LOW: f64 = 8.0;9const HIGHEST: f64 = 21.0;10const PER_OCTAVE: usize = 32;11const ZEROS: usize = 138;1213/// The novelty meter: the totients sieved once, the smoothed and the sharp novelty error read on a log grid of y, and the first zeros with the coefficients of their waves.14#[wasm_bindgen]15pub struct Novelty {16    heights: Vec<f64>,17    smooth: Vec<f64>,18    rough: Vec<f64>,19    gammas: Vec<f64>,20    coef: Vec<Complex>,21    sieve: usize,22}2324#[wasm_bindgen]25impl Novelty {26    /// Builds the meter on y = 2^-j for j from eight to the height at the given samples per octave, the totients sieved to two to the height plus one, with the given count of zeros, at most 138.27    #[wasm_bindgen(constructor)]28    pub fn new(high: f64, per_octave: usize, zeros: usize) -> Result<Novelty, Fault> {29        if !(LOW + 1.0..=HIGHEST).contains(&high) {30            return Err(Fault::new(format!(31                "the height runs from {} to {HIGHEST}.",32                LOW + 1.033            )));34        }35        if !(1..=PER_OCTAVE).contains(&per_octave) {36            return Err(Fault::new(format!(37                "an octave takes 1 to {PER_OCTAVE} samples."38            )));39        }40        if zeros > ZEROS {41            return Err(Fault::new(format!("the meter holds {ZEROS} zeros.")));42        }43        let sieve = 2.0f64.powf(high + 1.0).ceil() as usize;44        let phi = totients(sieve);45        let mut prefix = vec![0u64; sieve + 1];46        for n in 1..=sieve {47            prefix[n] = prefix[n - 1] + phi[n];48        }49        let main = novelty_main();50        let samples = ((high - LOW) * per_octave as f64).round() as usize + 1;51        let heights: Vec<f64> = (0..samples)52            .map(|k| LOW + k as f64 / per_octave as f64)53            .collect();54        let smooth = heights55            .iter()56            .map(|&j| {57                let y = 2.0f64.powf(-j);58                smoothed_novelty(&phi, y, main) / y.powf(1.5)59            })60            .collect();61        let rough = heights62            .iter()63            .map(|&j| {64                let y = 2.0f64.powf(-j);65                sharp_novelty(&prefix, y) / y66            })67            .collect();68        let line = Line::new();69        let gammas = line.zeros(zeros);70        let coef = line.novelty_coefficients(&gammas);71        Ok(Novelty {72            heights,73            smooth,74            rough,75            gammas,76            coef,77            sieve,78        })79    }80    fn first(&self, count: usize) -> Result<usize, Fault> {81        if count > self.gammas.len() {82            return Err(Fault::new(format!(83                "the meter holds {} zeros.",84                self.gammas.len()85            )));86        }87        Ok(count)88    }89    /// Returns j at every sample, the log base two of one over y.90    pub fn heights(&self) -> Vec<f64> {91        self.heights.clone()92    }93    /// Returns the error at every sample: the smoothed error over y to the three halves, or the sharp window's error over y.94    pub fn dots(&self, sharp: bool) -> Vec<f64> {95        if sharp {96            self.rough.clone()97        } else {98            self.smooth.clone()99        }100    }101    /// Returns the ordinate of every zero the meter holds.102    pub fn gammas(&self) -> Vec<f64> {103        self.gammas.clone()104    }105    /// Returns the modulus of every zero's wave coefficient.106    pub fn amplitudes(&self) -> Vec<f64> {107        self.coef.iter().map(|c| c.abs()).collect()108    }109    /// Sums the waves of the first zeros at the given heights j, on the scale of the dots: over y to the three halves, or over y when sharp, where the waves shrink by the root of y.110    pub fn wave(&self, count: usize, sharp: bool, heights: &[f64]) -> Result<Vec<f64>, Fault> {111        let count = self.first(count)?;112        let log2 = 2.0f64.ln();113        Ok(heights114            .iter()115            .map(|&j| {116                let sum = novelty_wave(&self.gammas[..count], &self.coef[..count], -j * log2);117                if sharp {118                    sum * 2.0f64.powf(-0.5 * j)119                } else {120                    sum121                }122            })123            .collect())124    }125    /// Returns the largest gap between the smoothed dots and the wave of the first zeros, over the largest dot.126    pub fn miss(&self, count: usize) -> Result<f64, Fault> {127        let wave = self.wave(count, false, &self.heights)?;128        let peak = self.smooth.iter().fold(0.0f64, |a, v| a.max(v.abs()));129        let gap = self130            .smooth131            .iter()132            .zip(&wave)133            .map(|(d, w)| (d - w).abs())134            .fold(0.0f64, f64::max);135        Ok(gap / peak)136    }137    /// Returns the reach of the totient sieve.138    pub fn sieve(&self) -> usize {139        self.sieve140    }141}