Recent searches
Search options
`expm1` and `log1p` don't exist in math libs purely because `exp(x) - 1` and `log(1 + x)` are common in computations (although they are) or to save a single add/subtract, but rather because they're common and evaluating them directly as `exp(x) - 1` or `log(1 + x)` respectively is very inaccurate when |x| << 1.
It's probably easiest to see with exp which has a power series that converges everywhere:
exp(x) = sum(k=0;inf) x^k/k! = 1 + x + x^2/2 + x^3/3! + x^4/4! + ...
Note the leading term is 1+. If |x| << 1, then "1 + x" will lose most of the digits of x, and "(1 + x) - 1" will therefore be quite inaccurate.
In this exact case, expm1 can conceptually compute the result using the first few terms (if it used the power series directly, which it probably doesn't, but: details!) and never add/subtract 1.
Much like physicists favorite small-angle approximation sin(x) ≈ x for small x, exp(x) - 1 ≈ x for small x (for pretty much the same reason too *waves hand at complex exp*), so if exp(x)-1 is what you want, your friendly neighborhood math library knows which values of x are in the ~danger zone~ and will do the right thing to give you a result with high accuracy.
log1p works out the same way (it's the inverse of expm1).
@rygorous See I think they missed a trick here. For small x, an IEEE-compliant exp(x) is much quicker to evaluate than exp(x)-1, so you can win more benchmarks.
@TomF math libraries know _all_ the tricks https://github.com/rutgers-apl/rlibm-all/blob/main/source/float/exp.c#L17
@rygorous <chef's kiss>
@TomF it's not just rlibm either
here's core-math expf32 https://gitlab.inria.fr/core-math/core-math/-/blob/master/src/binary32/exp/expf.c#L62
musl https://git.musl-libc.org/cgit/musl/tree/src/math/exp.c#n80
glibc https://github.com/bminor/glibc/blob/master/sysdeps/ieee754/dbl-64/e_exp.c#L108
@rygorous The fact that these functions exist is news to me, thank you for sharing!