Mailman 3 July 2010 - NumPy-Discussion

Re: [Numpy-discussion] [Matplotlib-users] Vectorization
by Benjamin Root 02 Jul '10

02 Jul '10

I am moving this over to numpy-discussion maillist... I don't have a firm answer for you, but I did notice one issue in your code. You call arange(len(dx) - 1) for your loops, but you probably really need arange(1, len(dx) - 1) because you are accessing elements both after *and* before the current index. An index of -1 is actually valid because that means the last element of the array, and may not be what you intended. Ben Root On Fri, Jul 2, 2010 at 1:15 PM, Nicolas Bigaouette <nbigaouette(a)gmail.com>wrote: > Hi all, > > I don't really know where to ask, so here it is. > > I was able to vectorize the normalization calculation in quantum mechanics: > <phi|phi>. Basically it's a volume integral of a scalar field. Using: > >> norm = 0.0 >> for i in numpy.arange(len(dx)-1): >> for j in numpy.arange(len(dy)-1): >> for k in numpy.arange(len(dz)-1): >> norm += psi[k,j,i]**2 * dx[i] * dy[j] * dz[k] >> > if dead slow. I replaced that with: > >> norm = (psi**2 * >> dx*dy[:,numpy.newaxis]*dz[:,numpy.newaxis,numpy.newaxis]).sum() >> > which is almost instantanious. > > I want to do the same for the calculation of the kinetic energy: > <phi|p^2|phi>/2m. There is a laplacian in the volume integral which > complicates things: > >> K = 0.0 >> for i in numpy.arange(len(dx)-1): >> for j in numpy.arange(len(dy)-1): >> for k in numpy.arange(len(dz)-1): >> K += -0.5 * m * phi[k,j,i] * ( >> (phi[k,j,i-1] - 2.0*phi[k,j,i] + phi[k,j,i+1]) / >> dx[i]**2 >> + (phi[k,j-1,i] - 2.0*phi[k,j,i] + phi[k,j+1,i]) / >> dy[j]**2 >> + (phi[k-1,j,i] - 2.0*phi[k,j,i] + phi[k+1,j,i]) / >> dz[k]**2 >> ) >> > > My question is, how would I vectorize such loops? I don't know how I would > manage the "numpy.newaxis" code-foo with neighbours dependency... Any idea? > > Thanx! > > > ------------------------------------------------------------------------------ > This SF.net email is sponsored by Sprint > What will you do first with EVO, the first 4G phone? > Visit sprint.com/first -- http://p.sf.net/sfu/sprint-com-first > _______________________________________________ > Matplotlib-users mailing list > Matplotlib-users(a)lists.sourceforge.net > https://lists.sourceforge.net/lists/listinfo/matplotlib-users > >

4 4

memory leak using numpy and cvxopt
by Tillmann Falck 02 Jul '10

02 Jul '10

Hi all, I am hitting a memory leak with the combination of numpy and cvxopt.matrix. As I am not where it occurs, I am cross posting. On my machine (Fedora 13, x86_64) this example quickly eats up all my memory. ----------- from cvxopt import matrix import numpy as np N = 2000 X = np.ones((N, N)) Y = matrix(0.0, (N, N)) while True: Y[:N, :N] = X ----------- I don't hit the leak if copy blocks of 1-d arrays. Regards, Tillmann

2 1

[ANN] Bento (ex-toydist) 0.0.3
by David 02 Jul '10

02 Jul '10

Hi, I am pleased to announce the release 0.0.3 for Bento, the pythonic packaging solution. Wherease the 0.0.2 release was mostly about getting the simplest-still-useful subset of distutils features, this new release adds quite a few significant features: - Add hooks to customize arbitrary stages in bento (there is a hackish example which shows how to use waf to build a simple C extension). The API for this is still in flux, though - Parallel and reliable build of C extensions through yaku build library. - One file distribution: no need for your users to install any new packages, just include one single file into your package to build with bento - Improved documentation - 2.4 -> 2.7 support, tested on linux/windows/mac os x You can download bento on github: http://github.com/cournape/Bento cheers, David

