ENTIER / Math.Floor bug

Post Reply
User avatar
Robert
Posts: 177
Joined: Sat Sep 28, 2013 11:04 am
Location: Edinburgh, Scotland

ENTIER / Math.Floor bug

Post by Robert »

Below I have a program and printout.

I expect ENTIER and Math.Floor to do the same thing numerically, but ENTIER gives the answer 94906267 and the Math.Floor gives 94906268.0. I don't think this has a reasonable explanation in terms of rounding or TYPE ranges. The Math Docu says they are identical, except for some small print which does not apply here .

Some analysis: (m − 1/2)^2 = m^2 − m + 1/4, which is less than n = m^2 − m + 1. Thus the square root of n is more than m − 1/2, so should round to m.
Obviously, using 64 bit REALs, one does not expect that Math.Sqrt is exact. According to my investigations it gives 9490626.5 exactly to 64 bits, but rather less in its internal 80 bit format, maybe 1 80-bit LSB less.

Summary: The Math library gives the correct answer, but only by making two errors. It thinks the square root is less than 94906267.5, when it is actually more, but rounds it up anyway.

Most CP REAL functions/operations work to the internal 80-bit accuracy. Math.Floor seems to be an exception, and probably should be fixed to align with ENTIER.

Any other opinions?

Code: Select all

PROCEDURE  EntierBug*;
  VAR
    m, n, v, w, y, z  :  REAL;
  BEGIN
    m  :=  94906268.;
    n  :=  m * m - m + 1.;

    v  :=  M.Sqrt (n);
    w  :=  M.Sqrt (n) - m;

    y  :=  ENTIER  (M.Sqrt (n) + 0.5);
    z  :=  M.Floor (M.Sqrt (n) + 0.5);

    HALT (60)
  END  EntierBug;
and a printout:
TRAP 60 (postcondition violated)

AaaKaraTest.EntierBug [00002AF1H]
.m REAL 94906268.0
.n REAL 9007199610781556.0
.v REAL 94906267.5
.w REAL -0.5000000013169483
.y REAL 94906267.0
.z REAL 94906268.0
Kernel.Call [000020BFH]
luowy
Posts: 87
Joined: Thu Dec 17, 2015 1:32 pm

Re: ENTIER / Math.Floor bug

Post by luowy »

this is a problem of precision loss in conversion from 80bit to 64bit;

1, the input:
in your example, actually,
v := M.Sqrt(n) + 0.5; (* 94906268.0 *)
but if you make a comparison , you will get:
v > M.Sqrt(n) + 0.5 is TRUE;
it means M.Sqrt(n) + 0.5 in 80bit is less the 94906268.0;


M.Sqrt(n) + 0.5 = 94906268.0 -delta, the delta is very small;
under the round nearest mode, it will convert to a 64 bit 94906268.0,


2,the procedure:
though Floor and ENTIER has same machine code, but the input parameter is not same:

in Floor(M.Sqrt(n) + 0.5), the input is 64bits,
as Floor() is a normal procedure, it need real64 paramter,
the result of input expression has to pop to stack from real80 to real64;

in ENTIER(M.Sqrt(n) + 0.5) , the input can be 80bits,
as ENTIER is buildin procedure, it can reuse the real80 reault of the input expresion directly;

3, the defference:
Floor() input a standard 94906268.0, after FRNDINT, the result is 94906268, less or equal 94906268.0, no need sub 1;
result is 94906268.0;

but, ENTIER input 94906268.0-delta, after FRNDINT, the result is 94906268, great than the input, need sub 1;
so the result is 94906268.0-1= 94906267.0;
User avatar
Robert
Posts: 177
Joined: Sat Sep 28, 2013 11:04 am
Location: Edinburgh, Scotland

Re: ENTIER / Math.Floor bug

Post by Robert »

Thanks for the helpful explanation.

My mistake was thinking:
Robert wrote:Most CP REAL functions/operations work to the internal 80-bit accuracy
In fact, as stated in Platform-Specific-Issues, REAL functions return 80-bit results, but they DO round their inputs, which I had forgotten.

Should we change the Docu for Math.Floor to read
This function is identical to ENTIER, except that it returns a REAL type value and it does not generate an integer overflow. Also ENTIER does not round its input to 64-bits, which can lead to different results.
Josef Templ
Posts: 262
Joined: Tue Sep 17, 2013 6:50 am

Re: ENTIER / Math.Floor bug

Post by Josef Templ »

could be added to the docu but it is very implementation specific. So I am not sure about that.

- Josef
User avatar
Robert
Posts: 177
Joined: Sat Sep 28, 2013 11:04 am
Location: Edinburgh, Scotland

Re: ENTIER / Math.Floor bug

Post by Robert »

On reflection I think a change is probably not required.

I was reading too much into the phrase "This function is identical to ENTIER ...". It is also slower!.
Post Reply