Opened 7 years ago

Closed 7 years ago

#16439 closed defect (fixed)

Allow inverse trig functions to all take complex input

Reported by: kcrisman Owned by:
Priority: minor Milestone: sage-6.3
Component: symbolics Keywords:
Cc: burcin Merged in:
Authors: Punarbasu Purkayastha Reviewers: Karl-Dieter Crisman, Ralf Stephan
Report Upstream: N/A Work issues:
Branch: 7c4c725 (Commits) Commit: 7c4c72596fe2a6de32e61ca578f2e671a51d7962
Dependencies: Stopgaps:

Description

Where by complex I mean complex Python type.

sage: arccos(1+2*i)
arccos(2*I + 1)
sage: arccos(CC(1,2))
1.14371774040242 - 1.52857091948100*I
sage: arccos(CDF(1,2))
1.1437177404 - 1.52857091948*I
sage: arccos(complex(1,2))
<snip>
TypeError: Unable to convert x (='(1+2j)') to real number.

Presumably there are others of this type. Related:

sage: arcsec(-0.992725195401)
NaN
sage: arcsec(float(-0.992725195401))
<snip>
    732         if parent is float:
    733             return math.acos(1/x)
--> 734         return (1/x).arccos()
    735 
    736     def _eval_numpy_(self, x):

AttributeError: 'complex' object has no attribute 'arccos'

but this may or may not be dealt with at #13246.

Change History (28)

comment:1 follow-up: Changed 7 years ago by ppurka

Why don't complex numbers have these attributes? This should be possible - for example wolfram alpha gives a value for acos(-0.9+i).

I think the computations should be performed using CDF and then coerced back to complex. Complex output from these functions should be handled by the plot function separately.

comment:2 in reply to: ↑ 1 ; follow-up: Changed 7 years ago by kcrisman

Why don't complex numbers have these attributes? This should be possible - for example wolfram alpha gives a value for acos(-0.9+i).

Complex number do have these attributes, just not the minimalist complex Python data type.

sage: a = complex(1,1)
sage: a.[tab]
a.conjugate  a.imag       a.real       
sage: a._[tab]
a.__abs__           a.__hash__          a.__rdivmod__
a.__add__           a.__init__          a.__reduce__
a.__class__         a.__int__           a.__reduce_ex__
a.__coerce__        a.__le__            a.__repr__
a.__delattr__       a.__long__          a.__rfloordiv__
a.__div__           a.__lt__            a.__rmod__
a.__divmod__        a.__mod__           a.__rmul__
a.__doc__           a.__mul__           a.__rpow__
a.__eq__            a.__ne__            a.__rsub__
a.__float__         a.__neg__           a.__rtruediv__
a.__floordiv__      a.__new__           a.__setattr__
a.__format__        a.__nonzero__       a.__sizeof__
a.__ge__            a.__pos__           a.__str__
a.__getattribute__  a.__pow__           a.__sub__
a.__getnewargs__    a.__radd__          a.__subclasshook__
a.__gt__            a.__rdiv__          a.__truediv__

And the built-in math module doesn't do arcsec

sage: math.a[tab]
math.acos   math.asin   math.atan   math.atanh  
math.acosh  math.asinh  math.atan2  

and even then only takes ones with real output

sage: math.acos(1)
0.0
sage: math.acos(2)
ValueError: math domain error
sage: math.acos(a)
TypeError: can't convert complex to float

I think the computations should be performed using CDF and then coerced back to complex. Complex output from these functions should be handled by the plot function separately.

Maybe. Is CDF basically the same as complex? Is there precedent for this?

comment:3 in reply to: ↑ 2 Changed 7 years ago by ppurka

Replying to kcrisman:

I think the computations should be performed using CDF and then coerced back to complex. Complex output from these functions should be handled by the plot function separately.

Maybe. Is CDF basically the same as complex? Is there precedent for this?

I do not know a precedent, but we will not be losing anything from doing this. The Python complex number comprises of two floating point numbers which are implemented as double in C - see this documentation. The CDF in Sage also use double precision floating point numbers.

sage: CDF.precision()
53

# Following the documentation at
# https://docs.python.org/2/library/sys.html#sys.float_info
sage: import sys
sage: sys.float_info.mant_dig
53

comment:4 Changed 7 years ago by ppurka

  • Branch set to u/ppurka/allow_inverse_trig_functions_to_all_take_complex_input

comment:5 follow-up: Changed 7 years ago by ppurka

  • Commit set to c4ba30d802f932b29205ec3f49f3eeeeb0942e76
  • Status changed from new to needs_review

The following functions are not fixed: atan, acos, asin because I am not very sure how they get evaluated. I have fixed the other functions to take complex type as input.

Now, on to the plot ticket to check the output *phew*


New commits:

b34ee39allow some trigonometric functions to take complex as input
c4ba30dadd couple of doctests for complex input

comment:6 in reply to: ↑ 5 ; follow-ups: Changed 7 years ago by kcrisman

