Is it safe to cascade hypot()?

Question

I would like to compute the norm (length) of three- and four-dimensional vectors. I'm using double-precision floating point numbers and want to be careful to avoid unnecessary overflow or underflow.

The C math library provides hypot(x,y) for computing the norm of two-dimensional vectors, being careful to avoid underflow/overflow in intermediate calculations.

My question: Is it safe to use hypot(x, hypot(y, z)) and hypot(hypot(w, x), hypot(y, z)) to compute the lengths of three- and four-dimensional vectors, respectively?

Answer 1

It's safe, but it's wasteful: you only need to compute sqrt() once, but when you cascade hypot(), you will call sqrt() for every call to hypot(). Ordinarily I might not be concerned about the performance, but this may also degrade the precision of the result. You could write your own:

double hypot3(double x, double y, double z) {
    return sqrt(x*x + y*y + z*z);
}

etc. This will be faster and more accurate. I don't think anyone would be confused when they see hypot3() in your code.

The standard library hypot() may have tricks to avoid overflow, but you may not be concerned about it. Ordinarily, hypot() is more accurate than sqrt(x*x + y*y). See e_hypot.c in the GLibC source code.

Answer 2

It safe (almost) to use hypot(x, hypot(y, z)) and hypot(hypot(w, x), hypot(y, z)) to compute the lengths of three- and four-dimensional vectors.

C does not strongly specify that hypot() must work for a double x, y that have a finite double answer. It has weasel words of "without undue overflow or underflow".

Yet given that hypot(x, y) works, a reasonable hypot() implementation will perform hypot(hypot(w, x), hypot(y, z)) as needed. There is only 1 increment (at the low end) /decrement (at the high end) of binary exponent range lost when with 4-D vs. 2-D.

Concerning speed, precision, and range, code profile against sqrtl((long double) w*w + (long double) x*x + (long double) y*y + (long double) z*z) as an alternative, but that seems only needed with select coding goals.