Mailman 3 January 2011 - NumPy-Discussion

ANN: Numexpr 1.4.2 released
by Francesc Alted 25 Jan '11

25 Jan '11

========================== Announcing Numexpr 1.4.2 ========================== Numexpr is a fast numerical expression evaluator for NumPy. With it, expressions that operate on arrays (like "3*a+4*b") are accelerated and use less memory than doing the same calculation in Python. What's new ========== This is a maintenance release. The most annying issues have been fixed (including the reduction bugs introduced in 1.4 series). Also, several performance enhancements are included too. In case you want to know more in detail what has changed in this version, see: http://code.google.com/p/numexpr/wiki/ReleaseNotes or have a look at RELEASE_NOTES.txt in the tarball. Where I can find Numexpr? ========================= The project is hosted at Google code in: http://code.google.com/p/numexpr/ You can get the packages from PyPI as well: http://pypi.python.org/pypi/numexpr Share your experience ===================== Let us know of any bugs, suggestions, gripes, kudos, etc. you may have. Enjoy! -- Francesc Alted

1 0

Numexpr giving randomized results on arrays larger than 2047 elements
by John Salvatier 25 Jan '11

25 Jan '11

Hello, I have discovered a strange bug with numexpr. numexpr.evaluate gives randomized results on arrays larger than 2047 elements. The following program demonstrates this: from numpy import * from numexpr import evaluate def func(x): return evaluate("sum(x, axis = 0)") x = zeros(2048)+.01 print evaluate("sum(x, axis = 0)") print evaluate("sum(x, axis = 0)") For me this prints different results each time, for example: 11.67 14.84 If we set the size to 2047 I get consistent results. 20.47 20.47 Interestingly, if I do not add .01 to x, it consistently sums to 0.

4 8

How to improve performance of slow tri*_indices calculations?
by eat 25 Jan '11

25 Jan '11

Hi, Running on: In []: np.__version__ Out[]: '1.5.1' In []: sys.version Out[]: '2.7.1 (r271:86832, Nov 27 2010, 18:30:46) [MSC v.1500 32 bit (Intel)]' For the reference: In []: X= randn(10, 125) In []: timeit dot(X.T, X) 10000 loops, best of 3: 170 us per loop In []: X= randn(10, 250) In []: timeit dot(X.T, X) 1000 loops, best of 3: 671 us per loop In []: X= randn(10, 500) In []: timeit dot(X.T, X) 100 loops, best of 3: 5.15 ms per loop In []: X= randn(10, 1000) In []: timeit dot(X.T, X) 100 loops, best of 3: 20 ms per loop In []: X= randn(10, 2000) In []: timeit dot(X.T, X) 10 loops, best of 3: 80.7 ms per loop Performance of triu_indices: In []: timeit triu_indices(125) 1000 loops, best of 3: 662 us per loop In []: timeit triu_indices(250) 100 loops, best of 3: 2.55 ms per loop In []: timeit triu_indices(500) 100 loops, best of 3: 15 ms per loop In []: timeit triu_indices(1000) 10 loops, best of 3: 59.8 ms per loop In []: timeit triu_indices(2000) 1 loops, best of 3: 239 ms per loop So the tri*_indices calculations seems to be unreasonable slow compared to for example calculations of inner products. Now, just to compare for a very naive implementation of triu indices. In []: def iut(n): ..: r= np.empty(n* (n+ 1)/ 2, dtype= int) ..: c= r.copy() ..: a= np.arange(n) ..: m= 0 ..: for i in xrange(n): ..: ni= n- i ..: mni= m+ ni ..: r[m: mni]= i ..: c[m: mni]= a[i: n] ..: m+= ni ..: return (r, c) ..: Are we really calculating the same thing? In []: triu_indices(5) Out[]: (array([0, 0, 0, 0, 0, 1, 1, 1, 1, 2, 2, 2, 3, 3, 4]), array([0, 1, 2, 3, 4, 1, 2, 3, 4, 2, 3, 4, 3, 4, 4])) In []: iut(5) Out[]: (array([0, 0, 0, 0, 0, 1, 1, 1, 1, 2, 2, 2, 3, 3, 4]), array([0, 1, 2, 3, 4, 1, 2, 3, 4, 2, 3, 4, 3, 4, 4])) Seems so, and then its performance: In []: timeit iut(125) 1000 loops, best of 3: 992 us per loop In []: timeit iut(250) 100 loops, best of 3: 2.03 ms per loop In []: timeit iut(500) 100 loops, best of 3: 5.3 ms per loop In []: timeit iut(1000) 100 loops, best of 3: 13.9 ms per loop In []: timeit iut(2000) 10 loops, best of 3: 39.8 ms per loop Even the naive implementation is very slow, but allready outperforms triu_indices, when n is > 250! So finally my question is how one could substantially improve the performance of indices calculations? Regards, eat

