Ask Your Question
1

compile a symbolic function and use it in cython

asked 2012-01-17 00:32:36 -0600

pang gravatar image

updated 2012-01-18 11:15:26 -0600

Jason Grout gravatar image

I have written code to create plots of Newton fractals. Cython was key, pure python was way too slow.

But I'd like to create an interact: write the polynomial and watch the fractal. The problem is the user should enter the polynomial in a text box as a symbolic expression, then I should send that to cython and the cython code should do the rest, calling that function when necessary.

I've managed to do so: I create a fast_callable on the expression, and pass the output of fast_callable to cython, but it is very slow. Unfortunately, I cannot make the different data types for complex numbers match, and there might be something else.

Maybe you can think of a different approach...

edit retag flag offensive close merge delete

2 answers

Sort by » oldest newest most voted
1

answered 2012-01-18 11:34:36 -0600

Jason Grout gravatar image

updated 2012-01-18 13:19:52 -0600

kcrisman gravatar image

You probably want to use the fast_callable's call_c function. I've posted some examples at http://sagenb.org/home/pub/4086/, but I repeat them below for convenience:

%cython

from sage.ext.interpreters.wrapper_rdf cimport Wrapper_rdf
from libc.math cimport sqrt

def benchmark(float n):
    cdef int i
    for i in range(1000):
        n = sqrt(n)
    return n

def iterate_c(Wrapper_rdf f,float n):
    cdef double *result = [n]
    cdef int i
    for i in range(1000):
        f.call_c(result, result)
    return result[0]

def iterate_py(f, float n):
    cdef int i
    for i in range(1000):
        n=f(n)
    return n

Then

g(x)=sqrt(x)
f=fast_callable(g,vars=[x],domain=float)

Timings:

sage: timeit('benchmark(1e23)')
625 loops, best of 3: 7.15 µs per loop
sage: timeit('iterate_c(f,1e23)')
625 loops, best of 3: 18.3 µs per loop
sage: timeit('iterate_py(f,1e23)')
625 loops, best of 3: 162 µs per loop
sage: timeit('iterate_py(g,1e23)')
5 loops, best of 3: 77.2 ms per loop
edit flag offensive delete link more

Comments

I can't for the life of me figure out what you updated, kcrisman. Is there a place I can see the history? (Edit: oh, I see; just "edit" and then I can see the revisions...you just updated a link. I was afraid I made a coding error.)

Jason Grout gravatar imageJason Grout ( 2012-01-18 15:17:48 -0600 )edit

Awesome! Repeating the trick with "double complex", "double_complex", and a rusty knowledge of C was not easy, though.

pang gravatar imagepang ( 2012-01-19 00:03:48 -0600 )edit

@pang: if you end up with something even better (better looking, faster, or both) for the Newton basins, please add it to #11837. Sage will probably always be slower on this sort of thing, but it doesn't have to be ''too'' slow.

kcrisman gravatar imagekcrisman ( 2012-01-19 02:05:01 -0600 )edit
1

answered 2012-01-17 04:16:03 -0600

kcrisman gravatar image

I'm not sure if this answers your actual question, but you might be interested in code that Simon King and I have at # 11837 for generating Newton basin fractals. You are right that it is slow. Hopefully not quite as slow as yours, though, maybe it's all you need.

Incidentally, even at David Joyce's awesome generator page he just takes the roots; trying to get that from a generic polynomial first will severely slow you down, I'm sure.

edit flag offensive delete link more

Comments

Well, asking for the roots instead of the polynomial does the trick for polynomials, but the question remains: how can I get a symbolic function (or some string representing a function) from an interact, compile it, and use it in cython with similar speed to a function that was written in cython.

pang gravatar imagepang ( 2012-01-17 22:02:24 -0600 )edit

That's why I wasn't sure if it answered your actual question; I wasn't sure if you wanted Newton basins specifically, or the cython. I don't know how to use symbolic things in Cython, though Ginac is C++ so maybe there is something possible.

kcrisman gravatar imagekcrisman ( 2012-01-18 02:15:33 -0600 )edit

Well, both things are great to know. Thanks!

pang gravatar imagepang ( 2012-01-19 00:04:50 -0600 )edit

Your Answer

Please start posting anonymously - your entry will be published after you log in or create a new account.

Add Answer

Question Tools

Stats

Asked: 2012-01-17 00:32:36 -0600

Seen: 458 times

Last updated: Jan 18 '12