Understanding numeric data types in functions

Question

This is a follow-up to the question on double-float that I posted earlier. I apologize for what's probably a fundamental Lisp concept, but I haven't grasped it yet.

For this problem, I am using GNU CLISP 2.49 (2010-07-07).

Suppose I have the following function, which simply determines square root via Newton's method:

(defun sr (n eps)
  (when (>= n 0)
    (do ((x (/ n 2.0) (/ (+ x (/ n x)) 2.0)))
      ((< (abs (- (* x x) n)) eps)
       x))))

I can call this as follows:

> (sr 2 0.00001)
1.4142157

It's giving me single precision float (the default). Makes sense. Due to the lack of precision, if I make eps too small, it doesn't function properly and goes into an infinite loop:

> (sr 2 0.00000001)
[just sits there...]

If I call it with double precision values, I still get single precision results:

> (sr 2.0d0 0.00001d0)
1.4142157
> (sr 2.0d0 0.00000001d0)
[just sits there...]

But if I redefine my function as follows:

(defun sr (n eps)
  (when (>= n 0)
    (do ((x (/ n 2.0d0) (/ (+ x (/ n x)) 2.0d0)))
      ((< (abs (- (* x x) n)) eps)
       x))))

I then get double precision no matter how I feed it:

> (sr 2 0.00001)
1.4142156862745097d0

And now feeding it a smaller eps works due to the increased precision:

> (sr 2 0.00000001)
1.4142135623746899d0

So my question is: is the precision applied by the function totally driven by the precision I specify in the constants it is using in the arithmetic expressions it contains? And if so, what if there were no constants anywhere in the function? What then determines the precision of the calculations and the result?

---

ADDENDUM

I just retested this on SBCL 1.0.57-1.fc17 and I get much more expected results, per the documentation that @JoshuaTaylor quoted in the comment.

Answer 1

Looks to me like this is an incompatibility of CLISP with the ANSI CL Standard.

The ANSI CL standard requires the result to be the largest type of all arguments:

CL-USER 20 > (/ 2.0d0 2.0)
1.0D0

CLISP gives:

[4]> (/ 2.0d0 2.0)
1.0

This should be a double float.

You can change that in CLISP to standard ANSI CL behavior:

[9]> CUSTOM:*FLOATING-POINT-CONTAGION-ANSI*
NIL
[10]> (setf CUSTOM:*FLOATING-POINT-CONTAGION-ANSI* t)
T
[11]> (sr 2.0d0 0.00001d0)
1.4142156862745097d0

Now it gives back the desired result.

See the CLISP manual: 12.2.4.1. Rule of Float Precision Contagion