mirror/guile
mirror/guile
1
Fork 0
mirror of https://git.savannah.gnu.org/git/guile.git synced 2025-05-28 16:00:22 +02:00
Code Activity

Fix bugs in numerical equality predicate.

Browse source 
* libguile/numbers.c (scm_num_eq_p): Fix bug comparing fractions to
  infinities (reported by Göran Weinholt <goran@weinholt.se>).  Fix
  erroneous comment describing the logic behind inum/flonum comparison.
  Use similar logic for inum/complex comparison to avoid rounding
  errors.  Make minor indentation fixes and simplifications.

* test-suite/tests/numbers.test (=): Add tests.
This commit is contained in:
Mark H Weaver 2013-07-16 00:18:40 -04:00
parent 4cc2e41cf7
commit 0132928891
2 changed files with 58 additions and 28 deletions
Download patch file Download diff file Expand all files Collapse all files

60
libguile/numbers.c
View file

@ -6542,9 +6542,11 @@ scm_num_eq_p (SCM x, SCM y)
to a double and compare. to a double and compare.
But on a 64-bit system an inum is bigger than a double and But on a 64-bit system an inum is bigger than a double and
casting it to a double (call that dxx) will round. dxx is at casting it to a double (call that dxx) will round.
worst 1 bigger or smaller than xx, so if dxx==yy we know yy is Although dxx will not in general be equal to xx, dxx will
an integer and fits a long. So we cast yy to a long and always be an integer and within a factor of 2 of xx, so if
dxx==yy, we know that yy is an integer and fits in
scm_t_signed_bits. So we cast yy to scm_t_signed_bits and
compare with plain xx. compare with plain xx.
An alternative (for any size system actually) would be to check An alternative (for any size system actually) would be to check
@ -6559,8 +6561,14 @@ scm_num_eq_p (SCM x, SCM y)
|| xx == (scm_t_signed_bits) yy)); || xx == (scm_t_signed_bits) yy));
} }
else if (SCM_COMPLEXP (y)) else if (SCM_COMPLEXP (y))
return scm_from_bool (((double) xx == SCM_COMPLEX_REAL (y)) {
&& (0.0 == SCM_COMPLEX_IMAG (y))); /* see comments with inum/real above */
double ry = SCM_COMPLEX_REAL (y);
return scm_from_bool ((double) xx == ry
&& 0.0 == SCM_COMPLEX_IMAG (y)
&& (DBL_MANT_DIG >= SCM_I_FIXNUM_BIT-1
|| xx == (scm_t_signed_bits) ry));
}
else if (SCM_FRACTIONP (y)) else if (SCM_FRACTIONP (y))
return SCM_BOOL_F; return SCM_BOOL_F;
else else
@ -6615,24 +6623,21 @@ scm_num_eq_p (SCM x, SCM y)
else if (SCM_BIGP (y)) else if (SCM_BIGP (y))
{ {
int cmp; int cmp;
if (isnan (SCM_REAL_VALUE (x))) if (isnan (xx))
return SCM_BOOL_F; return SCM_BOOL_F;
cmp = xmpz_cmp_d (SCM_I_BIG_MPZ (y), SCM_REAL_VALUE (x)); cmp = xmpz_cmp_d (SCM_I_BIG_MPZ (y), xx);
scm_remember_upto_here_1 (y); scm_remember_upto_here_1 (y);
return scm_from_bool (0 == cmp); return scm_from_bool (0 == cmp);
} }
else if (SCM_REALP (y)) else if (SCM_REALP (y))
return scm_from_bool (SCM_REAL_VALUE (x) == SCM_REAL_VALUE (y)); return scm_from_bool (xx == SCM_REAL_VALUE (y));
else if (SCM_COMPLEXP (y)) else if (SCM_COMPLEXP (y))
return scm_from_bool ((SCM_REAL_VALUE (x) == SCM_COMPLEX_REAL (y)) return scm_from_bool ((xx == SCM_COMPLEX_REAL (y))
&& (0.0 == SCM_COMPLEX_IMAG (y))); && (0.0 == SCM_COMPLEX_IMAG (y)));
else if (SCM_FRACTIONP (y)) else if (SCM_FRACTIONP (y))
{ {
double xx = SCM_REAL_VALUE (x); if (isnan (xx) || isinf (xx))
if (isnan (xx))
return SCM_BOOL_F; return SCM_BOOL_F;
if (isinf (xx))
return scm_from_bool (xx < 0.0);
x = scm_inexact_to_exact (x); /* with x as frac or int */ x = scm_inexact_to_exact (x); /* with x as frac or int */
goto again; goto again;
} }
@ -6642,8 +6647,15 @@ scm_num_eq_p (SCM x, SCM y)
else if (SCM_COMPLEXP (x)) else if (SCM_COMPLEXP (x))
{ {
if (SCM_I_INUMP (y)) if (SCM_I_INUMP (y))
return scm_from_bool ((SCM_COMPLEX_REAL (x) == (double) SCM_I_INUM (y)) {
&& (SCM_COMPLEX_IMAG (x) == 0.0)); /* see comments with inum/real above */
double rx = SCM_COMPLEX_REAL (x);
scm_t_signed_bits yy = SCM_I_INUM (y);
return scm_from_bool (rx == (double) yy
&& 0.0 == SCM_COMPLEX_IMAG (x)
&& (DBL_MANT_DIG >= SCM_I_FIXNUM_BIT-1
|| (scm_t_signed_bits) rx == yy));
}
else if (SCM_BIGP (y)) else if (SCM_BIGP (y))
{ {
int cmp; int cmp;
@ -6657,20 +6669,18 @@ scm_num_eq_p (SCM x, SCM y)
} }
else if (SCM_REALP (y)) else if (SCM_REALP (y))
return scm_from_bool ((SCM_COMPLEX_REAL (x) == SCM_REAL_VALUE (y)) return scm_from_bool ((SCM_COMPLEX_REAL (x) == SCM_REAL_VALUE (y))
&& (SCM_COMPLEX_IMAG (x) == 0.0)); && (SCM_COMPLEX_IMAG (x) == 0.0));
else if (SCM_COMPLEXP (y)) else if (SCM_COMPLEXP (y))
return scm_from_bool ((SCM_COMPLEX_REAL (x) == SCM_COMPLEX_REAL (y)) return scm_from_bool ((SCM_COMPLEX_REAL (x) == SCM_COMPLEX_REAL (y))
&& (SCM_COMPLEX_IMAG (x) == SCM_COMPLEX_IMAG (y))); && (SCM_COMPLEX_IMAG (x) == SCM_COMPLEX_IMAG (y)));
else if (SCM_FRACTIONP (y)) else if (SCM_FRACTIONP (y))
{ {
double xx; double xx;
if (SCM_COMPLEX_IMAG (x) != 0.0) if (SCM_COMPLEX_IMAG (x) != 0.0)
return SCM_BOOL_F; return SCM_BOOL_F;
xx = SCM_COMPLEX_REAL (x); xx = SCM_COMPLEX_REAL (x);
if (isnan (xx)) if (isnan (xx) || isinf (xx))
return SCM_BOOL_F; return SCM_BOOL_F;
if (isinf (xx))
return scm_from_bool (xx < 0.0);
x = scm_inexact_to_exact (x); /* with x as frac or int */ x = scm_inexact_to_exact (x); /* with x as frac or int */
goto again; goto again;
} }
@ -6686,10 +6696,8 @@ scm_num_eq_p (SCM x, SCM y)
else if (SCM_REALP (y)) else if (SCM_REALP (y))
{ {
double yy = SCM_REAL_VALUE (y); double yy = SCM_REAL_VALUE (y);
if (isnan (yy)) if (isnan (yy) || isinf (yy))
return SCM_BOOL_F; return SCM_BOOL_F;
if (isinf (yy))
return scm_from_bool (0.0 < yy);
y = scm_inexact_to_exact (y); /* with y as frac or int */ y = scm_inexact_to_exact (y); /* with y as frac or int */
goto again; goto again;
} }
@ -6699,10 +6707,8 @@ scm_num_eq_p (SCM x, SCM y)
if (SCM_COMPLEX_IMAG (y) != 0.0) if (SCM_COMPLEX_IMAG (y) != 0.0)
return SCM_BOOL_F; return SCM_BOOL_F;
yy = SCM_COMPLEX_REAL (y); yy = SCM_COMPLEX_REAL (y);
if (isnan (yy)) if (isnan (yy) || isinf(yy))
return SCM_BOOL_F; return SCM_BOOL_F;
if (isinf (yy))
return scm_from_bool (0.0 < yy);
y = scm_inexact_to_exact (y); /* with y as frac or int */ y = scm_inexact_to_exact (y); /* with y as frac or int */
goto again; goto again;
} }

