statrs should route all float ops through `num_traits` which means we'll rely on libm or std based. For a brief period, I thought, we should consider using `libm` since that's a growing default for the no_std environment, but then how complete does a math function backend need to be for consistency? libm implements special functions that aren't available in std and there are other crates that maintain special functions, this exposes cases where statrs and other functions aren't the same. The question I have on this is, > What footguns can we mitigate for users who need some statistics functions that wouldn't already be mitigating if they had a crash course in floating point? I think the community convergence onto a special functions crate could reduce some of this, since the number of places needed to make consistent choices could be reduced. Right now `num_traits` has this role for the overlap between `std` and `libm`. Perhaps a crate would be able to implement sane trait-level defaults while exposing specific algorithms as functions where the backend selection would only control the trait-level would be useful. Then we can contribute upstream to broader numerics needs as well. [`special`](https://docs.rs/special/0.14.1/special/) exposes its functions as traits on float and double atop `libm` [`puruspe`](https://crates.io/crates/puruspe) is also in this space, author is that of `peroxide` (which to me, it feels like Matlab in Rust). I don't know what direction to take yet, but some input on if this would be an issue someone actually stumbles on would be nice to know. If this is fruitful, would also try to get input from stainless-steel and Axect. --- Here's a summary on what we actually use and whether its available elsewhere. Note that no_std environment will typically imply `libm` version of it could be consumed by a user of `statrs`. | Function statrs needs | Where used | Hand-rolled? | In `libm`? | In `puruspe`? | In `special`? | |---|---|---|---|---|---| | `erf` / `erfc` | `function/erf.rs` | Boost-derived rational/continued-fraction (`erf_impl`) | `libm::erf`, `libm::erfc` | `erf`, `erfc` | `Error::error`, `Error::compl_error` | | `erf_inv` / `erfc_inv` | `function/erf.rs` | Boost-derived (`erf_inv_impl`) | - | `inverf`, `inverfc` | `Error::inv_error` only (no `erfc_inv`) | | `ln_gamma` | `function/gamma.rs` | Lanczos approx (Pugh, 2004) | `libm::lgamma` / `lgamma_r` | `ln_gamma` | `Gamma::ln_gamma` | | `gamma` (Γ) | `function/gamma.rs` | Lanczos approx (Pugh, 2004) | `libm::tgamma` | `gamma` | `Gamma::gamma` | | `digamma` / `inv_digamma` | `function/gamma.rs` | Algorithm AS 103 | - | - | `Gamma::digamma` only (no `inv_digamma`; also adds `trigamma`) | | incomplete gamma (`gamma_ui/li/ur/lr`) | `function/gamma.rs` | continued-fraction (Cephes-style, uncited) | - | `gammp`, `gammq`, `invgammp` | `Gamma::inc_gamma` | | `ln_beta`, `beta` | `function/beta.rs` | built from `ln_gamma` | - | `beta` | `Beta::ln_beta` only (no plain `beta`) | | incomplete/regularized beta (`beta_inc`, `beta_reg`, `inv_beta_reg`) | `function/beta.rs` | Algorithm AS 63/64/109 | - | `betai`, `invbetai` | `Beta::inc_beta`, `Beta::inv_inc_beta` | | `logistic`, `logit` | `function/logistic.rs` | trivial, built from `exp`/`ln` | - | - | - | | `factorial`, `ln_factorial`, `binomial`, `multinomial` | `function/factorial.rs` | lookup table + `ln_gamma` | - | - | - | | `harmonic`, `gen_harmonic` | `function/harmonic.rs` | definition | - | - | - |