Is Abs[z]^2 a bad way to calculate the square modulus of z?

Question

For a numerical quantity z, Abs[z] returns the square root of the sum of the squares of the real and imaginary parts of z. If we then square it, do we lose numerical accuracy because we are squaring after taking a square root?
Is it better to use Re[z]^2+Im[z]^2?

Answer 1

Let's find out.

Here is a function to compare the results of evaluating f on a machine approximation of an exact argument z to its exact value:

test[f_] := Function[{z}, f[N[z]] -  f[z]];

We expect to find problems around the 52nd significant binary digit or so. Let's magnify these differences, then, by about $2^{60}$ and use logarithmic axes to study them. (In order to look at the logs, we neglect the signs of the errors and also add a small starting value to them to avoid taking the logarithm of $0$. Even so, the y-axes more or less indicate the number of binary digits lost out of $60$. Thus, we should expect log errors to be typically around $8$ when z has a modulus near $1$.) The errors will likely vary with the argument of the complex number z and will tend to be directly proportional to its modulus.

The following code--which is readily extensible to alternative procedures by modifying its first line--explores three ways of evaluating a squared modulus, using Norm, Abs, and the sum of squares of real and imaginary parts. To do this, it plots the errors at equally spaced points along semi circles of radii $3/4$, $1$, and $4/3$ (colored blue, red, and gold, respectively).

functions = {Norm[#]^2 &, Abs[#]^2 &, (Re[#]^2 + Im[#]^2) &};
TableForm[
 {Table[ListPlot[
     Table[Log2[1 + 2^60 Abs[test[g][\[Rho] Exp[I Pi \[Theta]]]]], 
         {\[Rho], {3/4, 1, 4/3}}, {\[Theta], 0, 1, 1/Prime[80]}], 
    PlotStyle -> PointSize[0.01], DataRange -> {0, Pi}, PlotRange -> {Full, {-1, 11}}, 
    ImageSize -> 300], 
   {g, functions}]},
 TableHeadings -> {{}, Style[TraditionalForm[#], 16] & /@ functions}, 
 TableAlignments -> Center]

As expected, the error tends to be directly proportional to the modulus of the argument, as indicated by the approximately equal spacings between the three colors.

When it errs, Abs (shown in the middle column) tends to make relatively large errors, but it appears to err less often than the other two methods, which have similar (but not identical) error patterns. Moreover, the other two procedures have strong angular patterns in their errors: they tend to err more frequently for directions oriented along the axes, but may have slightly smaller errors there.

Abs probably uses a square root function that is optimized for accuracy: it works most of the time, but not all the time; and when it fails, it's just as bad--and sometimes slightly worse--than the alternatives.

Abs[#]^2 & behaves much like # Conjugate[#] & (not shown). This is a puzzle: the latter clearly does not optimize a square root, because it never uses one. This suggests that in the multiplication, Mathematica is exploiting the polar representation used to create points along the circles. Tracing the evaluation confirms this. Thus we learn that the accuracy of a calculation may depend not only on the value of the argument, but also on its form, even when the argument is only going to be evaluated to machine precision.