[Top][All Lists]
[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]
[Bug-gsl] GSL fails to honor zero's sign at branch cuts
From: |
Richard B. Kreckel |
Subject: |
[Bug-gsl] GSL fails to honor zero's sign at branch cuts |
Date: |
Sun, 13 Nov 2011 00:23:32 +0100 |
User-agent: |
Mozilla/5.0 (X11; U; Linux x86_64; en-US; rv:1.9.2.20) Gecko/20110826 Icedove/3.1.12 |
Systems that support signed zero should be careful about the value on
the branch cuts ("slits"). To quote William Kahan [0]:
"Since all the slits in question lie on either the real or the imaginary
axis, every point z on a slit is represented in at least two ways, at
least once with a +0 and at least once with a -0 for whichever of the
real and imaginary parts of z vanishes. Benignly, ambiguity in z at a
discontinuity of f(z) permits f(z) to be defined formally continuously,
except possibly at the ends of some slits, by continuation from inside
the principal domain. This continuity goes beyond mere formality. By
analytic continuation, the domain of each of our nine elementary
functions f(z) extends until it fills out a Riemann Surface; think of
this surface as a multiple covering wrapped like a bandage around the
Riemann Sphere and mapped onto it continuously by f. To construct f's
principal domain, cut the bandage along the slit(s) and discard all but
one layer covering the sphere. That layer is a closed surface mapped by
f continuously onto a subset of the sphere. The shadow of that layer
projected down upon the sphere is the principal domain; it consists of
the whole sphere, but with slit(s) covered twice. That is why we wish to
represent slits ambiguously."
The basic convention (adopted in C99, CLTL2, etc. and implemented in
GCC, GLibC, CLisp, etc.) for branch cuts with only one finite endpoint
is counter-clockwise continuity (CCC): "the function is continuous as
the cut is approached coming around the finite endpoint of the cut in a
counter clock wise direction" (wording from [1]).
GSL falls into the category of systems that should do this, since it
uses IEEE 754 double precision for constructing complex numbers.
However, GSL does not distinguish between +0 and -0 on the branch cuts!
I propose to fix that for two reasons:
1) Distinguishing between +0 and -0 is a superior concept.
2) Fixing this harmonizes GSL with other systems that support both
signed zero and complex numbers.
I attach a patch into that direction. It fixes sqrt, arcsin, arccos,
arctan, arcsinh, arccosh, and arctanh.
With that patch, some checks fail. Where they fail it should be
sufficient to flip the sign of the zero real or imaginary part in the
argument. (Of course, if one wishes, one can add new test vectors for
the opposite signed zero.)
Thanks for applying.
-richy.
[0] William Kahan: "Branch Cuts for Complex Elementary Functions, or
Much Ado About Nothing's Sign Bit"
[1] ISO/IEC 9899:1999 (E): "Programming Languages: C"
--
Richard B. Kreckel
<http://www.ginac.de/~kreckel/>
gsl.patch
Description: Text Data
- [Bug-gsl] GSL fails to honor zero's sign at branch cuts,
Richard B. Kreckel <=