yz-curvep-exs/src/ex.rs

/// explicit parametrizations of some curves
/// (t is a parameter in [0,1] with f(0) = f(1) periodic)

pub fn circle(t: f64) -> [f64; 2] {
    // t * τ
    let tt2pi = t * core::f64::consts::TAU;
    let (ttsin, ttcos) = tt2pi.sin_cos();
    [ttcos, ttsin]
}

pub fn circle_deriv(t: f64) -> [f64; 2] {
    use core::f64::consts::TAU;
    let tt2pi = t * TAU;
    let (ttsin, ttcos) = tt2pi.sin_cos();
    [-TAU * ttsin, TAU * ttcos]
}

pub fn ellipse(a: f64, b: f64) -> impl Fn(f64) -> [f64; 2] {
    move |t: f64| {
        let cdat = circle(t);
        [a * cdat[0], b * cdat[1]]
    }
}

pub fn intersect_eight(t: f64) -> [f64; 2] {
    //  __  __
    // /  \/  \
    // \  /\  /
    //  ''  ''

    let circle_data = circle(2.0 * t + 0.5);
    let l = circle_data[0] + 1.0;
    let l2 = if t >= 0.5 { -1.0 } else { 1.0 };
    [l * l2, circle_data[1]]
}

pub fn intersect_8pp(t: f64) -> [f64; 2] {
    // like `intersec_eight`, but we use 3/4 circles on the outer ends,
    // and a 90° intersect in the middle

    let sqrt2_16: f64 = 8.0 * (2.0_f64).sqrt();
    let sqrt2m1: f64 = (2.0_f64).sqrt() - 1.0;
    let afterh: bool = t >= 0.5;

    // circle1 : 1/16..7/16
    // circle2 : 9/16..15/16
    if (t >= 0.0625 && t < 0.4375) || (t >= 0.5625 && t < 0.9375) {
        let x8p = intersect_eight(t);
        let newx = x8p[0] + sqrt2m1 * if t >= 0.5 { -1.0 } else { 1.0 };
        [newx, x8p[1]]
    } else if (t < 0.0625) || (t >= 0.9375) {
        // /
        let tx = sqrt2_16 * if afterh { t - 1.0 } else { t };
        [tx, -tx]
    } else {
        // \
        let tx = sqrt2_16 * (t - 0.5);
        [-tx, -tx]
    }
}

// doesn't work yet...
pub fn stadion(t: f64, gap: f64) -> [f64; 2] {
    // circle1 : 1/8..3/8
    // circle2 : 5/8..7/8
    todo!();
    if t < 0.125 {
        [1.0, t * 4.0 * gap - 1.0]
    } else if t >= 0.125 && t < 0.375 {
        circle(2.0 * t - 0.25)
    } else if t >= 0.375 && t < 0.625 {
        [-1.0, 1.0 - t * 4.0 * gap]
    } else if t >= 0.625 && t < 0.875 {
        let cdat = circle(2.0 * t - 0.75);
        [cdat[0], cdat[1] - 2.0 * gap]
    } else {
        [1.0, (t - 1.0) * 4.0 * gap - 1.0]
    }
}

pub fn wobbly(t: f64, wob: u32) -> [f64; 2] {
    let tt2pi = t * core::f64::consts::TAU;
    let (ttsin, ttcos) = tt2pi.sin_cos();
    let (ttxsin, ttxcos) = (tt2pi * (wob as f64)).sin_cos();
    let wobinv = 1.0 / (wob as f64);
    [ttcos + wobinv * ttxcos, ttsin + wobinv * ttxsin]
}

pub fn lemniskate(t: f64, a: f64) -> [f64; 2] {
    let tt2pi = t * core::f64::consts::TAU;
    let (ttsin, ttcos) = tt2pi.sin_cos();
    [a * ttcos, a * ttcos * ttsin]
}

pub fn cassini_oval(a: f64, c: f64) -> impl Fn(f64) -> [f64; 2] {
    let c2 = c * c;
    let c4 = c2 * c2;
    let a2 = a * a;
    let a4 = a2 * a2;
    let cahyp = c.hypot(a);
    move |t: f64| {
        let tt2pi = t * core::f64::consts::TAU;
        let (ttsin, ttcos) = tt2pi.sin_cos();
        let denom = c2 + a2 * ttsin * ttsin;
        [
            (c2 * cahyp * ttcos) / denom,
            (cahyp * ttsin * (c4 - a4 * ttsin * ttsin).sqrt()) / denom,
        ]
    }
}

pub fn cassini_wobbly(a: f64, c: f64, wob: u32) -> impl Fn(f64) -> [f64; 2] {
    let cassini = cassini_oval(a, c);
    let wobx = (wob as f64) * core::f64::consts::TAU;
    let wobinv = 1.0 / (wob as f64);
    move |t: f64| {
        let ttx2pi = t * wobx;
        let (ttxsin, ttxcos) = (ttx2pi).sin_cos();
        let cas = cassini(t);
        [
            ttxcos.mul_add(wobinv, cas[0]),
            ttxsin.mul_add(wobinv, cas[1]),
        ]
    }
}

pub fn astroid_pedal(t: f64, a: f64) -> [f64; 2] {
    let tt2pi = t * core::f64::consts::TAU;
    let (ttsin, ttcos) = tt2pi.sin_cos();
    [a * ttcos * ttsin * ttsin, a * ttcos * ttcos * ttsin]
}

pub fn cassini_astroid(ca: f64, cc: f64, aa: f64) -> impl Fn(f64) -> [f64; 2] {
    let cassini = cassini_oval(ca, cc);
    move |t: f64| {
        let ap = astroid_pedal(t, aa);
        let cs = cassini(t);
        [ap[0] + cs[0], ap[1] + cs[1]]
    }
}

pub fn spiral(k: u32) -> impl Fn(f64) -> [f64; 2] {
    let x2pi = (k as f64) * core::f64::consts::TAU;
    move |t: f64| {
        let ttcos = (t * core::f64::consts::TAU).cos();
        let (ttxsin, ttxcos) = (t * x2pi).sin_cos();
        [ttcos * ttxcos, ttcos * ttxsin]
    }
}

pub fn spiral_1d(k: u32) -> impl Fn(f64) -> [f64; 2] {
    let k = k as f64;
    let x2pi = k * core::f64::consts::TAU;
    move |t: f64| {
        let ttcos = k * (t * core::f64::consts::TAU).cos();
        let ttxsin = (t * x2pi).sin();
        [ttcos, ttxsin]
    }
}