wiki-euler-mascheroni-constant.rs

1use mrlycore::errors::Result;2use mrlyfig::{ink, plot, save, Board};3use mrlynum::series::{harmonic, EULER};45const TOP: usize = 30;6const SAMPLES: usize = 1600;78fn main() -> Result<()> {9    let mut board = Board::square();10    let frame = board.frame(0.08);11    let left = 1.0;12    let right = (TOP + 1) as f64;13    let peak = harmonic(TOP);14    let across = |x: f64| frame.x + frame.w * (x - left) / (right - left);15    let up = |v: f64| frame.y + frame.h * (1.0 - v / peak);1617    let strip = ink::fade(ink::blue(), 0.28);18    let width = frame.w / SAMPLES as f64;19    for k in 0..SAMPLES {20        let x = left + (right - left) * (k as f64 + 0.5) / SAMPLES as f64;21        let crest = up(harmonic(x.floor() as usize));22        let foot = up(x.ln());23        board.rect(across(x) - width, crest, width * 2.0, foot - crest, strip);24    }2526    let mut curve = Vec::with_capacity(SAMPLES + 1);27    for k in 0..=SAMPLES {28        let x = left + (right - left) * k as f64 / SAMPLES as f64;29        curve.push((across(x), up(x.ln())));30    }31    board.polyline(&curve, 3.4, ink::blue());3233    let mut stair = Vec::with_capacity(2 * TOP);34    for n in 1..=TOP {35        let v = harmonic(n);36        stair.push((across(n as f64), up(v)));37        stair.push((across((n + 1) as f64), up(v)));38    }39    board.polyline(&stair, 3.4, ink::yellow());4041    plot::axis(&mut board, frame, ink::line());42    assert_eq!(stair.len(), 2 * TOP);43    assert!((peak - (TOP as f64).ln() - EULER).abs() < 0.02);44    save("wiki-euler-mascheroni-constant", &board)?;45    Ok(())46}