Slow CDF fast_callable powers

[robertwb]

Before

sage: var('x,y')
(x, y)
sage: f = (sum(x^k for k in [-10..10]) * sum(y^k for k in [-10..10])).expand()
sage: ff = fast_callable(f, CDF, (x,y))
sage: %timeit ff(1+2j, 3+4j)
125 loops, best of 3: 4.53 ms per loop
sage: len(ff.python_calls())
756
sage: ff.python_calls()[:10]
[(^10), (^10), (^10), (^9), (^9), (^10), (^10), (^8), (^9), (^9)]

which is still better than

sage: %timeit f(x=1+2j, y=3+4j)
25 loops, best of 3: 12 ms per loop

but nowhere close to

sage: var('x,y')
(x, y)
sage: f = (sum(x^k for k in [-10..10]) * sum(y^k for k in [-10..10])).expand()
sage: ff = fast_callable(f, CDF, (x,y))
sage: %timeit ff(1+2j, 3+4j)
625 loops, best of 3: 127 µs per loop
sage: len(ff.python_calls())
0
sage: 4.53 / .127
35.6692913385827

Apply attachment: 12513-cdf-pow-5.4.3.patch and attachment: 12513_solaris.patch

CC: @jasongrout

Component: numerical

Author: Robert Bradshaw

Reviewer: Timo Kluck, Jeroen Demeyer

Merged: sage-5.5.beta1

Issue created by migration from https://trac.sagemath.org/ticket/12513

[robertwb]

comment:2

Apply only 12513-cdf-pow-5.4.patch

[jasongrout]

comment:3

default is misspelled in the patch in one of the case statements (it's spelled as defualt).

[robertwb]

comment:4

Oops. Crazily enough, it interprets that as a label (so no compile error) and falls through to the base case (which also give the correct result, albeit slower).

[robertwb]

Attachment: 12513-cdf-pow-5.4.3.patch.gz

apply only this patch

[robertwb]

comment:5

Added some explicit tests of this codepath.

[jasongrout]

comment:6

Ah, I was wondering how it could compile. I guess you could have added a goto just to round it out :).

How does it know that the defualt: starts the base case? Are you sure it falls through to that? According to http://msdn.microsoft.com/en-us/library/66k51h7a(v=vs.80).aspx (Microsoft, I know...), for example, it seems that the default code would not be executed, "If the default statement is omitted, and no case match is found, none of the statements in the switch body are executed."

[robertwb]

comment:7

The "defualt" statement was not executed, so it falls all the way down to the last line return cpow(z, exp).

[jasongrout]

comment:8

Oh, right. I thought you meant the base case for the switch statement.

[tkluck]

comment:9

This looks excellent. I haven't tried the patch yet, but I noticed that you hard-code the integer range:

# supported for exponents that fit in an int 
self.ipow_range = (int(-2**31), int(2**31-1)) 

I'm not sure how platform independent that is. I think there's a preprocessor macro like MAX_INT or something. That's probably worth looking into.

[robertwb]

comment:10

Thanks for taking a look. Yes, INT_MAX is defined in stdint.c, but this is Python. (Note that sys.maxint is the wrong thing to use here, as that's the max value of a long.) The C standard mandates at least 16 bits, but I don't know of any (modern) processors that have less than 32-bit ints (including mobile processors like ARM).

Really, it's just the same bounds that are being used elsewhere in this file.

[tkluck]

Reviewer: tkluck

[tkluck]

comment:11

If there's no way in Cython to get the integer range, then I guess we should just leave it like this. I've tested the patch and it works as advertised. The doctests work, too. I'll give this a positive review.

[jdemeyer]

Changed reviewer from tkluck to Timo Kluck

[jdemeyer]

comment:13

Using the C99 constant I isn't as portable as it should be, it doesn't work on Solaris and OpenSolaris. I'll try a few things...

[jdemeyer]

Changed reviewer from Timo Kluck to Timo Kluck, Jeroen Demeyer

[jdemeyer]

comment:15

Attachment: 12513_solaris.patch.gz

[robertwb]

comment:16

"isn't as portable as it should be" yeah, it's a pretty basic part of the standard. It's sad when standard compliant code != portable code.

I assume you've tested this on solaris? If so, looks good to me. (Probably wouldn't hurt to have a solaris patchbot running...)

[jdemeyer]

Merged: sage-5.5.beta1