This code is taken from Sussman and Wisdom's Structure and Interpretation of Classical Mechanics, its purpose is to derive (close to) the smallest positive floating point the host machine supports. https://github.com/hnarayanan/sicm/blob/e37f011db68f8efc51ae309cd61bf497b90970da/scmutils/src/kernel/numeric.scm
Running it in DrRacket results in 2.220446049250313e-016 on my machine.
My question, what causes this to even return a value? This code is tail recursive, and it makes sense at some point the computer can no longer divide by 2. Why does it not throw?
(define *machine-epsilon*
(let loop ((e 1.0))
(if (= 1.0 (+ e 1.0))
(* 2 e)
(loop (/ e 2)))))
*machine-epsilon*
This code is tail recursive, and it makes sense at some point the computer can no longer divide by 2. Why does it not throw?
No, the idea is different: at some point the computer still can divide by 2, but the result (e
) becomes indistinguishable from 0 [upd: in the context of floating-point addition only - very good point mentioned in the comment] (e + 1.0 = 1.0
, this is exactly what if
clause is checking). We know for sure that the previous e
was still greater than zero "from the machine point of view" (otherwise we wouldn't get to the current execution point), so we simply return e*2
.