chickadee » math » fpexpm1

x

flonum

Like (fplog (+ 1.0 x)) and (- (fpexp x) 1.0), but accurate when x is small (within 1 ulp - see fpulp).

For example, one difficult input for (fplog (+ 1.0 x)) and (- (fpexp x) 1.0) is x = 1e-14, which fplog1p and fpexpm1 compute correctly:

> (fplog (+ 1.0 1e-14))
9.992007221626358e-15
> (fplog1p 1e-14)
9.99999999999995e-15
> (- (fpexp 1e-14) 1.0)
9.992007221626409e-15
> (fpexpm1 1e-14)
1.0000000000000049e-14

These functions are mutual inverses:

Notice that both graphs pass through the origin. Thus, inputs close to 0.0, around which flonums are particularly dense, result in outputs that are also close to 0.0. Further, both functions are approximately the identity function near 0.0, so the output density is approximately the same.

Many flonum functions defined in terms of fplog and fpexp become much more accurate when their defining expressions are put in terms of fplog1p and fpexpm1. The functions exported by this module and by math/special-functions use them extensively.

One notorious culprit is (fpexpt (- 1.0 x) y), when x is near 0.0. Computing it directly too often results in the wrong answer:

> (fpexpt (- 1.0 1e-20) 1e+20)
1.0

We should expect that multiplying a number just less than 1.0 by itself that many times would result in something less than 1.0. The problem comes from subtracting such a small number from 1.0 in the first place:

> (- 1.0 1e-20)
1.0

Fortunately, we can compute this correctly by putting the expression in terms of fplog1p, which avoids the error-prone subtraction:

> (fpexp (* 1e+20 (fplog1p (- 1e-20))))
0.36787944117144233

See fpexpt1p, which is more accurate still.