Bug in nexttoward/2 and/or float/1

Help!
1

I am trying to dig deeper in the problems here on the windows
platform. It seems there is a bug in nexttoward/2 and/or float/1.
I find that in SWI-Prolog 8.5.17:

/*  SWI-Prolog 8.5.17, Windows 10 */
?- X is 51.0**10.
X = 1.1904242382761301e+17.

?- X is float(51^10).
X = 1.19042423827613e+17.

So the nexttoward of 1.19042423827613e+17 should
be 1.1904242382761301e+17. But I get an instable result:

/* SWI-Prolog 8.5.17, Windows 10 */
?- X is 51^10, Z is float(X), Y is nexttoward(Z, 1E300).
X = 119042423827613001,
Z = 1.1904242382761301e+17,
Y = 1.1904242382761302e+17.

?- Z is float(51^10), Y is nexttoward(Z, 1E300).
Z = 1.19042423827613e+17,
Y = 1.1904242382761301e+17.

In the above there are two different behaviours of float/1,
although the same argument value is used. The first float/1
call yields 1.1904242382761301e+17, and the second float/1

call yields 1.1904242382761302e+17, so maybe nexttoward/2
is not to blame, rather float/1 and the way arithmetic expressions
are processed in SWI-Prolog?

2

The different behaviour of float/1 can be also
observed on Unix, not only on Windows. I get:

/* SWI-Prolog 8.4.2, WSL2 Unix */
?- X is 51^10, Z is float(X).
X = 119042423827613001,
Z = 1.1904242382761301e+17.

?- Z is float(51^10).
Z = 1.19042423827613e+17.

The same effect doesn’t happen in Scryer Prolog and Trealla Prolog,
but they both do wrong float/1 rounding (1.19042423827613e+17
instead 1.1904242382761301E17):

/* Scryer Prolog 0.9.0 */
?- X is 51^10, Z is float(X).
   X = 119042423827613001, Z = 1.19042423827613e17.
?- Z is float(51^10).
   Z = 1.19042423827613e17.

But Trealla Prolog has a problem with (^)/2, wrong result of
bigint pow() (i.e. 119042423827613008 instead
119042423827613001):

/* Trealla Prolog 2.2.10 */
?- X is 51^10, Z is float(X).
   X = 119042423827613008, Z = 1.19042423827613e+17.
?- Z is float(51^10).
   Z = 1.19042423827613e+17.
3

As a lot of this has recently been changed or is being changed, these are old results. The current bigint-to-float branch does all this correct, both on Windows and Linux.

4

Ok, will do some retesting if a new release arrives. Maybe
my share here was contributing some test cases? Also figured
out which float result is the correct one. Take these two results:

/*  SWI-Prolog 8.5.17, Windows 10 */
?- X is float(51^10).
X = 1.19042423827613e+17.

/* ECLiPSe Prolog 7.0.61 */
?- X is float(51^10).
X = 1.1904242382761301e+17

Which one is correct? One needs to think in binary
representation. We find, using a precise float to
bigint conversion:

?- X is integer(1.1904242382761301E17).
X = 119042423827613008.

?- X is integer(1.19042423827613E17).
X = 119042423827612992.

Now which one is closer?

?- X is 119042423827613001-119042423827613008.
X = -7.

?- X is 119042423827613001-119042423827612992.
X = 9.

So 1.1904242382761301E17 is the correct float. It even
doesn’t need a HALF_EVEN, or HALF_UP or HALF_DOWN
rule. Its basically the rule of the nearest. Yes or no?

5

Seems there are very few where this happens. Maybe this
helps tracking down the cause? Or maybe this is a known
artefact of something? I started with a simulation of mpz_get_d/2:

mpz_get_d(X, Y) :-
    M is msb(X),
    mpz_get_d(M, X, Y).

mpz_get_d(M, X, Y) :-
    M < 53, !, Y is float(X).
mpz_get_d(M, X, Y) :-
    Y is float(X>>(M-52))*(1<<(M-52)).

And then I checked whether mpz_get_d/2 behaves as float/1
in SWI-Prolog 8.5.17. And it depends how I cast it, whether I
compute N^10 inside the float argument or outside:

?- aggregate_all(count, (between(1,1000000,N),
    A is N^10, mpz_get_d(A,B), B =\= float(N^10)), C).
C = 0.

?- aggregate_all(count, (between(1,1000000,N),
    A is N^10, mpz_get_d(A,B), B =\= float(A)), C).
C = 8.