Phil Steitz wrote: > Luc Maisonobe wrote: >> Phil Steitz a écrit : >> >>> Ted Dunning wrote: >>> >>>> I just check on how sbcl handles this. The result was mixed, but >>>> instructive as predicted. >>>> >>>> Positive and negative infinity behave as expected. >>>> >>>> ** (/ 1 0.0) >>>> #.SB-EXT:SINGLE-FLOAT-POSITIVE-INFINITY >>>> >>>> * (/ -1 0.0) >>>> * >>>> ***#.SB-EXT:SINGLE-FLOAT-NEGATIVE-INFINITY >>>> >>>> * >>>> Square roots of negative infinity work as for option 3. >>>> >>>> ** (sqrt (/ -1 .0)) >>>> #C(0.0 #.SB-EXT:SINGLE-FLOAT-POSITIVE-INFINITY) >>>> >>>> * >>>> Fourth roots are reasonable as well. >>>> ** (sqrt (sqrt (/ -1 0.0))) >>>> #C(#.SB-EXT:SINGLE-FLOAT-POSITIVE-INFINITY >>>> #.SB-EXT:SINGLE-FLOAT-POSITIVE-INFINITY) >>>> * >>>> But, the exponentiation operator is inconsistent with these results >>>> and >>>> behaves like option (1) (without documentation as far as I can tell). >>>> >>>> ** (expt (/ -1 .0) 0.25) >>>> #.SB-EXT:SINGLE-FLOAT-POSITIVE-INFINITY >>>> >>>> * (expt (/ -1 .0) 0.5) >>>> #.SB-EXT:SINGLE-FLOAT-POSITIVE-INFINITY >>>> >>> Thanks, Ted. This is instructive. R also handles square roots in a >>> sort of 3)-like way, but inconsistently and without maintaining z = >>> sqrt(z) * sqrt(z). Arg works as expected, consistent with atan2 (what >>> we use now); but Arg(sqrt(z)) = Arg(z)/2 fails in the lower half-plane. >>> >>> >>>> z <- complex(real=Inf, imaginary=0) >>>> Arg(z) >>>> >>> [1] 0 >>> >>>> sqrt(z) >>>> >>> [1] Inf+0i >>> >>>> Arg(sqrt(z)) >>>> >>> [1] 0 >>> >>>> sqrt(z)*sqrt(z) >>>> >>> [1] Inf+NaNi >>> >>>> z <- complex(real=0, imaginary=Inf) >>>> Arg(z) >>>> >>> [1] 1.570796 >>> >>>> sqrt(z) >>>> >>> [1] Inf+Infi >>> >>>> Arg(sqrt(z)) >>>> >>> [1] 0.7853982 >>> >>>> sqrt(z)*sqrt(z) >>>> >>> [1] NaN+Infi >>> >>>> z <- complex(real=-Inf, imaginary=Inf) >>>> sqrt(z) >>>> >>> [1] Inf+Infi >>> >>>> Arg(z) >>>> >>> [1] 2.356194 >>> >>>> Arg(sqrt(z)) >>>> >>> [1] 0.7853982 >>> >>>> sqrt(z)*sqrt(z) >>>> >>> [1] NaN+Infi >>> >>>> z <- complex(real=-Inf, imaginary=0) >>>> Arg(z) >>>> >>> [1] 3.141593 >>> >>>> sqrt(z) >>>> >>> [1] 0+Infi >>> >>>> Arg(sqrt(z)) >>>> >>> [1] 1.570796 >>> >>>> sqrt(z)*sqrt(z) >>>> >>> [1] -Inf+NaNi >>> >>>> z <- complex(real=-Inf, imaginary=-Inf) >>>> Arg(z) >>>> >>> [1] -2.356194 >>> >>>> sqrt(z) >>>> >>> [1] Inf-Infi >>> >>>> Arg(sqrt(z)) >>>> >>> [1] -0.7853982 >>> >>>> sqrt(z)*sqrt(z) >>>> >>> [1] NaN-Infi >>> >>> The problem with 3) is that it unfortunately *only* works for sqrt and >>> 4th roots. Doh! The problem is that pi/4 splits are the only args >>> that >>> you can "reasonably" represent with real and imaginary parts from {0, >>> +-Inf}. The implementation of atan2 collapses all args for pairs with >>> an infinite element into these, which makes sense from a limit >>> perspective. >>> >>> So I think we are left with either 1) or 2) or finding a "standard" >>> impl >>> or spec to imitate. Unfortunatly, the C99x standard is not public >>> and I >>> don't recall complex nth roots being included there in any case. If >>> anyone knows of a reasonable public standard to implement here, that >>> would be great to at least look at. >>> >> >> Solution 2 seems a good compromise to me. It is reasonably logical and >> could probably be implemented simply. Solution 1 looks clumsy and we >> could probably not withstand it on the long run, which would imply >> breaking compatibility some day to go back to a more definitive >> solution. >> >> >>> An improvement to 1), lets call it 1.5), that would make behavior a >>> little less strange would be to avoid the NaNs from Inf * 0 and the >>> "spurious" Infs from Inf * cos(phi) or Inf * sin(phi) where rounding >>> error has moved phi off an analytical 0. You can see examples of both >>> of these in first two infinite test cases in ComplexTest.java as it is >>> today. If the code was changed to work around these cases, the >>> square >>> roots of Inf + 0i, for example, would be {Inf + 0i, -Inf + 0i}. >>> I have not tested this, but it looks to me like the NaN in the first >>> root is from Inf * 0 and the -Inf in the imaginary part of the second >>> root is from Inf * sin(phi) with phi from 2 * pi / 2 just slightly >>> greater than pi. >>> >> >> Is this simpler or more difficult to implement than 2 ? If simpler, then >> this could be a reasonable solution too. >> > 2) is easiest, since all you have to do is check isInfinite(). The > 1.5) idea is more complicated and unfortunately the "avoid spurious > infs from args near analytical zeros" bit requires an arbitrary choice > of what "near" means and also still does not guarantee that the roots > solve the definitional equation. So my vote is for 2). If there are > no objections, I will implement that with the other changes. > > Phil Committed in r731822. I noticed one more thing that I just commented for now in the code. It might be faster and/or more accurate to use a solver to compute the nth root of the modulus. Should we consider doing this and exposing the associated real function in MathUtils? If yes, what solver with what parameters? Phil >> Luc >> >> >>> Phil >>> >>> >>>> * >>>> >>>> On Fri, Jan 2, 2009 at 8:19 PM, Ted Dunning >>>> wrote: >>>> >>>> >>>> >>>>> I think that (2) is the easiest for a user to understand. Obviously, >>>>> the >>>>> documentation aspect of (1) should be brought along. >>>>> >>>>> The behavior of Common Lisp is always instructive in these >>>>> regards. The >>>>> complex arithmetic was generally very well thought out. >>>>> >>>>> On Fri, Jan 2, 2009 at 7:14 PM, Phil Steitz >>>>> wrote: >>>>> >>>>> >>>>>> ... I noticed another thing. Behavior for complex numbers with >>>>>> NaN and >>>>>> >>>>>>> infinite parts is, well, "interesting." It is consistent with what >>>>>>> we do >>>>>>> elsewhere in this class to just apply computational formulas and >>>>>>> document >>>>>>> behavior for infinite and NaN; but the current implementation is >>>>>>> hard to >>>>>>> document correctly and the results are a little strange for >>>>>>> infinite values. >>>>>>> >>>>>>> I see three reasonable choices here, and am interested in others' >>>>>>> opinions. Please select from the following or suggest something >>>>>>> else. >>>>>>> >>>>>>> 1) Leave as is and just point out that the computational formula >>>>>>> will >>>>>>> behave strangely for infinite values. >>>>>>> >>>>>>> 2) return {Complex.NaN} or {Complex.INF} respectively when the >>>>>>> number has >>>>>>> a NaN or infinite part. >>>>>>> >>>>>>> 3) return {Complex.NaN} when either part is NaN; but for infinite >>>>>>> values, >>>>>>> compute the argument using getArgument (atan2), generate the >>>>>>> arguments for >>>>>>> the roots from this and select the real/im parts of the roots from >>>>>>> {0, inf, >>>>>>> -inf} to match the argument (this will always be possible because >>>>>>> atan2 >>>>>>> always returns a multiple of pi/4 for infinite arguments). For >>>>>>> example, the >>>>>>> 4th roots of real positive infinity would be {inf + 0i, 0 + infi, >>>>>>> -inf + 0i, >>>>>>> 0 + -infi} >>>>>>> >>>>>>> 2) is probably the most defensible mathematically; but 3) is >>>>>>> closer to >>>>>>> the spirit of C99x. Unfortunately, since our implementation of >>>>>>> multiply is >>>>>>> 2)-like, 3) does not guarantee that nth roots actually solve the >>>>>>> equation >>>>>>> r^n = z. >>>>>>> >>>>>>> >>>> >>> --------------------------------------------------------------------- >>> To unsubscribe, e-mail: dev-unsubscribe@commons.apache.org >>> For additional commands, e-mail: dev-help@commons.apache.org >>> >>> >> >> >> --------------------------------------------------------------------- >> To unsubscribe, e-mail: dev-unsubscribe@commons.apache.org >> For additional commands, e-mail: dev-help@commons.apache.org >> >> > --------------------------------------------------------------------- To unsubscribe, e-mail: dev-unsubscribe@commons.apache.org For additional commands, e-mail: dev-help@commons.apache.org