3 3

need a better way to fill a grid
by John 24 Jan '11

24 Jan '11

Hello, I'm trying to cycle over some vectors (lat,lon,emissions) of irregularly spaced lat/lon spots, and values. I need to sum the values each contributing to grid on a regular lat lon grid. This is what I have presently, but it is too slow. Is there a more efficient way to do this? I would prefer not to create an external module (f2py, cython) unless there is really no way to make this more efficient... it's the looping through the grid I guess that takes so long. Thanks, john def grid_emissions(lon,lat,emissions,grid.dx, grid.dy, grid.outlat0, grid.outlon0, grid.nxmax, grid.nymax): """ sample the emissions into a grid to fold into model output """ dx = grid.dxout dy = grid.dyout # Generate a regular grid to fill with the sum of emissions xi = np.linspace(grid.outlon0, grid.outlon0+(grid.nxmax*grid.d), grid.nxmax) yi = np.linspace(grid.outlat0, grid.outlat0+(grid.nymax*grid.dy), grid.nymax) X, Y = np.meshgrid(yi, xi) Z = np.zeros(X.shape) for i,x in enumerate(xi): for j,y in enumerate(yi): Z[i,j] = np.sum( emissions[\ np.where(((lat>y-dy) & (lat<y+dy)) & ((lon>x-dx) & (lon<x+dx)))[0]]) return Z

3 5

OpenMP-fication of loops.
by Vladimir Voznesensky 23 Jan '11

23 Jan '11

Hello. I've hacked loops.c.src to use OpenMP. Is anybody here interested in my hacks? Cheers. - VV

3 8

64 bit Windows installers for NumPy?
by Peter 22 Jan '11

22 Jan '11

Hi all, Are there plans to provide official 64bit Windows installers for NumPy? I'm aware that Christoph Gohlke had been able to do this, since he offers unofficial plain builds and MKL builds for NumPy here: http://www.lfd.uci.edu/~gohlke/pythonlibs/ Regards, Peter

5 8

Cython workshop
by Dag Sverre Seljebotn 20 Jan '11

20 Jan '11

As we have funding for it, we're talking about organizing a Cython workshop sometimes this year (possibly in Munich, Germany, though it's not decided yet). It's still not clear how user-centric vs. developer-centric the workshop will be, or how strong a role numerical computation will have vs. more general language features. We're just getting in touch with people potentially interested in joining the workshop, and then we'll take it from there. Respond on the wiki page or on cython-dev. http://wiki.cython.org/workshop1 <http://wiki.cython.org/workshop1#preview> Dag Sverre

1 0

[Ann] EuroScipy 2011 - Call for papers
by Tiziano Zito 20 Jan '11

20 Jan '11

