java - ¿Por qué la función hypot() es tan lenta?

c++ hypotenuse

Question

Hice algunas pruebas con C++ hypot() y Java Math.hypot . Ambos parecen ser significativamente más lentos que sqrt(a*a + b*b) . ¿Es eso por una mejor precisión? ¿Qué método para calcular una hipotenusa utiliza la función de hypot ? Sorprendentemente, no pude encontrar ningún indicio de bajo rendimiento en la documentación.

Answer 1

Aquí hay una implementación más rápida, cuyos resultados también están más cerca de java.lang.Math.hypot: (NB: para la implementación de Delorie, es necesario agregar el manejo de las entradas NaN y + -Infinity)

private static final double TWO_POW_450 = Double.longBitsToDouble(0x5C10000000000000L);
private static final double TWO_POW_N450 = Double.longBitsToDouble(0x23D0000000000000L);
private static final double TWO_POW_750 = Double.longBitsToDouble(0x6ED0000000000000L);
private static final double TWO_POW_N750 = Double.longBitsToDouble(0x1110000000000000L);
public static double hypot(double x, double y) {
x = Math.abs(x);
y = Math.abs(y);
if (y < x) {
double a = x;
x = y;
y = a;
} else if (!(y >= x)) { // Testing if we have some NaN.
if ((x == Double.POSITIVE_INFINITY) || (y == Double.POSITIVE_INFINITY)) {
return Double.POSITIVE_INFINITY;
} else {
return Double.NaN;
}
}
if (y-x == y) { // x too small to substract from y
return y;
} else {
double factor;
if (x > TWO_POW_450) { // 2^450 < x < y
x *= TWO_POW_N750;
y *= TWO_POW_N750;
factor = TWO_POW_750;
} else if (y < TWO_POW_N450) { // x < y < 2^-450
x *= TWO_POW_750;
y *= TWO_POW_750;
factor = TWO_POW_N750;
} else {
factor = 1.0;
}
return factor * Math.sqrt(x*x+y*y);
}
}

Answer 2

Encontré Math.hypot () para ser extremadamente lento. Descubrí que podía codificar una versión java rápida utilizando el mismo algoritmo que produce resultados idénticos. Esto puede ser usado para reemplazarlo.

/**
* <b>hypot</b>
* @param x
* @param y
* @return sqrt(x*x +y*y) without intermediate overflow or underflow.
* @Note {@link Math#hypot} is unnecessarily slow. This returns the identical result to
* Math.hypot with reasonable run times (~40 nsec vs. 800 nsec).
* <p>The logic for computing z is copied from "Freely Distributable Math Library"
* fdlibm''s e_hypot.c. This minimizes rounding error to provide 1 ulb accuracy.
*/
public static double hypot(double x, double y) {
if (Double.isInfinite(x) || Double.isInfinite(y)) return Double.POSITIVE_INFINITY;
if (Double.isNaN(x) || Double.isNaN(y)) return Double.NaN;
x = Math.abs(x);
y = Math.abs(y);
if (x < y) {
double d = x;
x = y;
y = d;
}
int xi = Math.getExponent(x);
int yi = Math.getExponent(y);
if (xi > yi + 27) return x;
int bias = 0;
if (xi > 510 || xi < -511) {
bias = xi;
x = Math.scalb(x, -bias);
y = Math.scalb(y, -bias);
}
// translated from "Freely Distributable Math Library" e_hypot.c to minimize rounding errors
double z = 0;
if (x > 2*y) {
double x1 = Double.longBitsToDouble(Double.doubleToLongBits(x) & 0xffffffff00000000L);
double x2 = x - x1;
z = Math.sqrt(x1*x1 + (y*y + x2*(x+x1)));
} else {
double t = 2 * x;
double t1 = Double.longBitsToDouble(Double.doubleToLongBits(t) & 0xffffffff00000000L);
double t2 = t - t1;
double y1 = Double.longBitsToDouble(Double.doubleToLongBits(y) & 0xffffffff00000000L);
double y2 = y - y1;
double x_y = x - y;
z = Math.sqrt(t1*y1 + (x_y*x_y + (t1*y2 + t2*y))); // Note: 2*x*y <= x*x + y*y
}
if (bias == 0) {
return z;
} else {
return Math.scalb(z, bias);
}
}

Answer 3

No es una simple función sqrt. Debe verificar este enlace para la implementación del algoritmo: http://www.koders.com/c/fid7D3C8841ADC384A5F8DE0D081C88331E3909BF3A.aspx

Tiene tiempo de bucle para verificar la convergencia.

/* Slower but safer algorithm due to Moler and Morrison. Never
produces any intermediate result greater than roughly the
larger of X and Y. Should converge to machine-precision
accuracy in 3 iterations. */
double r = ratio*ratio, t, s, p = abig, q = asmall;
do {
t = 4. + r;
if (t == 4.)
break;
s = r / t;
p += 2. * s * p;
q *= s;
r = (q / p) * (q / p);
} while (1);

EDITAR (Actualización de JM):

Here está el documento original de Moler-Morrison, y here hay un buen seguimiento debido a Dubrulle.

Answer 4

http://steve.hollasch.net/cgindex/math/pythag-root.txt

sugiere que la implementación más rápida utilizando sqrt () es cuadrática frente a cúbica para Moler & Morrison, con aproximadamente las mismas características de desbordamiento.