Raymond Toy pushed to branch rtoy-setexception-inexact at cmucl / cmucl

Commits:

4 changed files:

Changes:

  • src/lisp/k_sin.c
    --- a/src/lisp/k_sin.c
    +++ b/src/lisp/k_sin.c
    @@ -67,8 +67,14 @@ S6  =  1.58969099521155010221e-10; /* 0x3DE5D93A, 0x5ACFD57C */
     
             ux.d = x;
     	ix = ux.i[HIWORD]&0x7fffffff;	/* high word of x */
    -	if(ix<0x3e400000)			/* |x| < 2**-27 */
    -	   {if((int)x==0) return x;}		/* generate inexact */
    +	if(ix<0x3e400000) {			/* |x| < 2**-27 */
    +            /* return x inexact except 0 */
    +            if (x != 0) {
    +                fdlibm_setexception(x, FDLIBM_INEXACT);
    +            }
    +
    +            return x;
    +        }
     	z	=  x*x;
     	v	=  z*x;
     	r	=  S2+z*(S3+z*(S4+z*(S5+z*S6)));
    

  • src/lisp/k_tan.c
    --- a/src/lisp/k_tan.c
    +++ b/src/lisp/k_tan.c
    @@ -78,31 +78,34 @@ __kernel_tan(double x, double y, int iy) {
     	hx = ux.i[HIWORD];		/* high word of x */
     	ix = hx & 0x7fffffff;			/* high word of |x| */
     	if (ix < 0x3e300000) {			/* x < 2**-28 */
    -		if ((int) x == 0) {		/* generate inexact */
    -			if (((ix | ux.i[LOWORD]) | (iy + 1)) == 0)
    -				return one / fabs(x);
    -			else {
    -				if (iy == 1)
    -					return x;
    -				else {	/* compute -1 / (x+y) carefully */
    -					double a, t;
    +            /* return x inexact except 0 */
    +            if (x != 0) {
    +                fdlibm_setexception(x, FDLIBM_INEXACT);
    +            }
     
    -					z = w = x + y;
    -                                        uz.d = z;
    -					uz.i[LOWORD] = 0;
    -                                        z = ux.d;
    +	    if (((ix | ux.i[LOWORD]) | (iy + 1)) == 0)
    +		return one / fabs(x);
    +	    else {
    +		if (iy == 1)
    +		    return x;
    +		else {	/* compute -1 / (x+y) carefully */
    +		    double a, t;
    +
    +		    z = w = x + y;
    +		    uz.d = z;
    +		    uz.i[LOWORD] = 0;
    +		    z = ux.d;
                                             
    -					v = y - (z - x);
    -					t = a = -one / w;
    -                                        uz.d = t;
    -                                        uz.i[LOWORD] = 0;
    -                                        t = uz.d;
    +		    v = y - (z - x);
    +		    t = a = -one / w;
    +		    uz.d = t;
    +		    uz.i[LOWORD] = 0;
    +		    t = uz.d;
                                             
    -					s = one + t * z;
    -					return t + a * (s + t * v);
    -				}
    -			}
    +		    s = one + t * z;
    +		    return t + a * (s + t * v);
     		}
    +	    }
     	}
     	if (ix >= 0x3FE59428) {	/* |x| >= 0.6744 */
     		if (hx < 0) {
    

  • src/lisp/s_atan.c
    --- a/src/lisp/s_atan.c
    +++ b/src/lisp/s_atan.c
    @@ -104,7 +104,12 @@ huge   = 1.0e300;
     	    else     return -atanhi[3]-atanlo[3];
     	} if (ix < 0x3fdc0000) {	/* |x| < 0.4375 */
     	    if (ix < 0x3e200000) {	/* |x| < 2^-29 */
    -		if(huge+x>one) return x;	/* raise inexact */
    +		/* return x inexact except 0 */
    +		if (x != 0) {
    +		    fdlibm_setexception(x, FDLIBM_INEXACT);
    +		}
    +
    +		return x;
     	    }
     	    id = -1;
     	} else {
    

  • tests/fdlibm.lisp
    --- a/tests/fdlibm.lisp
    +++ b/tests/fdlibm.lisp
    @@ -258,7 +258,19 @@
     		(kernel:%atan *snan*))
       (assert-true (ext:float-nan-p (kernel:%atan *qnan*)))
       (kernel::with-float-traps-masked (:invalid)
    -    (assert-true (ext:float-nan-p (kernel:%atan *snan*)))))
    +    (assert-true (ext:float-nan-p (kernel:%atan *snan*))))
    +  ;; atan(x) = x for |x| < 2^-29, signaling inexact.
    +  (let ((x (scale-float 1d0 -30))
    +	(x0 0d0))
    +    (with-inexact-exception-enabled
    +	;; This must not throw an inexact exception because the result
    +	;; is exact when the arg is 0.
    +	(assert-eql 0d0 (kernel:%atan x0)))
    +    (with-inexact-exception-enabled
    +	;; This must throw an inexact exception for non-zero x even
    +	;; though the result is exactly x.
    +	(assert-error 'floating-point-inexact
    +		      (kernel:%atan x)))))
     
     (define-test %log10.exceptions
       (:tag :fdlibm)
    @@ -622,3 +634,34 @@
     	;; though the result is exactly x.
     	(assert-error 'floating-point-inexact
     		      (kernel:%cos x)))))
    +
    +(define-test %sin.exceptions
    +    (:tag :fdlibm)
    +  ;; sin(x) = x for |x| < 2^-27.  Signal inexact unless x = 0
    +  (let ((x (scale-float 1d0 -28))
    +	(x0 0d0))
    +    (with-inexact-exception-enabled
    +	;; This must not throw an inexact exception because the result
    +	;; is exact when the arg is 0.
    +	(assert-eql 0d0 (kernel:%sin x0)))
    +    (with-inexact-exception-enabled
    +	;; This must throw an inexact exception for non-zero x even
    +	;; though the result is exactly x.
    +	(assert-error 'floating-point-inexact
    +		      (kernel:%sin x)))))
    +
    +(define-test %tan.exceptions
    +    (:tag :fdlibm)
    +  ;; tan(x) = x for |x| < 2^-28.  Signal inexact unless x = 0
    +  (let ((x (scale-float 1d0 -29))
    +	(x0 0d0))
    +    (with-inexact-exception-enabled
    +	;; This must not throw an inexact exception because the result
    +	;; is exact when the arg is 0.
    +	(assert-eql 0d0 (kernel:%tan x0)))
    +    (with-inexact-exception-enabled
    +	;; This must throw an inexact exception for non-zero x even
    +	;; though the result is exactly x.
    +	(assert-error 'floating-point-inexact
    +		      (kernel:%tan x)))))
    +