The following functions are not fixed: atan, acos, asin because I am not very sure how they get evaluated. I have fixed the other functions to take complex type as input.

You'll notice that here

arctan = atan = Function_arctan()

so atan is exactly the same as arctan - one fix works for both!

comment:7 in reply to: ↑ 6 Changed 7 years ago by kcrisman

so atan is exactly the same as arctan - one fix works for both!

(Though I haven't actually tried the branch.)


I understand and approve of the changes for trig.py, pending verification and messing around. But I don't quite get the changes in function.pyx. In particular, it seems like there is too much stuff. If org_parent is float then return org_parent(RDF(res)) does nothing (or at least nothing that wouldn't have happened with float, I hope). Why would that be better than just converting to float? I'd go straight to trying the CDF piece. Is there some specific case I'm not considering that would make me think that?

comment:8 in reply to: ↑ 6 Changed 7 years ago by kcrisman

The following functions are not fixed: atan, acos, asin because I am not very sure how they get evaluated. I have fixed the other functions to take complex type as input.

I see what you mean now, you also mean arctan, arccos, and arcsin. Sorry for my confusion. I expect they are evaluated in Ginac. I will investigate that a bit. But at least those guys plot fine.

comment:9 follow-up: Changed 7 years ago by kcrisman

Got it!

