Menu ▾ ▴

#5281 Float atan and acot take the other side of their branch cuts from bigfloat and exact evaluation

None
closed
nobody
5
5 days ago
6 days ago
No

By Claude:

On part of the imaginary axis, float atan and acot differ by %pi from bfloat and rectform of the same value: atan below -%i, and acot on its whole cut between -%i and %i. The other ten inverse trigonometric and hyperbolic functions agree on their cuts.

(%i1) [atan(-2.0*%i), rectform(atan(-2*%i))];
(%o1) [1.5707963267948966-0.5493061443340549*%i,-((%i*log(3))/2)-%pi/2]
(%i2) [acot(0.5*%i), rectform(acot(%i/2))];
(%o2) [1.5707963267948966-0.5493061443340549*%i,-((%i*log(3))/2)-%pi/2]

The exact and bigfloat values keep atan odd and equal to -%i*atanh(%i*x), which is how the simplifier rewrites atan(-2*%i). The float value does neither: atan(2.0*%i) is 1.5708+0.5493*%i, so float atan takes the right half-plane's side on both halves of the cut. As a result a simplified form and its float evaluation disagree: cos(atan(-2.0*%i)) gives 0.5774*%i, but subst(x = -2.0*%i, cos(atan(x))) gives -0.5774*%i.

Discussion

  • Raymond Toy

    Raymond Toy - 6 days ago

    Using a bulid from Sep 22, I get:

    (%i3) atan(-2.0*%i);
    (%o3)             1.5707963267948966 - 0.5493061443340549 %i
    

    But clisp, cmucl, and sbcl say (atan #c(0 -2.0)) is #C(1.5707964 -0.54930615). It's going to be really confusing if we make Maxima return different values.

     
    • David Scherfgen

      David Scherfgen - 6 days ago

      Regarding CL, I think it's down to signed zero. (atan #c(0d0 -2d0)) gives #C(1.5707963267948966d0 -0.5493061443340549d0), while (atan (complex -0d0 -2d0)) gives #C(-1.5707963267948966d0 -0.5493061443340549d0). The value depends on which side the zero says we came from. A Maxima float like -2.0*%i has no signed zero, and passing it on as #c(0.0 -2.0) quietly picks the right half-plane.

      Maxima already declines to follow the Lisp in the analogous case. SBCL's (atanh 2d0) is 0.549+1.571i, but Maxima's atanh(2.0) is 0.549-1.571i, because maxima-branch-atanh computes cut values itself ("Some Lisp implementations goof up branch cuts for ASIN, ACOS, and/or ATANH..."). atan is the one that goes straight to cl:atan, so Maxima's own floats now contradict each other:

      (%i1) [atan(-2.0*%i), expand(-%i*atanh(2.0))];
      (%o1) [1.5707963267948966-0.5493061443340549*%i,-(0.5493061443340549*%i)-1.5707963267948966]
      

      Everything else in Maxima gives the second value: rectform, bfloat, the logarithmic form from logarc, and the simplifier's own rewrite atan(-2*%i) → -%i*atanh(2).

      A maxima-branch-atan could fix this the same way maxima-branch-atanh does. It would use -%i*maxima-branch-atanh(%i*x) exactly on the cut (real part 0.0, abs(imagpart) > 1) and cl:atan elsewhere; acot, which goes through atan(1/x), would follow. The cost is differing from (atan #c(0.0 -2.0)) on the lower half of the cut.

      Keeping the CL value would mean changing the exact and bigfloat side instead. That gives up atan(-x) = -atan(x) on the cut, which the simplifier's reflection rule relies on everywhere.

       
      • Raymond Toy

        Raymond Toy - 6 days ago

        clisp doesn't have signed zeros and returns #c(1.57 -0.549). cmucl returns a different value for atanh(2): #C(0.54930615 -1.5707964). This is because long ago I carefully derived the value from first principals and using Kahan's convention of not spuriously adding a signed real part unless absolutely necessary. (Or maybe that was for acot? I forget now.)

        Anyway, dealing with branch cuts is quite complicated with or without signed-zeroes, and I've generally followed Kahan's conventions as shown in his paper "Much Ado about (nothing's) sign". I think his paper says that many expected identities won't hold on the branch cuts. Or if they do, other things break.

         
  • David Scherfgen

    David Scherfgen - 5 days ago
    • status: open --> closed
     
  • David Scherfgen

    David Scherfgen - 5 days ago

    Fixed by commit [62d4a5].
    After a few further minor bugfixes that I'll commit next, all the inverse trigonometric/hyperbolic functions agree on 5 paths: simplifier (exact values, reflection rules), float evaluation, bigfloat evaluation, rectform and logarc.

     

    Related

    Commit: [62d4a5]


Log in to post a comment.