========================= Announcing EuroScipy 2011 ========================= --------------------------------------------- The 4th European meeting on Python in Science --------------------------------------------- **Paris, Ecole Normale Supérieure, August 25-28 2011** We are happy to announce the 4th EuroScipy meeting, in Paris, August 2011. The EuroSciPy meeting is a cross-disciplinary gathering focused on the use and development of the Python language in scientific research. This event strives to bring together both users and developers of scientific tools, as well as academic research and state of the art industry. Main topics =========== - Presentations of scientific tools and libraries using the Python language, including but not limited to: - vector and array manipulation - parallel computing - scientific visualization - scientific data flow and persistence - algorithms implemented or exposed in Python - web applications and portals for science and engineering. - Reports on the use of Python in scientific achievements or ongoing projects. - General-purpose Python tools that can be of special interest to the scientific community. Tutorials ========= There will be two tutorial tracks at the conference, an introductory one, to bring up to speed with the Python language as a scientific tool, and an advanced track, during which experts of the field will lecture on specific advanced topics such as advanced use of numpy, scientific visualization, software engineering... Keynote Speaker: Fernando Perez =============================== We are excited to welcome Fernando Perez (UC Berkeley, Helen Wills Neuroscience Institute, USA) as our keynote speaker. Fernando Perez is the original author of the enhanced interactive python shell IPython and a very active contributor to the Python for Science ecosystem. Important dates =============== Talk submission deadline: Sunday May 8 Program announced: Sunday May 29 Tutorials tracks: Thursday August 25 - Friday August 26 Conference track: Saturday August 27 - Sunday August 28 Call for papers =============== We are soliciting talks that discuss topics related to scientific computing using Python. These include applications, teaching, future development directions, and research. We welcome contributions from the industry as well as the academic world. Indeed, industrial research and development as well academic research face the challenge of mastering IT tools for exploration, modeling and analysis. We look forward to hearing your recent breakthroughs using Python! Submission guidelines ===================== - We solicit talk proposals in the form of a one-page long abstract. - Submissions whose main purpose is to promote a commercial product or service will be refused. - All accepted proposals must be presented at the EuroSciPy conference by at least one author. The one-page long abstracts are for conference planing and selection purposes only. We will later select papers for publication of post-proceedings in a peer-reviewed journal. How to submit an abstract ========================= To submit a talk to the EuroScipy conference follow the instructions here: http://www.euroscipy.org/card/euroscipy2011_call_for_papers Organizers ========== Chairs: - Gaël Varoquaux (INSERM, Unicog team, and INRIA, Parietal team) - Nicolas Chauvat (Logilab) Local organization committee: - Emmanuelle Gouillart (Saint-Gobain Recherche) - Jean-Philippe Chauvat (Logilab) Tutorial chair: - Valentin Haenel (MKP, Technische Universität Berlin) Program committee: - Chair: Tiziano Zito (MKP, Technische Universität Berlin) - Romain Brette (ENS Paris, DEC) - Emmanuelle Gouillart (Saint-Gobain Recherche) - Eric Lebigot (Laboratoire Kastler Brossel, Université Pierre et Marie Curie) - Konrad Hinsen (Soleil Synchrotron, CNRS) - Hans Petter Langtangen (Simula laboratories) - Jarrod Millman (UC Berkeley, Helen Wills NeuroScience institute) - Mike Müller (Python Academy) - Didrik Pinte (Enthought Inc) - Marc Poinot (ONERA) - Christophe Pradal (CIRAD/INRIA, Virtual Plantes team) - Andreas Schreiber (DLR) - Stéfan van der Walt (University of Stellenbosch) Website ======= http://www.euroscipy.org/conference/euroscipy_2011

1 0

Prime size FFT: bluestein transform vs general chirp/z transform ?
by David Cournapeau 19 Jan '11

19 Jan '11

Hi, I finally took the time to clean up my code to speed up prime-size FFT (which use a O(N^2) algo in both numpy and scipy). The code is there: https://github.com/cournape/numpy/tree/bluestein (most of the code is tests, because numpy.fft had almost none). Bottom line: it is used only for prime numbers, and is faster than the current code for complex transforms > 500. Because of python + inherent bluestein overhead, this is mostly useful for "long" fft (where the speed up is significant - already 100x speed up for prime size ~ 50000). Several comments: - the overhead is pretty significant (on my machine, bluestein transfrom is slower for prime size < 500) - it could be used as such for real transforms, but the overhead would be even more significant (there is no bluestein transform for real transforms, so one needs to re-rexpress real transforms in term of complex ones, multiplying the overhead by 2x). There are several alternatives to make things faster (Rader-like transform, as used by fftw), but I think this would be quite hard to do in python without significant slowdown, because the code cannot be vectorized. - one could also decide to provide a chirp-z transform, of which Bluestein transform is a special case. Maybe this is more adapted to scipy ? - more generic code will require a few simple (but not trivial) arithmetic-like functions (find prime factors, find generator of Z/nZ groups with n prime, etc...). Where should I put those ? cheers, David

3 4

adding numexpr user-provided ufunc evaluation
by John Salvatier 19 Jan '11

19 Jan '11

I am thinking about building a patch to implement user-provided ufunc evaluation for numexpr as was discussed in this thread http://www.mail-archive.com/numpy-discussion@scipy.org/msg26292.html. Does anyone have any advice for learning/understanding the numexpr source and/or suggestions about implementation? All advice is helpful :) Francesc Alted mentioned that Mark Wiebe's other NEP https://github.com/m-paradox/numpy/blob/mw_neps/doc/neps/deferred-ufunc-eva…, would provide this feature (along with many other benefits). It looks really amazing, but I don't have a good sense of how developed it is. Best Regards, John Salvatier

1 0