Function normal_logcdf 

pub fn normal_logcdf(x: f64) -> f64

Numerically stable ln Φ(x) for the standard normal CDF. For x ≥ 0, evaluates ln(1 - 0.5 erfc(x/sqrt(2))) with ln_1p, retaining the small negative result after Φ(x) itself rounds to one. For x < 0, rewrites ln Φ(x) = −u² + ln(½·erfcx(u)), u = −x/√2, which preserves digits throughout the representable left tail without a probability floor. Returns the corresponding IEEE limit at infinities and propagates NaN.

Examples found in repository

examples/special_audit_dump.rs (line 80)

fn main() {
    // ---- Bessel channels -------------------------------------------------
    for eta in geometric_grid(1e-6, 1e12, 400) {
        let (centered, ratio, d1) = special::bessel_i0_centered_terms(eta);
        emit("bessel_centered_log", eta, centered);
        emit("bessel_ratio", eta, ratio);
        emit("bessel_d1", eta, d1);
        emit(
            "bessel_d2",
            eta,
            special::bessel_i0_centered_second_log_derivative_from_log_abs(eta.ln()),
        );
    }
    // Dense sweep across the ascending/asymptotic seam at 20.
    for eta in linear_grid(15.0, 25.0, 201) {
        let (centered, ratio, d1) = special::bessel_i0_centered_terms(eta);
        emit("bessel_centered_log", eta, centered);
        emit("bessel_ratio", eta, ratio);
        emit("bessel_d1", eta, d1);
        emit(
            "bessel_d2",
            eta,
            special::bessel_i0_centered_second_log_derivative_from_log_abs(eta.ln()),
        );
    }

    // ---- Polygamma family ------------------------------------------------
    for x in geometric_grid(1e-8, 1e10, 400) {
        emit("digamma", x, special::digamma(x));
        emit("trigamma", x, special::trigamma(x));
        emit("tetragamma", x, special::tetragamma(x));
        emit("pentagamma", x, special::pentagamma(x));
    }
    for x in linear_grid(0.5, 30.0, 300) {
        emit("digamma", x, special::digamma(x));
        emit("trigamma", x, special::trigamma(x));
        emit("tetragamma", x, special::tetragamma(x));
        emit("pentagamma", x, special::pentagamma(x));
    }

    // ---- Normal-distribution channels ------------------------------------
    let mut normal_args: Vec<f64> = linear_grid(-40.0, 40.0, 801);
    normal_args.extend(geometric_grid(1e-8, 1e2, 200));
    normal_args.extend(geometric_grid(1e-8, 1e2, 200).into_iter().map(|v| -v));
    for x in normal_args {
        emit("normal_pdf", x, prob::normal_pdf(x));
        emit("normal_cdf", x, prob::normal_cdf(x));
        emit("normal_logcdf", x, prob::normal_logcdf(x));
        emit("normal_logsf", x, prob::normal_logsf(x));
        let (log_cdf, mills) = prob::signed_probit_logcdf_and_mills_ratio(x);
        emit("probit_logcdf", x, log_cdf);
        emit("probit_mills", x, mills);
        let derivatives = prob::normal_logcdf_derivatives(x);
        for (order, value) in derivatives.iter().enumerate() {
            emit(&format!("logcdf_d{order}"), x, *value);
        }
        if x >= 0.0 {
            emit("erfcx", x, prob::erfcx_nonnegative(x));
        }
    }
    for x in geometric_grid(1e-3, 1e6, 400) {
        emit("erfcx", x, prob::erfcx_nonnegative(x));
    }

    // ---- log1mexp --------------------------------------------------------
    for a in geometric_grid(1e-14, 50.0, 300) {
        emit("log1mexp", a, prob::log1mexp_positive(a));
    }

    // ---- Normal quantile -------------------------------------------------
    for p in geometric_grid(1e-300, 0.5, 400) {
        if let Ok(q) = prob::standard_normal_quantile(p) {
            emit("normal_quantile", p, q);
        }
    }
    for p in linear_grid(0.001, 0.999, 400) {
        if let Ok(q) = prob::standard_normal_quantile(p) {
            emit("normal_quantile", p, q);
        }
    }
    for log_p in linear_grid(-700.0, -1e-6, 400) {
        if let Ok(q) = prob::standard_normal_quantile_from_log_cdf(log_p) {
            emit("normal_quantile_from_log", log_p, q);
        }
    }

    // ---- Gauss-Legendre --------------------------------------------------
    for n in [
        3usize, 4, 5, 7, 8, 12, 15, 16, 20, 24, 31, 32, 40, 48, 63, 64, 80, 96, 100, 127, 128, 160,
        200, 256,
    ] {
        let (nodes, weights) = special::gauss_legendre(n);
        for (i, (node, weight)) in nodes.iter().zip(weights.iter()).enumerate() {
            println!("gl_node\t{n}\t{i}\t{node:e}");
            println!("gl_weight\t{n}\t{i}\t{weight:e}");
        }
    }

    // ---- Binomial coefficient --------------------------------------------
    for n in 0usize..=60 {
        for k in 0..=n {
            println!(
                "binomial\t{n}\t{k}\t{:e}",
                special::binomial_coefficient_f64(n, k)
            );
        }
    }
}