Re: [Numpy-discussion] Optimizing recursive loop - updated example
Thanks for the hint. I thought about Cython myself but I was unable to get even the slightest speed gain out of it. Here is the equivalent Cython code with the timing and setup.py. I typed (I think). Am I missing something obvious? Cheers Robert On 25.07.2011 01:38, Joon Ro wrote: For those cases where you cannot vectorize the operation, numpy is usually does not help much. Try using Cython. You will be able to compile the part of the code and the loop will be much faster (can be more than 100 times faster). http://docs.cython.org/ -Joon On Sun, 24 Jul 2011 18:10:14 -0500, Robert Elsner ml...@re-factory.de wrote: Boiled it down a bit more to only include code that actually takes time. First time around I found the other variant more instructive because it shows the discrepancy between the DCT and the loop but might be confusing. Thus here the bare minimum that correctly calculates the coefficients of the first derivative from the coefficients of the Chebyshev polynomials. Cheers Robert On 25.07.2011 00:43, Robert Elsner wrote: Hey Everybody, I am approximating the derivative of nonperiodic functions on [-1,1] using Chebyshev polynomials. The implementation is straightforward and works well but is painfully slow. The routine wastes most of its time on a trivial operation (see comment in the code) Unfortunately the spectral coefficients are calculated recursively and thus I haven't found a way to speed up the code. I am aware of list comprehensions, the speed advantage of calling native numpy functions etc. But no luck yet finding an alternate formulation that speeds up the calculation. I attached a stripped down variant of the code. I am very grateful for tips or directions where to look for a solution (well I know writing a small C extension might be the way to go but Id love to have speedy Python) Note: The DCT takes almost no time compared to the statement in the loop. Cheers Robert ___ NumPy-Discussion mailing list NumPy-Discussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpy-discussion from optimization import chebyshev_derivative_gc import timeit # Timing code t = timeit.Timer( chebyshev_derivative_gc ) print t.repeat( 3, 1000 ) import numpy as np cimport numpy as np # Fixing datatypes DTYPE = np.float ctypedef np.float_t DTYPE_t def chebyshev_derivative_gc(): Calculate the first derivative of a function cdef int m = 6 cdef int N = 64 cdef np.ndarray[DTYPE_t, ndim = 2] a = np.random.random( ( m, N ) ) cdef np.ndarray[DTYPE_t, ndim = 2] b = np.zeros( (m, N), dtype = DTYPE ) cdef int l,g,h,j l = 2 g = 1 h = 0 b[:, N - l] = l * ( N - g ) * a[:, N - g] for j in range( N - l, h, -g ): # The next line is the offender - try by commenting out. b[:, j - g] = ( l * j ) * a[:, j] + b[:, j + g] Created on 25.07.2011 @author: robert from distutils.core import setup from distutils.extension import Extension from Cython.Distutils import build_ext ext_modules = [Extension( optimization, [optimization.pyx] )] setup( name = 'Testing Cython speedups', cmdclass = {'build_ext': build_ext}, ext_modules = ext_modules ) ___ NumPy-Discussion mailing list NumPy-Discussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpy-discussion
Re: [Numpy-discussion] Optimizing recursive loop - updated example
Yes I did. Slicing and Cython do not mix too well. Using an explicit loop fixes the problem. In case anybody is interested the code is attached. Thanks for your help Robert On 25.07.2011 12:30, Robert Elsner wrote: Thanks for the hint. I thought about Cython myself but I was unable to get even the slightest speed gain out of it. Here is the equivalent Cython code with the timing and setup.py. I typed (I think). Am I missing something obvious? Cheers Robert On 25.07.2011 01:38, Joon Ro wrote: For those cases where you cannot vectorize the operation, numpy is usually does not help much. Try using Cython. You will be able to compile the part of the code and the loop will be much faster (can be more than 100 times faster). http://docs.cython.org/ -Joon On Sun, 24 Jul 2011 18:10:14 -0500, Robert Elsner ml...@re-factory.de wrote: Boiled it down a bit more to only include code that actually takes time. First time around I found the other variant more instructive because it shows the discrepancy between the DCT and the loop but might be confusing. Thus here the bare minimum that correctly calculates the coefficients of the first derivative from the coefficients of the Chebyshev polynomials. Cheers Robert On 25.07.2011 00:43, Robert Elsner wrote: Hey Everybody, I am approximating the derivative of nonperiodic functions on [-1,1] using Chebyshev polynomials. The implementation is straightforward and works well but is painfully slow. The routine wastes most of its time on a trivial operation (see comment in the code) Unfortunately the spectral coefficients are calculated recursively and thus I haven't found a way to speed up the code. I am aware of list comprehensions, the speed advantage of calling native numpy functions etc. But no luck yet finding an alternate formulation that speeds up the calculation. I attached a stripped down variant of the code. I am very grateful for tips or directions where to look for a solution (well I know writing a small C extension might be the way to go but Id love to have speedy Python) Note: The DCT takes almost no time compared to the statement in the loop. Cheers Robert ___ NumPy-Discussion mailing list NumPy-Discussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpy-discussion ___ NumPy-Discussion mailing list NumPy-Discussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpy-discussion import numpy as np cimport numpy as np # Fixing datatypes DTYPE = np.float ctypedef np.float_t DTYPE_t def chebyshev_derivative_gc(): Calculate the first derivative of a function cdef int m = 6 cdef int N = 64 cdef np.ndarray[DTYPE_t, ndim = 2] a = np.random.random( ( m, N ) ) cdef np.ndarray[DTYPE_t, ndim = 2] b = np.zeros( (m, N), dtype = DTYPE ) cdef int i, j for i in xrange(m): b[i, N - 2] = 2 * ( N - 1 ) * a[i, N - 1] for i in xrange(m): for j in xrange( N - 2, 0, -1 ): # The next line is the offender - try by commenting out. b[i, j - 1] = ( 2 * j ) * a[i, j] + b[i, j + 1] ___ NumPy-Discussion mailing list NumPy-Discussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpy-discussion
Re: [Numpy-discussion] Optimizing recursive loop - updated example
Mon, 25 Jul 2011 12:30:34 +0200, Robert Elsner wrote: Thanks for the hint. I thought about Cython myself but I was unable to get even the slightest speed gain out of it. Here is the equivalent Cython code with the timing and setup.py. I typed (I think). Am I missing something obvious? Cython doesn't vectorize the slice operations such as b[:,j-g] but falls back to Numpy on them. You'll need to convert the slice notation to a loop to get speed gains. -- Pauli Virtanen ___ NumPy-Discussion mailing list NumPy-Discussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpy-discussion
Re: [Numpy-discussion] Optimizing recursive loop - updated example
On Sun, Jul 24, 2011 at 5:10 PM, Robert Elsner ml...@re-factory.de wrote: Boiled it down a bit more to only include code that actually takes time. First time around I found the other variant more instructive because it shows the discrepancy between the DCT and the loop but might be confusing. Thus here the bare minimum that correctly calculates the coefficients of the first derivative from the coefficients of the Chebyshev polynomials. Have you tried using an (inverse) discrete sine transform to get the derivative? dT_n/dx = n*U_{n-1}, where U_n is the Chebyshev polynomial of the second kind, sin((n+1)\theta)/sin(theta) where cos(\theta) = x. I don't believe the discrete sine transform is part of scipy, but you can just use the inverse fft instead. snip Chuck ___ NumPy-Discussion mailing list NumPy-Discussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpy-discussion
Re: [Numpy-discussion] Optimizing recursive loop - updated example
I didn't look into that but it definitely sounds interesting. Especially as the coefficient manipulation is mildly unstable for higher derivatives. Need to work out the math first though ;). Thanks for the hint. On 25.07.2011 15:59, Charles R Harris wrote: On Sun, Jul 24, 2011 at 5:10 PM, Robert Elsner ml...@re-factory.de wrote: Boiled it down a bit more to only include code that actually takes time. First time around I found the other variant more instructive because it shows the discrepancy between the DCT and the loop but might be confusing. Thus here the bare minimum that correctly calculates the coefficients of the first derivative from the coefficients of the Chebyshev polynomials. Have you tried using an (inverse) discrete sine transform to get the derivative? dT_n/dx = n*U_{n-1}, where U_n is the Chebyshev polynomial of the second kind, sin((n+1)\theta)/sin(theta) where cos(\theta) = x. I don't believe the discrete sine transform is part of scipy, but you can just use the inverse fft instead. snip Chuck ___ NumPy-Discussion mailing list NumPy-Discussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpy-discussion ___ NumPy-Discussion mailing list NumPy-Discussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpy-discussion
Re: [Numpy-discussion] code review/build test for datetime business day API
Hey all, On Tue, Jun 14, 2011 at 4:34 PM, Mark Wiebe mwwi...@gmail.com wrote: These functions are now fully implemented and documented. As always, code reviews are welcome here: https://github.com/numpy/numpy/pull/87 I haven't been keeping up with the datetime developments, but I noticed the introduction of more names into the root numpy namespace. About a year (or two?) ago at SciPy, there were discussions about organising the NumPy namespace for 2.0, and halting the introduction of new functions into the root namespace. What is the status quo? Regards Stéfan ___ NumPy-Discussion mailing list NumPy-Discussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpy-discussion
Re: [Numpy-discussion] code review/build test for datetime business day API
2011/7/25 Stéfan van der Walt ste...@sun.ac.za Hey all, On Tue, Jun 14, 2011 at 4:34 PM, Mark Wiebe mwwi...@gmail.com wrote: These functions are now fully implemented and documented. As always, code reviews are welcome here: https://github.com/numpy/numpy/pull/87 I haven't been keeping up with the datetime developments, but I noticed the introduction of more names into the root numpy namespace. About a year (or two?) ago at SciPy, there were discussions about organising the NumPy namespace for 2.0, and halting the introduction of new functions into the root namespace. What is the status quo? I'm trying to make things fit into the existing system as naturally as possible. The discussion you're talking about ideally should have resulted in some guideline documentation about namespaces, but I don't recall seeing something like that in a prominent place anywhere. -Mark Regards Stéfan ___ NumPy-Discussion mailing list NumPy-Discussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpy-discussion ___ NumPy-Discussion mailing list NumPy-Discussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpy-discussion
Re: [Numpy-discussion] code review/build test for datetime business day API
On Mon, Jul 25, 2011 at 12:35 PM, Mark Wiebe mwwi...@gmail.com wrote: I'm trying to make things fit into the existing system as naturally as possible. The discussion you're talking about ideally should have resulted in some guideline documentation about namespaces, but I don't recall seeing something like that in a prominent place anywhere. Probably should have! Either way, it's something to consider: if we introduce those functions now, people will start to use them where they are (np.xyz), introducing another change in usage comes 2.0 (or 3.0 or whichever). Regards Stéfan ___ NumPy-Discussion mailing list NumPy-Discussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpy-discussion
[Numpy-discussion] Problems with numpy binary for Python2.7 + OS X 10.6
As best I can tell, I have Python 2.7.2 for my system Python: [ijstokes@moose ~]$ python -V Python 2.7.2 [ijstokes@moose ~]$ which python /Library/Frameworks/Python.framework/Versions/2.7/bin/python however when I attempt to install the recent numpy binary python-2.7.2-macosx10.6.dmg I get stopped at the first stage of the install procedure with the error: numpy 1.6.1 can't be installed on this disk. numpy requires System Python 2.7 to install. Any idea what I might be doing wrong? Is it looking for /usr/bin/python2.7? For that, I only have up to 2.6 available. (and 2.5) Cheers, Ian ___ NumPy-Discussion mailing list NumPy-Discussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpy-discussion
Re: [Numpy-discussion] code review/build test for datetime business day API
2011/7/25 Stéfan van der Walt ste...@sun.ac.za On Mon, Jul 25, 2011 at 12:35 PM, Mark Wiebe mwwi...@gmail.com wrote: I'm trying to make things fit into the existing system as naturally as possible. The discussion you're talking about ideally should have resulted in some guideline documentation about namespaces, but I don't recall seeing something like that in a prominent place anywhere. Probably should have! Either way, it's something to consider: if we introduce those functions now, people will start to use them where they are (np.xyz), introducing another change in usage comes 2.0 (or 3.0 or whichever). Absolutely, do you have any suggestions about organizing the datetime functionality? It's as good a place as any to try applying a good namespace convention. -Mark Regards Stéfan ___ NumPy-Discussion mailing list NumPy-Discussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpy-discussion ___ NumPy-Discussion mailing list NumPy-Discussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpy-discussion
Re: [Numpy-discussion] code review/build test for datetime business day API
2011/7/25 Stéfan van der Walt ste...@sun.ac.za Hey all, On Tue, Jun 14, 2011 at 4:34 PM, Mark Wiebe mwwi...@gmail.com wrote: These functions are now fully implemented and documented. As always, code reviews are welcome here: https://github.com/numpy/numpy/pull/87 I haven't been keeping up with the datetime developments, but I noticed the introduction of more names into the root numpy namespace. About a year (or two?) ago at SciPy, there were discussions about organising the NumPy namespace for 2.0, and halting the introduction of new functions into the root namespace. What is the status quo? Datetime is now a numpy type, so to that extent is in the base namespace. One could maybe argue that the calender functions belong in another namespace, and some of what is in numpy/lib (poly1d, financial functions) should have been in separate namespaces to begin with. I'm not sure what else in datetime might belong in its own namespace. Chuck ___ NumPy-Discussion mailing list NumPy-Discussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpy-discussion
Re: [Numpy-discussion] code review/build test for datetime business day API
On Mon, Jul 25, 2011 at 1:02 PM, Mark Wiebe mwwi...@gmail.com wrote: Probably should have! Either way, it's something to consider: if we introduce those functions now, people will start to use them where they are (np.xyz), introducing another change in usage comes 2.0 (or 3.0 or whichever). Absolutely, do you have any suggestions about organizing the datetime functionality? It's as good a place as any to try applying a good namespace convention. The first thought that comes to mind is simply to keep them in a submodule, so that users can do something like: from numpy.datetime import some_date_func That convention should be very easy to support across the restructuring. The important thing then is to document clearly that these functions are available. Regards Stéfan ___ NumPy-Discussion mailing list NumPy-Discussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpy-discussion
Re: [Numpy-discussion] code review/build test for datetime business day API
2011/7/25 Stéfan van der Walt ste...@sun.ac.za On Mon, Jul 25, 2011 at 1:02 PM, Mark Wiebe mwwi...@gmail.com wrote: Probably should have! Either way, it's something to consider: if we introduce those functions now, people will start to use them where they are (np.xyz), introducing another change in usage comes 2.0 (or 3.0 or whichever). Absolutely, do you have any suggestions about organizing the datetime functionality? It's as good a place as any to try applying a good namespace convention. The first thought that comes to mind is simply to keep them in a submodule, so that users can do something like: from numpy.datetime import some_date_func That convention should be very easy to support across the restructuring. The important thing then is to document clearly that these functions are available. Can't use numpy.datetime, since that conflicts with Python's datetime library, especially in pylab. Can't use numpy.datetime64, since that's already the name of the scalar type. I don't like numpy.dt, that name belongs to a delta t variable in Python. I'm not sure what a good name for the namespace is, actually. -Mark Regards Stéfan ___ NumPy-Discussion mailing list NumPy-Discussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpy-discussion ___ NumPy-Discussion mailing list NumPy-Discussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpy-discussion
Re: [Numpy-discussion] code review/build test for datetime business day API
On Mon, Jul 25, 2011 at 03:52:48PM -0500, Mark Wiebe wrote: Can't use numpy.datetime, since that conflicts with Python's datetime library, especially in pylab. I don't understand that: isn't the point of namespaces to avoid those naming conflicts. To me that's just like saying that numpy.sum shouldn't be named sum because it would conflict with the sum builtin. My 2 euro cents (a endangered currency) Gael ___ NumPy-Discussion mailing list NumPy-Discussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpy-discussion
Re: [Numpy-discussion] code review/build test for datetime business day API
On Mon, Jul 25, 2011 at 3:00 PM, Gael Varoquaux gael.varoqu...@normalesup.org wrote: On Mon, Jul 25, 2011 at 03:52:48PM -0500, Mark Wiebe wrote: Can't use numpy.datetime, since that conflicts with Python's datetime library, especially in pylab. I don't understand that: isn't the point of namespaces to avoid those naming conflicts. To me that's just like saying that numpy.sum shouldn't be named sum because it would conflict with the sum builtin. It's just asking for import problems and general confusion to shadow a Python module, that's why we renamed io to npyio. I'm curious as to what would be in the module. My 2 euro cents (a endangered currency) Aren't they all? Chuck ___ NumPy-Discussion mailing list NumPy-Discussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpy-discussion
Re: [Numpy-discussion] code review/build test for datetime business day API
On Mon, Jul 25, 2011 at 2:11 PM, Charles R Harris charlesr.har...@gmail.com wrote: It's just asking for import problems and general confusion to shadow a Python module, that's why we renamed io to npyio. Why? Users can simply do import numpy.io as npyio ? Stéfan ___ NumPy-Discussion mailing list NumPy-Discussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpy-discussion
Re: [Numpy-discussion] code review/build test for datetime business day API
On Mon, Jul 25, 2011 at 1:52 PM, Mark Wiebe mwwi...@gmail.com wrote: Can't use numpy.datetime, since that conflicts with Python's datetime library, especially in pylab. We're allowed to name the modules under numpy whatever we like--people know that doing from numpy import * can (and already does) cause havoc. But maybe numpy.time would suffice as a grouping. Regards Stéfan ___ NumPy-Discussion mailing list NumPy-Discussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpy-discussion
Re: [Numpy-discussion] code review/build test for datetime business day API
On Monday, July 25, 2011, Gael Varoquaux gael.varoqu...@normalesup.org wrote: On Mon, Jul 25, 2011 at 03:52:48PM -0500, Mark Wiebe wrote: Can't use numpy.datetime, since that conflicts with Python's datetime library, especially in pylab. I don't understand that: isn't the point of namespaces to avoid those naming conflicts. To me that's just like saying that numpy.sum shouldn't be named sum because it would conflict with the sum builtin. My 2 euro cents (a endangered currency) Gael I think the problem is that numpy's datetime might not be a drop-in replacement for python's datetime. For operations in python where one would use python's sum, numpy's sum would produce effectively identical results. We even have a few bug reports in our tracker for the use of all and any being ever-so-slightly different than python's, and it can cause some confusion in pylab mode. Admittedly, though, I can't come up with a better name, either. My two cents (also endangered)... Ben Root ___ NumPy-Discussion mailing list NumPy-Discussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpy-discussion
Re: [Numpy-discussion] code review/build test for datetime business day API
2011/7/25 Stéfan van der Walt ste...@sun.ac.za: On Mon, Jul 25, 2011 at 2:11 PM, Charles R Harris charlesr.har...@gmail.com wrote: It's just asking for import problems and general confusion to shadow a Python module, that's why we renamed io to npyio. Why? Users can simply do import numpy.io as npyio ? IIRC this was changed because of a (now fixed) bug in 2to3. Skipper ___ NumPy-Discussion mailing list NumPy-Discussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpy-discussion
Re: [Numpy-discussion] Problems with numpy binary for Python2.7 + OS X 10.6
In article 4e2dcb72.3070...@hkl.hms.harvard.edu, Ian Stokes-Rees ijsto...@hkl.hms.harvard.edu wrote: As best I can tell, I have Python 2.7.2 for my system Python: [ijstokes@moose ~]$ python -V Python 2.7.2 [ijstokes@moose ~]$ which python /Library/Frameworks/Python.framework/Versions/2.7/bin/python however when I attempt to install the recent numpy binary python-2.7.2-macosx10.6.dmg I get stopped at the first stage of the install procedure with the error: numpy 1.6.1 can't be installed on this disk. numpy requires System Python 2.7 to install. Any idea what I might be doing wrong? Is it looking for /usr/bin/python2.7? For that, I only have up to 2.6 available. (and 2.5) Cheers, Ian I believe the error message is misleading (a known bug). From the path you are probably running python.org python (though it could be ActiveState or built from source). Assuming it really is python.org, the next question is: which of the two flavors of python.org Python do you have: the 10.3 version (which is 32-bit only, but very backward compatible), or the 10.6 version (which includes 64-bit support but requires MacOS X 10.6 or later)? There is a separate numpy installer for each (and unfortunately they are not listed near each other in the file list). Maybe you got that match wrong? If in doubt you could reinstall python from python.org. -- Russell ___ NumPy-Discussion mailing list NumPy-Discussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpy-discussion
Re: [Numpy-discussion] code review/build test for datetime business day API
2011/7/25 Stéfan van der Walt ste...@sun.ac.za On Mon, Jul 25, 2011 at 2:11 PM, Charles R Harris charlesr.har...@gmail.com wrote: It's just asking for import problems and general confusion to shadow a Python module, that's why we renamed io to npyio. Why? Users can simply do import numpy.io as npyio ? It caused problems with 2to3 for one thing because it was getting imported as io in the package. It is just a bad idea to shadow python modules and best avoided. Chuck ___ NumPy-Discussion mailing list NumPy-Discussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpy-discussion
[Numpy-discussion] Get a 1D slice of a 3D data set?
Hi all, I have a 3D orthogonal and non-uniform grid representing a scalar field. I'm using matplotlib.image.NonUniformImage() to plot it similarly to imshow(). What I'd like to do is plot the values of the scalar field across a specific line (say, from point A to B). Any suggestion? Thanks! ___ NumPy-Discussion mailing list NumPy-Discussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpy-discussion