From: Robert Dodier <robert.dodier@gm...>  20060222 20:19:36

Hello, So far as I can tell, SBCL 0.9.9 doesn't conform to CLHS (to judge by the test script appended to this message) for branch cuts of the functions asin, acos, and atanh. I didn't check any other functions. (1) acos and asin: CLHS says (http://www.lispworks.com/documentation/HyperSpec/Body/f_asin_.htm) that the branch cut on [1, \infty) is continuous w/ 4th quadrant and on (\infty, 1] is continuous w/ 2nd quadrant. SBCL: continuous w/ 1st and 2nd quadrants respectively (2) atanh: CLHS says (http://www.lispworks.com/documentation/HyperSpec/Body/f_sinh_.htm) that the branch cut on [1, \infty) is continuous w/ 1st quadrant and on (\infty, 1] is continuous w/ 3rd quadrant. SBCL: continuous w/ 1st and 2nd respectively Thanks for your attention to this topic. Robert Dodier PS. Here's a test script: (cl:asin #C(2.5 0.01)) (cl:asin 2.5) (cl:asin #C(2.5 0.01)) (cl:asin #C(2.5 0.01)) (cl:asin 2.5) (cl:asin #C(2.5 0.01)) (cl:acos #C(2.5 0.01)) (cl:acos 2.5) (cl:acos #C(2.5 0.01)) (cl:acos #C(2.5 0.01)) (cl:acos 2.5) (cl:acos #C(2.5 0.01)) (cl:atanh #C(2.5 0.01)) (cl:atanh 2.5) (cl:atanh #C(2.5 0.01)) (cl:atanh #C(2.5 0.01)) (cl:atanh 2.5) (cl:atanh #C(2.5 0.01)) 
From: Christophe Rhodes <csr21@ca...>  20060222 22:04:29

"Robert Dodier" <robert.dodier@...> writes: > So far as I can tell, SBCL 0.9.9 doesn't conform to CLHS > (to judge by the test script appended to this message) for branch cuts of > the functions asin, acos, and atanh. I didn't check any other functions. > > (1) acos and asin: CLHS says > (http://www.lispworks.com/documentation/HyperSpec/Body/f_asin_.htm) > that the branch cut on [1, \infty) is continuous w/ 4th quadrant > and on (\infty, 1] is continuous w/ 2nd quadrant. > > SBCL: continuous w/ 1st and 2nd quadrants respectively > > (2) atanh: CLHS says > (http://www.lispworks.com/documentation/HyperSpec/Body/f_sinh_.htm) > that the branch cut on [1, \infty) is continuous w/ 1st quadrant > and on (\infty, 1] is continuous w/ 3rd quadrant. > > SBCL: continuous w/ 1st and 2nd respectively > > Thanks for your attention to this topic. This is nasty. I think that the essential problem, despite the best efforts of Steele (and Kahan) is that CL isn't quite right in this area. The question raised by your test script is: what does 2.5 mean in the complex plane? Is it 2.5+0.0i, 2.50.0i or 2.5+<exactly 0>i? Does 2.5 mean 2.5+0.0i, 2.50.0i or 2.5+<exactly 0>i? I suppose that the simple view would be that these real floats are on the axis, but unfortunately that's not what the rest of the spec says; (COMPLEX 2.5) is #c(2.5 0.0), but while (complex 2.5) is #c(2.5 0.0), according to the spec (IMAGPART 2.5) is 0.0. I think that the first thing to do is check whether or not the various trig functions obey the branch cuts when given arguments of type complex: that is, with an explicit positive or negative zero imaginary part. After that, we can discuss what the most useful interpretation of the real numbers in terms of the branch cuts are. (The reason I say CL doesn't get this quite right is that there's no pure imaginary floating point quantity, so composition of branch cuts can't work right no matter what choices are made.) You might also be interested in the ieeefptests project on commonlisp.net. Cheers, Christophe 
From: Robert Dodier <robert.dodier@gm...>  20060223 17:14:38

Hello Harald, I'm investigating some stuff, in the meantime here is a quick reply. > In any case, I think anybody thinking of rewriting this stuff for SBCL > ought to have Kahan's paper as a handy reference. Heck, I might want > to have a look myself, only just not this week. I'm aware of Kahan's work on this topic. However I'm not sure it's relevant for resolving questions about CL. Let me back up a little. I am working on the Maxima algebra system. Maxima attempts to be supported (i.e. compiles and passes its test suite) on several CL implementations. From Maxima's POV, we need to determine what is unambiguously specified by CLHS and what is not. For the former, Maxima can assume behavior according to spec and punt on test suite errors ("it's somebody else's problem") or work around them. For the latter, Maxima has to either propagate the ambiguity into its own spec, or arrange to call functions only when the result is unambiguous. There are various nonnormative documents auxilliary to CLHS (including Kahan, CLtL2, and the "clean up" documents). But since claimedtobeconforming CL implementations aren't required to pay attention to them, Maxima can't count on anything derived from them. I'm assuming here that CLHS is the sole spec for CL. If that's not the case you'll let me know. For what it's worth, Robert Dodier 
From: Harald HancheOlsen <hanche@ma...>  20060223 17:34:17

+ "Robert Dodier" <robert.dodier@...>:  I'm aware of Kahan's work on this topic.  However I'm not sure it's relevant for resolving questions about CL. I'm sorry  I was too terse, didn't have time to write more clearly. I didn't mean to use Kahan to resolve problems with the CLHS. But he has useful things to say about implementation issues assuming (as I believe to be the case) that his ideas of branch cuts coincides with the HyperSpec's. ... oh, but now when I browse the CLHS I do see there are some issues  and references to Kahan, too  I should probably shut my mouth until I have had a closer look. It's been too long ...  I'm assuming here that CLHS is the sole spec for CL. Indeed, and I had not intended to say otherwise. (Except, strictly speaking, I suppose the ANSI spec is THE authority.)  Harald 
From: Robert Dodier <robert.dodier@gm...>  20060227 14:07:10

On 2/22/06, Christophe Rhodes <csr21@...> wrote: > The question raised by your test script is: what does 2.5 mean in the > complex plane? Is it 2.5+0.0i, 2.50.0i or 2.5+<exactly 0>i? Does > 2.5 mean 2.5+0.0i, 2.50.0i or 2.5+<exactly 0>i? I suppose that > the simple view would be that these real floats are on the axis, but > unfortunately that's not what the rest of the spec says; (COMPLEX 2.5) > is #c(2.5 0.0), but while (complex 2.5) is #c(2.5 0.0), according to > the spec (IMAGPART 2.5) is 0.0. Well, the discrepancy between (imagpart (complex 2.5)) and (imagpart 2.5) is certainly disheartening ... From a broader perspective, that the CLHS allows implementations to return different results depending on the presence or absence of signed zero makes it very burdensome for a programmer to write portable code. Too bad about that I guess. > I think that the first thing to do is check whether or not the various > trig functions obey the branch cuts when given arguments of type > complex: that is, with an explicit positive or negative zero imaginary > part. After that, we can discuss what the most useful interpretation > of the real numbers in terms of the branch cuts are. OK, I've studied this problem some more and here's what I've come to. I believe the CLHS is OK wrt asin and acos. About atanh, the CLHS specifies behavior according to the formula atanh x =3D (log(1 + x)  log(1  x))/2 and then goes on to give a description of branch cuts which, to the best of my knowledge, is inconsistent with the formula. I'm assuming the formula is the spec and the stuff about branch cuts is just commentary. So (I've concluded) the branch cuts for asin, acos, and atanh should all (\infty, 1) continuous w/ 2nd quadrant, and (1, \infty) continuous w/ 4th quadrant. For arguments which are complex numbers, it appears that the SBCL asin, acos, and atanh have branch cuts which are continuous with 2nd and 1st quadrants, respectively. See the PS below. For the record, other CL implementations that I know about (GCL, Clisp, CMUCL) make differing choices of branch cuts, and no two of the four make the same choices on all three functions asin, acos, and atanh. For what it's worth, Robert Dodier PS. test script: (format t "~S~%~S~%~S~%~%" (cl:asin #C(2.5 0.01)) (cl:asin #C(2.5 0.0)) (cl:asin #C(2.5 0.01))) (format t "~S~%~S~%~S~%~%" (cl:asin #C(2.5 0.01)) (cl:asin #C(2.5 0.0)) (cl:asin #C(2.5 0.01))) (format t "~S~%~S~%~S~%~%" (cl:acos #C(2.5 0.01)) (cl:acos #C(2.5 0.0)) (cl:acos #C(2.5 0.01))) (format t "~S~%~S~%~S~%~%" (cl:acos #C(2.5 0.01)) (cl:acos #C(2.5 0.0)) (cl:acos #C(2.5 0.01))) (format t "~S~%~S~%~S~%~%" (cl:atanh #C(2.5 0.01)) (cl:atanh #C(2.5 0.0)) (cl:atanh #C(2.5 0.01))) (format t "~S~%~S~%~S~%~%" (cl:atanh #C(2.5 0.01)) (cl:atanh #C(2.5 0.0)) (cl:atanh #C(2.5 0.01))) output: #C(1.566432 1.5668098) #C(1.5707964 1.5667993) #C(1.566432 1.5668098) #C(1.566432 1.5668098) #C(1.5707964 1.5667993) #C(1.566432 1.5668098) #C(0.004364322 1.5668098) #C(0.0 1.5667993) #C(0.004364322 1.5668098) #C(3.1372283 1.5668098) #C(3.1415927 1.5667993) #C(3.1372283 1.5668098) #C(0.42363986 1.5688916) #C(0.42364892 1.5707964) #C(0.42363986 1.5688916) #C(0.42363986 1.5688916) #C(0.42364892 1.5707964) #C(0.42363986 1.5688916) 
From: Harald HancheOlsen <hanche@ma...>  20060223 07:28:22

+ Christophe Rhodes <csr21@...>:  This is nasty. I think that the essential problem, despite the best  efforts of Steele (and Kahan) is that CL isn't quite right in this  area. Some years ago, I was looking into similar branch cut problems in python (the language, not the cmucl compiler). (I was in fact planning to rewrite pythons cmathmodule to get it right, but then I never got around to it. But I digress.) At the time, I was looking at a paper by Kahan [1] as well as the HyperSpec. If I remember correctly, I concluded that the HyperSpec and Kahan were in agreement, and thinking that python ought to do the same thing. In any case, I think anybody thinking of rewriting this stuff for SBCL ought to have Kahan's paper as a handy reference. Heck, I might want to have a look myself, only just not this week.  Harald [1] W. Kahan: Branch Cuts for Complex Elementary Functions, or Much Ado About Nothing's Sign Bit. In A. Iserles, M. J. D. Powell (eds.): Proceedings of the Joint IMA/SIAM Conference on The State of the Art in Numerical Analysis held at the University of Birmingham, 1418 April 1986. Clarendon Press, Oxford 1987. ISBN 0198536143. 