When GinacFunction goes up to BuiltinFunction it first tries doing it in Function. That coerces your complex to SR. (Incorrectly, perhaps, because

sage: SR(float(1))
1.0
sage: SR(complex(1,1))
(1+1j)

and I'm not sure what's up with the parentheses, but I don't really know much about coercion in Sage.)

Anyway, then it goes to the arcsin() method in sage.symbolic.expression.Expression (!) because now we have a symbolic expression for the previously complex input. And I'm pretty sure that this is where things go wrong, because here it looks like it tries to coerce our now-symbolic thing to RR, which of course we can't do.

Now you'll notice that in e.g. py_tan in the same file it then would go on to try coercing to CC, but it doesn't here. In fact, in a whole bunch of places this happens, and

sage: sinh(complex(1,1))
<snip>
TypeError: Unable to convert x (='(1+1j)') to real number.

so I am pretty sure this is where things end.

The reason your fix for the other arc trig functions works is because they are not defined in Ginac/Pynac?. But these functions don't work with complex input, presumably because here it doesn't handle complex stuff, or maybe because when we replaced Ginac's internal numbers with Sage ones something happened.

I'm cc:ing Burcin, who will probably have a suggestion as to the best place to fix this. But my druthers would be to just copycat what happens in sage/symbolic/pynac.pyx for other complex input, like py_exp or py_acosh.

Though as you say if we can't get it done quickly enough here, we can do this stuff and open a new ticket for it. But I really don't like that coercion to SR of complexes.

comment:10 Changed 7 years ago by kcrisman

  • Cc burcin added

Oops, cc:ing Burcin as promised.

Burcin, if you have a moment, I think there is a very easy question to answer about where to put a fix for something involving Ginac not doing complex input for arc trig functions. Thanks!

comment:11 in reply to: ↑ 9 Changed 7 years ago by ppurka

Replying to kcrisman:

Anyway, then it goes to the arcsin() method in sage.symbolic.expression.Expression (!) because now we have a symbolic expression for the previously complex input. And I'm pretty sure that this is where things go wrong, because here it looks like it tries to coerce our now-symbolic thing to RR, which of course we can't do.

Ah. I didn't figure out this step. I understood that it goes into Function, but lost track from there.

The reason your fix for the other arc trig functions works is because they are not defined in Ginac/Pynac?. But these functions don't work with complex input, presumably because here it doesn't handle complex stuff, or maybe because when we replaced Ginac's internal numbers with Sage ones something happened.

Yes. I know. These functions have their own _evalf_ method which gets executed finally. What is weird is that there are some other functions like sec, csc, etc, which also have an _evalf_ method, but that method probably never gets called! Instead the _eval_ method gets called. And, somehow that _eval_ method takes care of complex type input.

Though as you say if we can't get it done quickly enough here, we can do this stuff and open a new ticket for it. But I really don't like that coercion to SR of complexes.

Yes, all this multiple coercion is unclear to me also. SR is slow; and we should avoid doing those kinds of coercion. The only time when that should be done is when hold=True is passed as a parameter.

comment:12 Changed 7 years ago by git

  • Commit changed from c4ba30d802f932b29205ec3f49f3eeeeb0942e76 to e776b188e54d567b2112aadb7cd8e13e6e6cbe59

Branch pushed to git repo; I updated commit sha1. New commits:

4a6ef0bNeed to convert to complex only
e776b18fix unneeded line number in doctest

comment:13 Changed 7 years ago by kcrisman

Just one point - almost certainly some of those digits are 'extra' from Python. Compare (after this):

sage: sec(1.+i)
0.498337030555187 + 0.591083841721045*I
sage: sec(complex(1,1))
(0.49833703055518686+0.5910838417210451j)
sage: sec(CDF(complex(1,1)))
0.498337030555 + 0.591083841721*I

Interesting.

Here is something else I don't know what to do with.

sage: arccot(1.+i)
arccot(1.00000000000000 + 1.00000000000000*I)  # in SR
sage: cot(1.+i)
0.217621561854403 - 0.868014142895925*I

One would think this should be the same type of output.

That said, SR is useful for some other cases of exact input.

I'm wondering whether we want to try to fix all of this here, though.

comment:14 Changed 7 years ago by ppurka

I suggest that we do the minimum here so that the plot ticket can be doctested.

Fixing this entirely looks like a lengthy procedure, which someone very familiar with the coercion framework should do. I suspect there should be one big fix that should happen only in Function or BuiltinFunction such that it fixes all the other trigonometric and other functions.

comment:15 Changed 7 years ago by rws

  • Reviewers set to Karl-Dieter Crisman, Ralf Stephan
  • Status changed from needs_review to positive_review

As you say we will not be losing anything from doing this. The fact that it's all a minimal change, and patchbot is happy (ignore the spurious plugin failure) leads me to accept this. I would just ask Karl-Dieter and Burcin to feel free to open a ticket on what came up additionally.

comment:16 Changed 7 years ago by vbraun

  • Status changed from positive_review to needs_work

author name

comment:17 Changed 7 years ago by ppurka

  • Authors set to Punarbasu Purkayastha
  • Status changed from needs_work to positive_review

comment:18 Changed 7 years ago by vbraun

  • Status changed from positive_review to needs_work

Precision issues on multiple platforms

sage -t --long src/sage/functions/trig.py
**********************************************************************
File "src/sage/functions/trig.py", line 179, in sage.functions.trig.Function_sec._evalf_
Failed example:
    sec(complex(1,1))
Expected:
    (0.49833703055518686+0.5910838417210451j)
Got:
    (0.4983370305551868+0.591083841721045j)
**********************************************************************
File "src/sage/functions/trig.py", line 721, in sage.functions.trig.Function_arccsc._evalf_
Failed example:
    arccsc(complex(1,1))
Expected:
    (0.45227844715119064-0.5306375309525178j)
Got:
    (0.45227844715119064-0.5306375309525179j)
**********************************************************************
File "src/sage/functions/trig.py", line 797, in sage.functions.trig.Function_arcsec._evalf_
Failed example:
    arcsec(complex(1,1))
Expected:
    (1.118517879643706+0.5306375309525178j)
Got:
    (1.118517879643706+0.5306375309525179j)
**********************************************************************

Use # rel tol.

comment:19 Changed 7 years ago by kcrisman

Well, at least this justifies my comment about "extra" digits...

comment:20 Changed 7 years ago by git

  • Commit changed from e776b188e54d567b2112aadb7cd8e13e6e6cbe59 to 33449c150096c7db29eb46f2810f3694bfaf47bf

Branch pushed to git repo; I updated commit sha1. New commits:

33449c1add "# rel tol" to the tests

comment:21 Changed 7 years ago by ppurka

  • Status changed from needs_work to needs_review

@vbraun: added the # rel tol. Need review, because I don't have access to exotic archs :)

comment:22 follow-up: Changed 7 years ago by kcrisman

Neither do I, but I do know you need to add something WITH rel tol for it to work. See the developer guide (granted, only one example...)

comment:23 Changed 7 years ago by git

  • Commit changed from 33449c150096c7db29eb46f2810f3694bfaf47bf to 3355773583dad6fd88a8c7019a9fa19e305b3350

Branch pushed to git repo; I updated commit sha1. New commits:

3355773I hope I got the correct rel tol value :-/

comment:24 in reply to: ↑ 22 Changed 7 years ago by ppurka

Replying to kcrisman:

Neither do I, but I do know you need to add something WITH rel tol for it to work. See the developer guide (granted, only one example...)

Thanks. I just put in 1e-12 - as I understand it, rel tol gives the maximum magnitude of the error that must be tolerated. The results seem correct up to almost the 15th digit, but I put in a more conservative figure.

comment:25 Changed 7 years ago by vbraun

Special functions should always come out to within a few ulp, anything else is a bug IMHO. If it were really 1e-12 then that would be a bug for sure. I'd suggest 1e-15.

comment:26 Changed 7 years ago by git

  • Commit changed from 3355773583dad6fd88a8c7019a9fa19e305b3350 to 7c4c72596fe2a6de32e61ca578f2e671a51d7962

Branch pushed to git repo; I updated commit sha1. New commits:

7c4c725set it to 1e-15 as per vbraun advice

comment:27 Changed 7 years ago by vbraun

  • Status changed from needs_review to positive_review

comment:28 Changed 7 years ago by vbraun

  • Branch changed from u/ppurka/allow_inverse_trig_functions_to_all_take_complex_input to 7c4c72596fe2a6de32e61ca578f2e671a51d7962
  • Resolution set to fixed
  • Status changed from positive_review to closed
Note: See TracTickets for help on using tickets.