[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: 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.


[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

Attachment: gsl.patch
Description: Text Data

reply via email to

[Prev in Thread] Current Thread [Next in Thread]