5 6

numpy.all docstring reality check
by David Goldsmith 02 Jul '10

02 Jul '10

Hi, folks. Under Parameters, the docstring for numpy.core.fromnumeric.all says: "out : ndarray, optionalAlternative output array in which to place the result. It must have the same shape as the expected output and *the type is preserved*." [emphasis added].I assume this is a copy-and-paste-from-another-docstring "typo" (shouldn't it be (possibly ndarray of) bool), but I just wanted to double check. DG -- Mathematician: noun, someone who disavows certainty when their uncertainty set is non-empty, even if that set has measure zero. Hope: noun, that delusive spirit which escaped Pandora's jar and, with her lies, prevents mankind from committing a general suicide. (As interpreted by Robert Graves)

3 3

where support full broadcasting, right?
by David Goldsmith 02 Jul '10

02 Jul '10

Hi. The docstring (in the wiki) for where states: x, y : array_like, optionalValues from which to choose. *x* and *y* need to have the same shape as *condition*.But: >>> x = np.eye(2) >>> np.where(x,2,3) array([[2, 3], [3, 2]]) So apparently where supports broadcasting of scalars at least; does it provide full broadcasting support? Thanks! DG

1 0

Determine slices in a sorted array
by Laszlo Nagy 01 Jul '10

01 Jul '10

Given an array with two axes, sorted by a column 'SLICE_BY', how can I extract slice indexes for rows with the same 'SLICE_BY' value? Here is an example program, demonstrating the problem: from numpy import * a = random.randint(0,100,(20,4)) SLICE_BY = 0 # Make slices of array 'a' by column SLICE_BY a.sort(SLICE_BY) slices = [] prev_val = None sidx = -1 for rowidx,row in enumerate(a): val = row[SLICE_BY] if val!=prev_val: if prev_val is None: prev_val = val sidx = rowidx else: slices.append((prev_val,sidx,rowidx)) sidx = rowidx prev_val = val if sidx<a.shape[0]-1: slices.append((val,sidx,a.shape[0])) print a print slices This program would print: [[ 1 0 8 1] [ 4 5 17 9] [ 4 11 19 23] [11 12 24 23] [13 16 28 23] [14 26 29 36] [15 33 32 37] [20 38 38 40] [28 47 47 45] [33 50 50 57] [45 55 52 65] [47 67 60 65] [56 76 71 68] [61 76 71 78] [70 83 82 83] [89 83 84 85] [91 84 85 87] [95 96 86 88] [98 96 89 88] [99 98 92 88]] [(1, 0, 1), (4, 1, 3), (11, 3, 4), (13, 4, 5), (14, 5, 6), (15, 6, 7), (20, 7, 8), (28, 8, 9), (33, 9, 10), (45, 10, 11), (47, 11, 12), (56, 12, 13), (61, 13, 14), (70, 14, 15), (89, 15, 16), (91, 16, 17), (95, 17, 18), (98, 18, 19)] Altough my demonstration program is functionally correct, it is not efficient. I need to do this with 10 million rows. Number of slices is relatively small (10 to 10000). Is is possible to construct my "slices" with pure numpy functions? E.g. anything that does not involve big number of python bytecode instructions, constucting Python objects, referencing/dereferencing 10 million times etc. Thanks, Laszlo

1 0

[ANN] PyTables 2.2 released: enter the multi-core age
by Francesc Alted 01 Jul '10

01 Jul '10

================================= Announcing PyTables 2.2 (final) ================================= I'm happy to announce PyTables 2.2 (final). After 18 months of continuous development and testing, this is, by far, the most powerful and well-tested release ever. I hope you like it too. What's new ========== The main new features in 2.2 series are: * A new compressor called Blosc, designed to read/write data to/from memory at speeds that can be faster than a system `memcpy()` call. With it, many internal PyTables operations that are currently bounded by CPU or I/O bandwith are speed-up. Some benchmarks: http://blosc.pytables.org/trac/wiki/SyntheticBenchmarks And a demonstration on how Blosc can improve PyTables performance: http://www.pytables.org/docs/manual/ch05.html#chunksizeFineTune * Support for HDF5 hard links, soft links and external links (kind of mounting external filesystems). A new tutorial about its usage has been added to the 'Tutorials' chapter of User's Manual. See: http://www.pytables.org/docs/manual/ch03.html#LinksTutorial * A new `tables.Expr` module (based on Numexpr) that allows to do persistent, on-disk computations on many algebraic operations. For a brief look on its performance, see: http://pytables.org/moin/ComputingKernel * Suport for 'fancy' indexing (i.e., à la NumPy) in all the data containers in PyTables. Backported from the implementation in the h5py project. Thanks to Andrew Collette for his fine work on this! * Binaries for both Windows 32-bit and 64-bit are provided now. As always, a large amount of bugs have been addressed and squashed too. In case you want to know more in detail what has changed in this version, have a look at: http://www.pytables.org/moin/ReleaseNotes/Release_2.2 You can download a source package with generated PDF and HTML docs, as well as binaries for Windows, from: http://www.pytables.org/download/preliminary For an on-line version of the manual, visit: http://www.pytables.org/docs/manual-2.2 What it is? =========== PyTables is a library for managing hierarchical datasets and designed to efficiently cope with extremely large amounts of data with support for full 64-bit file addressing. PyTables runs on top of the HDF5 library and NumPy package for achieving maximum throughput and convenient use. Resources ========= About PyTables: http://www.pytables.org About the HDF5 library: http://hdfgroup.org/HDF5/ About NumPy: http://numpy.scipy.org/ Acknowledgments =============== Thanks to many users who provided feature improvements, patches, bug reports, support and suggestions. See the ``THANKS`` file in the distribution package for a (incomplete) list of contributors. Most specially, a lot of kudos go to the HDF5 and NumPy (and numarray!) makers. Without them, PyTables simply would not exist. Share your experience ===================== Let us know of any bugs, suggestions, gripes, kudos, etc. you may have. ---- **Enjoy data!** -- The PyTables Team -- Francesc Alted

1 1

Ticket #1223...
by David Goldsmith 01 Jul '10

01 Jul '10

...concerns the behavior of numpy.random.multivariate_normal; if that's of interest to you, I urge you to take a look at the comments (esp. mine :-) ); otherwise, please ignore the noise. Thanks! DG

4 7

sum up to a specific value
by Renato Fabbri 01 Jul '10

01 Jul '10

hi, i need to find which elements of an array sums up to an specific value any idea of how to do this? best, rf -- GNU/Linux User #479299 skype: fabbri.renato

3 5

NumPy re-factoring project
by Jason McCampbell 01 Jul '10

01 Jul '10

Hi everyone, This is a follow-up to Travis's message on the re-factoring project from May 25th and the subsequent discussion. For background, I am a developer at Enthought working on the NumPy re-factoring project with Travis and Scott. The immediate goal from our perspective is to re-factor the core of NumPy into two architectural layers: a true core that is CPython-independent and an interface layer that binds the core to CPython. A design proposal is now posted on the NumPy developer wiki: http://projects.scipy.org/numpy/wiki/NumPyRefactoring The write-up is a fairly high-level description of what we think the split will look like and how to deal with issues such as memory management. There are also placeholders listed as 'TBD' where more investigation is still needed and will be filled in over time. At the end of the page there is a section on the C API with a link to a function-by-function breakdown of the C API and whether the function belongs in the interface layer, the core, or need to be split between the two. All functions listed as 'core' will continue to have an interface-level wrapper with the same name to ensure source-compatibility. All of this, particularly the interface/core function designations, is a first analysis and in flux. The goal is to get the information out and elicit discussion and feedback from the community. Best regards, Jason Jason McCampbell Enthought, Inc.

15 61