Several default methods for floats are wrong

[fingolfin]

Several default methods for floats in the GAP library are wrong, and can e.g. lead to mathematically inconsistent implementations of floateans. Possible actions include:

• removing these default operations
• introducing additional filters (e.g. to differentiate "real" from "complex" floats generically, and also "interval" floats from "non-interval" floats.

Here are some of the functions I consider wrong:

InstallMethod( AbsoluteValue, "for floats", [ IsFloat ], -1,
        function ( x )
    if x < Zero(x) then return -x; else return x; fi;
end );

InstallMethod( Norm, "for floats", [ IsFloat ], -1,
        function ( x )
    return x*x;
end );

InstallMethod( Argument, "for floats", [ IsFloat ], -1,
        function ( x )
    return Zero(x);
end );

These are wrong if x is a complex float. Since there doesn't seem a way to distinguish complex floats from real ones (e.g. no filter for that, as far as I can tell), this can easily lead to mathematically inconsistent implementations of floateans.

InstallMethod( SignFloat, "for floats", [ IsFloat ], -1,
        function ( x )
    if x < Zero(x) then return -1; elif IsZero(x) then return 0; else return 1; fi;
end );

This leads to incorrect results:

gap> SignFloat(-0.0);
0
gap> nan:=0.0/0.0;
nan
gap> SignFloat(nan);
1

For several operations, I also wonder whether it is wise to provide numerically bad default implementations. E.g. rounding is a very sensitive thing, and while e.g. Floor(x+MakeFloat(x,1/2)) as default for Round may seem compelling at a first glance, it gives incorrect results for a large number of floats. Example:

gap> a:=2^52+1;; x:=a*1.;; Int(x) = a; y:=Floor(x+MakeFloat(x,1/2));; x=y;
true
false
gap> x-y;
-1.

So an exact integral float is rounded incorrectly. There are more problems, for details see e.g. http://blog.frama-c.com/index.php?post/2013/05/02/nearbyintf1

A more extensive test suite would be helpful in all this, too.

[laurentbartholdi]

Dear @Max: I agree with you that this could be done better. The idea was to have at least something, but of course it must not give wrong results. There are different issues, as you sorted them out:

• there should be a subcategory "IsRealFloat" on which < should be defined and preserve the field operations. Norm, Argument, AbsoluteValue should have a default only on these

• the sign of -0 should be 0; but the sign of NaN is undefined in IEEE754, so what's wrong with GAP returning 1?

• I agree that losing an ulp in rounding is unfortunate; but there's no easy code that rounds correctly. It is best not to have any default? I'm happy with either. The link you give does not point to directly usable code, either: it doesn't handle negative values, and it relies on 52 bits of mantissa, which is precisely what I want to avoid in a generic default method.

[fingolfin]

• Re IsRealFloat -- that sounds good. I guess there'd also be an IsComplexFloat. To go along. And, independently, perhaps an IsIntervalFloat.

• I'd expect SignFloat(NaN) to give an error, but of course you are right, it is also valid to return 1, or fail, or anything else.

• The link on rounding was not meant as a suggestion for how to fix this, but rather as a source of additional examples (and I am not even claiming it's the best for that, I just googled for something with a few more examples and that's the first I found which quickly got to the point). Personally, I'd rather not provide a default method than to provide one which is sometimes incorrect.