[Tcl-bugs] [ tcl-Bugs-2143217 ] isqrt 'make check' failure

[SourceForge.net <no...@so...>]

Bugs item #2143217, was opened at 2008-10-02 22:43
Message generated for change (Tracker Item Submitted) made by Item Submitter
You can respond by visiting: 
https://sourceforge.net/tracker/?func=detail&atid=110894&aid=2143217&group_id=10894

Please note that this message will contain a full copy of the comment thread,
including the initial issue submission, for this request,
not just the latest update.
Category: 56. LibTomMath
Group: current: 8.6a2
Status: Open
Resolution: None
Priority: 5
Private: No
Submitted By: Nobody/Anonymous (nobody)
Assigned to: Kevin B KENNY (kennykb)
Summary: isqrt 'make check'  failure

Initial Comment:
'make check' produces -

==== expr-47.12 isqrt of various sizes of integer FAILED
==== Contents of test case:

    set faults 0
    for {set i 0} {$faults < 10 && $i <= 1024} {incr i} {
        set root [expr {1 << $i}]
        set rm1 [expr {$root - 1}]
        set arg [expr {1 << (2 * $i)}]
        set tval [expr {isqrt($arg-1)}]
        if {$tval != $rm1} {
            append trouble "i = " $i ": isqrt(" $arg "-1) == " $tval  " != " $rm1 "\n"
            incr faults
        }
        set tval [expr {isqrt($arg)}]
        if {$tval != $root} {
            append trouble "i = " $i ": isqrt(" $arg ") == " $tval  " != " $root "\n"
            incr faults
        }
        set tval [expr {isqrt($arg+1)}]
        if {$tval != $root} {
            append trouble "i = " $i ": isqrt(" $arg "+1) == " $tval  " != " $root "\n"
            incr faults
        }
    }
    set trouble

---- Result was:
i = 28: isqrt(72057594037927936-1) == 279620267 != 268435455
i = 28: isqrt(72057594037927936+1) == 279620267 != 268435456
i = 29: isqrt(288230376151711744-1) == 545818760 != 536870911
i = 29: isqrt(288230376151711744+1) == 545818760 != 536870912
i = 30: isqrt(1152921504606846976-1) == 1079707056 != 1073741823
i = 30: isqrt(1152921504606846976+1) == 1079707056 != 1073741824
i = 31: isqrt(4611686018427387904-1) == 2150992608 != 2147483647
i = 31: isqrt(4611686018427387904+1) == 2150992608 != 2147483648
i = 32: isqrt(18446744073709551616-1) == 4296881274 != 4294967295
i = 41: isqrt(4835703278458516698824704-1) == 2199023259647 != 2199023255551

---- Result should have been (exact matching):

==== expr-47.12 FAILED

This is with gcc-4.2.1, -O2.  Without -O2, the test passes.

I believe the cause is in bn_mp_sqrt.c; the following patch fixes it.

--- tcl8.6a2/libtommath/bn_mp_sqrt.c    2006-12-01 05:48:23.000000000 +0000
+++ tcl8.6a2-patched/libtommath/bn_mp_sqrt.c    2008-10-01 00:02:55.000000000 +0100
@@ -69,7 +69,7 @@
   if (dig) {
       t1.used = i+2;
       d -= ldexp((double) dig, DIGIT_BIT);
-      if (d != 0.0) {
+      if (d >= 1.0) {
          t1.dp[i+1] = dig;
          t1.dp[i] = ((mp_digit) d) - 1;
       } else {

The reason is that d can become small, much less than 1.0 but not exactly 0.0 - in that case t1 is not set to a valid representation.
t1.dp[i] is -1, outside the range of a 28bit mp_digit.
Without -O2, d is exactly 0.0 for the first failing case, hence t1 is valid and the final result is correct.
This could be due to intermediate results being rounded to 64 bit double, but with -O2 everything is kept as 80 bit long double.

With attached code extracted from bn_mp_isqrt.c, I get
% make CC=/usr/local/gcc-4.2.1/bin/gcc CFLAGS=-O2 LDFLAGS=-lm isqrt && ./isqrt
/usr/local/gcc-4.2.1/bin/gcc -O2  -lm  isqrt.c   -o isqrt
d=-1.86265e-09=-0x1p-29 dig=1


System is
Linux piglet 2.4.36 #1 Tue Jan 22 01:12:23 GMT 2008 i686 unknown unknown GNU/Linux


Malc.
mal...@vi...

----------------------------------------------------------------------

You can respond by visiting: 
https://sourceforge.net/tracker/?func=detail&atid=110894&aid=2143217&group_id=10894