Re: Slaying the HYPOTamus

On Mon, Aug 24, 2009 at 6:52 PM, David Fetter<david(at)fetter(dot)org> wrote:
> double hypot( double x, double y )
> {
>    double yx;
>
>    if( isinf(x) || isinf(y) )
>    return get_float8_infinity();
>
>    if( isnan(x) || isnan(y) )
>    return get_float8_nan();

For what it's worth though, check out the code in float.c to see how
postgres handles floating point overflows. I'm not sure whether it's
forced by the standard, we discussed and settled on this behaviour, or
the person who wrote it just thought it was a good idea but we may as
well be consistent.

The behaviour we have now is that if you pass Inf or -Inf then we
happily return Inf or -Inf (or whatever the result is) as appropriate.
But if the operation overflows despite being passed reasonable values
then it throws an error.

Afaict hypot() can still overflow even with this new coding if you
have, say, hypot(MAX_FLOAT, MAX_FLOAT) which will return MAX_FLOAT *
sqrt(2). In that case we should throw an error unless either input was
Inf.

Also, the question arises what should be returned for hypot(Inf,NaN)
which your code returns Inf for. Empirically, it seems the normal
floating point behaviour is to return NaN so I think the NaN test
should be first.

Lastly I find the swap kind of confusing and also think it might
perform less well than just having one branch and simple logic in each
branch. This is just a style question that you could see either way
though, the performance difference is probably immeasurable if even if
my guess is right.

so I would suggest:

if (isnan(x) || isnan(y)
return float8_get_nan();
else if (isinf(x) || isinf(y))
return float8_get_inf();
else if (x == 0.0 && y == 0.0)
return 0.0;

x = fabs(x);
y = fabs(y);

if (x > y)
retval = x * sqrt((y/x)*(y/x) + 1);
else
retval = y * sqrt((x/y)*(x/y) + 1);

if (isinf(retval))
ereport(overflow...)

return retval;
}

--
greg
http://mit.edu/~gsstark/resume.pdf