26
test-suite/tests/numbers.test
View file

@ -2034,7 +2034,28 @@
(pass-if (not (= (ash-flo 1.0 58) (1- (ash 1 58))))) (pass-if (not (= (ash-flo 1.0 58) (1- (ash 1 58)))))
(pass-if (= (ash 1 58) (ash-flo 1.0 58))) (pass-if (= (ash 1 58) (ash-flo 1.0 58)))
(pass-if (not (= (1+ (ash 1 58)) (ash-flo 1.0 58)))) (pass-if (not (= (1+ (ash 1 58)) (ash-flo 1.0 58))))
(pass-if (not (= (1- (ash 1 58)) (ash-flo 1.0 58))))) (pass-if (not (= (1- (ash 1 58)) (ash-flo 1.0 58))))
;; prior to guile 2.0.10, inum/complex comparisons were done just by
;; converting the inum to a double, which on a 64-bit would round making
;; say inexact 2^58 appear equal to exact 2^58+1
(pass-if (= (+ +0.0i (ash-flo 1.0 58)) (ash 1 58)))
(pass-if (not (= (+ +0.0i (ash-flo 1.0 58)) (1+ (ash 1 58)))))
(pass-if (not (= (+ +0.0i (ash-flo 1.0 58)) (1- (ash 1 58)))))
(pass-if (= (ash 1 58) (+ +0.0i (ash-flo 1.0 58))))
(pass-if (not (= (1+ (ash 1 58)) (+ +0.0i (ash-flo 1.0 58)))))
(pass-if (not (= (1- (ash 1 58)) (+ +0.0i (ash-flo 1.0 58)))))
;; prior to guile 2.0.10, fraction/flonum and fraction/complex
;; comparisons mishandled infinities.
(pass-if (not (= 1/2 +inf.0)))
(pass-if (not (= 1/2 -inf.0)))
(pass-if (not (= +inf.0 1/2)))
(pass-if (not (= -inf.0 1/2)))
(pass-if (not (= 1/2 +inf.0+0.0i)))
(pass-if (not (= 1/2 -inf.0+0.0i)))
(pass-if (not (= +inf.0+0.0i 1/2)))
(pass-if (not (= -inf.0+0.0i 1/2))))
;;; ;;;
;;; < ;;; <
@ -2085,6 +2106,9 @@
(pass-if "n = 0.0" (pass-if "n = 0.0"
(not (< 0.0 0.0))) (not (< 0.0 0.0)))
(pass-if "n = -0.0"
(not (< 0.0 -0.0)))
(pass-if "n = 1" (pass-if "n = 1"
(< 0.0 1)) (< 0.0 1))