chickadee » math » lg+

logx

flonum

logy

flonum

Like (fplog (+ (fpexp logx) (fpexp logy))) and (fplog (- (fpexp logx) (fpexp logy))), respectively, but more accurate and less prone to overflow and underflow.

When logy > logx, lg- returns +nan.0. Both functions correctly treat -inf.0 as log-space 0.0.

To add more than two log-space numbers with the same guarantees, use lgsum.

Examples:

> (lg+ (fplog 0.5) (fplog 0.2))
-0.35667494393873234
> (fpexp (lg+ (fplog 0.5) (fplog 0.2)))
0.7000000000000001
> (lg- (fplog 0.5) (fplog 0.2))
-1.203972804325936
> (fpexp (lg- (fplog 0.5) (fplog 0.2)))
0.30000000000000004
> (lg- (fplog 0.2) (fplog 0.5))
+nan.0

Though more accurate than a naive implementation, both functions are prone to catastrophic cancellation (see fpulp-error) in regions where they output a value close to 0.0 (or log-space 1.0). While these outputs have high relative error, their absolute error is very low, and when exponentiated, nearly have just rounding error. Further, catastrophic cancellation is unavoidable when logx and logy themselves have error, which is by far the common case.

These are, of course, excuses—but for floating-point research generally. There are currently no reasonably fast algorithms for computing lg+ and lg- with low relative error. For now, if you need that kind of accuracy, use (math bigfloat) (Note: still unimplemented in CHICKEN).