Robert wrote:Josef Templ wrote:Very nice example. To me it looks like NaN is not the problem here, right?
In fact the problem above is actually caused by the '+' function. I'll post an explanation later.
I'm writing this mainly for my own benefit.
What are these mysterious signed zeros?
I've seen explanations that suggest they are used to
approximate non-zero numbers too small to otherwise represent. I don't like that explanation; numbers are just numbers, with neither resolution or accuracy. But, of course, computer arithmetic is different to mathematical arithmetic as in includes rounding and other artifacts.
Number
formats have an associated resolution. Numbers, may, in a given context have an associated accuracy specific to that context. The numbers INTEGER 3, SHORTREAL 3.0, and REAL 3.0 all are the same; they don't have different accuracies, three is simply three.
The way to think about +0.0 & −0.0 are as packages that both contain the number 0.0, and an additional bit of information (or dis-information) that the application programmer may, optionally, exploit.
What are they used for?
When we use a function we are often also interested in closely related functions such as its inverse or reciprocal. Sometimes one of these functions has a discontinuity at x = 0. The limit of F(x) as x tends to 0 from below may exist, call it F−. The limit of F(x) as x tends to 0 from above may also exist, call it F+. How do we define F(0)? The mathematician has to make a choice; common choices are to leave F(0) undefined, define it as F−, define it as F+, define it as (F− + F+) / 2, or ...
But the programmer has a further choice; their subroutine can return F− when asked for the limit from below, and return F+ when asked for the limit from above. These requests can be signalled by calling F(x) with x = −0 or +0 respectively.
A simple example; 1 / (−0) returns −INF as that is the limit of 1 / x as x tends up to 0, and the same for 1 / (+0) as x tends to 0 from above.
−0 & 0 are not different numbers, they are ways of asking for different limits at the same point, namely at zero.
Examples of this process with real functions seem to be rare; it is much more common with complex functions. But, in programming, complex functions are built with real functions, so it has to be built into the real arithmetic system.
A (simple) complex example
Consider squaring and then square-rooting z, with z somewhere in the right half plane (ie with positive real part). (Complex) squaring involves squaring the magnitude and doubling the argument (angle) of z; square-rooting involves square-rooting the magnitude and halving the angle; so it brings us back to z.
But for z on the imaginary axis doubling the angle gets us onto the negative real axis, and essentially we do not know if the square-root of −1 is +i or −i. But if we know we are approaching the axis from above (real part of z is +0) continuity says we want the answer +i, and if we are approaching from below (real part of z is −0) continuity says we want the answer −i. This enables us to extend the squaring & square-root functions to include the imaginary axis, not just the strictly right half plane. This is both useful and convenient as the Borda example shows.
Radial flow - pre-Borda
In this picture imagine water in the right half plane, air in the left, both under some pressure, and a tiny hole at the origin. Both fluids flow radially into this sink hole - the water flow lines are shown. Now insert two barriers (shown in blue) with air on the left and water on the right. They do not affect the flow.
(This picture is easy to plot!)−
![Radial flow](./download/file.php?id=226&sid=fad47fd180a11fcbee1621fdf49a9872)
- Radial flow
- PreBorda.png (4.69 KiB) Viewed 8775 times
Flow lines, Borda's mouthpiece
Distort this picture by replacing z by
F(z) ≡ 1 + g(z) + Ln (g(z)) where g(z) ≡ z² + z √(1 + z²).
![Borda's mouthpiece](./download/file.php?id=227&sid=fad47fd180a11fcbee1621fdf49a9872)
- Borda's mouthpiece
- Borda.png (5.25 KiB) Viewed 8775 times
This rotates the barriers to be horizontal (the edges of a pipe extending infinitely to the left), and moves the sink hole also to ( −∞ + 0 i).
Since F is a
conformal map (angle preserving) the distorted picture also represents flow lines, now the flow coming out of the pipe into an infinite reservoir.
It is extremely convenient to draw not just the flow lines strictly within the flow, but also those defining the edges. These are the vertical grey lines in the first picture, those on the imaginary axis. These require √(z²) to be computed correctly for both z = +i & z = −i, as in the
simple complex example. Thus signed zeros are used to indicate the direction we need the √ function to be continuous.
Why do most libraries, and some languages, fail for this example?
The explanation lies in about the simplest part of G, namely the expression √(1 + z²).
Consider, as examples, the red and green crosses in the pictures.
For them we (should) have (red then green)
Code: Select all
z: 0 − 2 i 0 + 2 i
z²: − 4 − 0 i − 4 + 0 i
(1 + z²): − 3 − 0 i − 3 + 0 i
√(1 + z²). 0 − 1.73205 i 0 + 1.73205 i
However many libraries cannot add the real number 1! Instead they type cast it into the complex number (1 + 0 i) first, then do a complex add. This is not only slower, it is
incorrect. Since 0 + (−0) = 0 + (+0) = 0 both spots have (1 + z²) computed as − 4 + 0 i, and the flow lines adjacent to the bottom pipe are misplotted on the top pipe.
Why does only zero have this continuity hint?
1 - The most common elementary functions in sub-routine libraries are things like Ln, ArcSin, ArcCos, ArcTan, ArcCsc, ArcSec, ArcCot, ArcSinh, ArcCosh, ArcTanh, ArcSch, ArcSech, ArcCoth, and these
all have their branch cuts (places of discontinuity) defined on the complex axes. So the direction of continuity requested can be signalled by the sign of a zero real part, or the sign of a zero imaginary part.
Functions with branch cuts elsewhere are simply
very rare.
2 - If one wished to exploit this technique on a function H (x) that had is discontiuity at, say, 5 one could simply define it in terms of another function H'(x) ≡ H(x − 5), and now H' does have its discontinuity at zero; problem solved.