Mailman 3 February 2013 - NumPy-Discussion

numpy.einsum bug?
by Jaakko Luttinen 22 Feb '13

22 Feb '13

Hi, Is this a bug in numpy.einsum? >>> np.einsum(3, [], 2, [], []) ValueError: If 'op_axes' or 'itershape' is not NULL in theiterator constructor, 'oa_ndim' must be greater than zero I think it should return 6 (i.e., 3*2). Regards, Jaakko

1 0

Seeking help and support for next-gen math modeling tools using Python
by Rob Clewley 20 Feb '13

20 Feb '13

Hi all, and apologies for a little cross-posting: First, thanks to those of you who have used and contributed to the PyDSTool math modeling environment [1]. This project has greatly benefitted from the underlying platform of numpy / scipy / matplotlib / ipython. Going forward I have three goals, for which I would like the urgent input and support of existing or potential users. (i) I have ideas for expanding PyDSTool with innovative tools in my research area, which is essentially for the reverse engineering of complex mechanisms in multi-scale dynamic systems [2]. These tools have already been prototyped and show promise, but they need a lot of work. (ii) I want to grow and develop the community of users who will help drive new ideas, provide feedback, and collaborate on writing and testing code for both the core and application aspects of PyDSTool. (iii) The first two goals will help me to expand the scientific / engineering applications and use cases of PyDSTool as well as further sustain the project in the long-term. I am applying for NSF funding to support these software and application goals over the next few years [3], but the proposal deadline is in just four weeks! If you are interested in helping in any way I would greatly appreciate your replies (off list) to either of the following queries: I need to better understand my existing and potential users, many of whom may not be registered on the sourceforge users list. Please tell me who you are and what you use PyDSTool for. If you are not using it yet but you’re interested in this area then please provide feedback regarding what you would like to see change. If you are interested in these future goals, even if you are not an existing user but may be in the future, please write a brief letter of support on a letterhead document that I will send in with the proposal as PDFs. I have sample text that I can send you, as well as my draft proposal’s introduction and specific aims. These letters can make a great deal of difference during review. Without funding, collaborators, user demand and community support, these more ambitious goals for PyDSTool will not happen, although I am committed to a basic level of maintenance. For instance, based on user feedback I am about to release an Ubuntu-based Live CD [4] that will allow users to try PyDSTool on any OS without having to install it. PyDSTool will also acquire an improved setup procedure and will be added to the NeuroDebian repository [5], among others. I am also finalizing an integrated interface to CUDA GPUs to perform fast parallel ODE solving [6]. Thanks for your time, Rob Clewley [1] http://pydstool.sourceforge.net [2] http://www.ni.gsu.edu/~rclewley/Research/index.html, and in particular http://www.ploscompbiol.org/article/info%3Adoi%2F10.1371%2Fjournal.pcbi.100… [3] NSF Software Infrastructure for Sustained Innovation (SI2-SSE) program solicitation: http://www.nsf.gov/pubs/2013/nsf13525/nsf13525.htm [4] http://help.ubuntu.com/community/LiveCD [5] http://neuro.debian.net/ [6] http://www.nvidia.com/object/cuda_home_new.html -- Robert Clewley, Ph.D. Assistant Professor Neuroscience Institute and Department of Mathematics and Statistics Georgia State University PO Box 5030 Atlanta, GA 30302, USA tel: 404-413-6420 fax: 404-413-5446 http://neuroscience.gsu.edu/rclewley.html

4 4

pull request: generalized ufunc signature fix and lineal algebra generalized ufuncs
by Oscar Villellas 20 Feb '13

20 Feb '13

Hello, At Continuum Analytics we've been working on a submodule implementing a set of lineal algebra operations as generalized ufuncs. This allows specifying arrays of lineal algebra problems to be computed with a single Python call, allowing broadcasting as well. As the vectorization is handled in the kernel, this gives a speed edge on the operations. We think this could be useful to the community and we want to share the work done. I've created a couple of pull-requests: The first one contains a fix for a bug in the handling of certain signatures in the gufuncs. This was found while building the submodule. The fix was done by Mark Wiebe, so credit should go to him :). https://github.com/numpy/numpy/pull/2953 The second pull request contains the submodule itself and builds on top of the previous fix. It contains a rst file that explains the submodule, enumerates the functions implemented and details some implementation bits. The entry point to the module is in written in Python and contains detailed docstrings. https://github.com/numpy/numpy/pull/2954 We are open to discussion and to make improvements to the code if needed, in order to adapt to NumPy standards. Thanks, Oscar.

3 4

savetxt trouble
by Neal Becker 20 Feb '13

20 Feb '13

I tried to save a vector as a csv, but it didn't work. The vector is: a[0,0] array([-0.70710678-0.70710678j, 0.70710678+0.70710678j, 0.70710678-0.70710678j, 0.70710678+0.70710678j, -0.70710678-0.70710678j, 0.70710678-0.70710678j, -0.70710678+0.70710678j, -0.70710678+0.70710678j, -0.70710678-0.70710678j, -0.70710678-0.70710678j, 0.70710678-0.70710678j, -0.70710678-0.70710678j, -0.70710678+0.70710678j, 0.70710678-0.70710678j, 0.70710678-0.70710678j, 0.70710678+0.70710678j, 0.70710678-0.70710678j, -0.70710678-0.70710678j, -0.70710678-0.70710678j, 0.70710678-0.70710678j, 0.70710678+0.70710678j, 0.70710678+0.70710678j, -0.70710678+0.70710678j, 0.70710678+0.70710678j, -0.70710678-0.70710678j, -0.70710678+0.70710678j, 0.70710678-0.70710678j, -0.70710678+0.70710678j, 0.70710678+0.70710678j, 0.70710678+0.70710678j, 0.70710678-0.70710678j, -0.70710678-0.70710678j, 0.70710678-0.70710678j, -0.70710678-0.70710678j, -0.70710678+0.70710678j, 0.70710678+0.70710678j, 0.70710678-0.70710678j, 0.70710678-0.70710678j, 0.70710678+0.70710678j, -0.70710678+0.70710678j]) np.savetxt ('test.out', a[0,0], delimiter=',') The saved txt file says: (-7.071067811865540120e-01+-7.071067811865411334e-01j) (7.071067811865535679e-01+7.071067811865415775e-01j) (7.071067811865422437e-01+-7.071067811865529018e-01j) (7.071067811865520136e-01+7.071067811865431318e-01j) ...

3 3

Array accumulation in numpy
by Tony Ladd 19 Feb '13

19 Feb '13

I want to accumulate elements of a vector (x) to an array (f) based on an index list (ind). For example: x=[1,2,3,4,5,6] ind=[1,3,9,3,4,1] f=np.zeros(10) What I want would be produced by the loop for i=range(6): f[ind[i]]=f[ind[i]]+x[i] The answer is f=array([ 0., 7., 0., 6., 5., 0., 0., 0., 0., 3.]) When I try to use implicit arguments f[ind]=f[ind]+x I get f=array([ 0., 6., 0., 4., 5., 0., 0., 0., 0., 3.]) So it takes the last value of x that is pointed to by ind and adds it to f, but its the wrong answer when there are repeats of the same entry in ind (e.g. 3 or 1) I realize my code is incorrect, but is there a way to make numpy accumulate without using loops? I would have thought so but I cannot find anything in the documentation. Would much appreciate any help - probably a really simple question. Thanks Tony -- Tony Ladd Chemical Engineering Department University of Florida Gainesville, Florida 32611-6005 USA Email: tladd-"(AT)"-che.ufl.edu Web http://ladd.che.ufl.edu Tel: (352)-392-6509 FAX: (352)-392-9514

4 4

Re: [Numpy-discussion] Array accumulation in numpy
by Alan G Isaac 19 Feb '13

19 Feb '13

x=[1,2,3,4,5,6] ind=[1,3,9,3,4,1] f=np.zeros(10) np.bincount(ind,x) hth, Alan Isaac

2 1

Suggestion for a speed improved accumarray receipe
by Michael Löffler 17 Feb '13

17 Feb '13

Hello numpy list! I had tried to contribute to the AccumarrayLike recipe cookbook, but access seems to be restricted and no registration allowed. The current recipe mimics the matlab version in most aspects, but concerning performance it's horribly slow. The following snippet handles only the most simple usecase, but gives about 16x speed improvement compared to the current receipe, so it would probably worth to also mention it in the cookbook: def accum_custom(accmap, a, func=np.sum): indices = np.where(np.ediff1d(accmap, to_begin=[1], to_end=[np.nan]))[0] vals = np.zeros(len(indices) - 1) for i in xrange(len(indices) - 1): vals[i] = func(a[indices[i]:indices[i+1]]) return vals accmap = np.repeat(np.arange(100000), 20) a = np.random.randn(accmap.size) %timeit accum(accmap, a, func=np.sum) >>> 1 loops, best of 3: 16.7 s per loop %timeit accum_custom(accmap, a, func=np.sum) >>> 1 loops, best of 3: 945 ms per loop Best regards, Michael

1 0

calling NumPy from Julia - a plea for fewer macros
by Steven G. Johnson 17 Feb '13

17 Feb '13

Dear NumPy developers, I've been working on a glue package that allows the Julia language (http://julialang.org/) to call Python routines easily https://github.com/stevengj/PyCall.jl and I'm using NumPy to pass multidimensional arrays between languages. Julia has the ability to call C functions directly (without writing C glue), and I've been exploiting this to write PyCall purely in Julia. (This is nice for a number of reasons; besides programming and linking convenience, it means that I can dynamically load different Python versions on the same machine, and don't need to recompile if e.g. NumPy is updated.) However, calling NumPy has been a challenge, because of NumPy's heavy reliance on macros in its C API. I wanted to make a couple of suggestions to keep in mind as you plan for NumPy 2.0: 1) Dynamically linking to NumPy's C API was challenging, to say the least. Assuming you stick with the PyArray_API lookup table of pointers, it would be much easier to call from other languages if you include e.g. a numpy.core.multiarray._ARRAY_API_NAMES variable in the Python module that is a list of strings giving the symbol names corresponding to the numpy.core.multiarray._ARRAY_API pointer. (Plus documentation, of course.) Currently I need to parse __multiarray_api.h to extract this information, which is somewhat hackish. 2) Please provide non-macro equivalents (exported in the _ARRAY_API symbol table or otherwise) of PyArray_NDIM etcetera to access PyArrayObject members. (e.g. call them PyArray_ndim etc. Note that inline functions are not enough, since they are not loadable dynamically.) Right now, the only ways[*] I can see to access this information are either to use C glue (which I want to avoid for the reasons above) or to call Python to access the __array_interface__ attribute (which is suboptimal from a performance standpoint). Thanks for all your efforts! Any feedback on PyCall would be welcome, too. --SGJ [*] A third way would be to parse ndarraytypes.h to extract the format of the PyArrayObject_fields structure, and use upcoming Julia support for accessing C struct types to read the fields. This is likely to require tracking NumPy releases carefully to avoid breakage, however, as well as involving some care with the PyObject_HEAD macro. PS. If you want to try out PyCall with NumPy, note that a patch to Julia is currently required for this to work: https://github.com/JuliaLang/julia/pull/2317

3 4

FFTW bindings now implement numpy.fft interface
by Henry Gomersall 17 Feb '13

17 Feb '13

Some of you may be interested in the latest release of my FFTW bindings. It can now serve as a drop in replacement* for numpy.fft and scipy.fftpack. This means you can get most of the speed-up of FFTW with a one line code change or monkey patch existing libraries. Lots of other goodness too of course. Source here: https://github.com/hgomersall/pyFFTW pypi here: http://pypi.python.org/pypi/pyFFTW docs here: http://hgomersall.github.com/pyFFTW/ It's GPL3 due to license restrictions on FFTW. Get in touch if you want a different license and I'm sure we can reach an agreement ;) Cheers, Henry *In the case where the input array is not a numpy array, it doesn't work. This was an oversight and will be fixed in the next release. That said, if you're converting from a list on every transform, you have better optimisations than using FFTW. One other small caveat in a corner case to do with repeated axes - described in the docs.

2 2

np.percentile docstring
by Andreas Hilboll 17 Feb '13

17 Feb '13

In my numpy 1.6.1 (from Ubuntu 12.04LTS repository), the docstring of np.percentile is wrong. I'm not just submitting a PR because I don't understand something. in the "Notes" and "Examples" sections, there seems to be some confusion if ``q`` should be in the range [0,1] or in [0,100]. The parameters section is very clear about this, it should be [0,100]. However, in the Examples section, it says >>> np.percentile(a, 0.5, axis=0) array([ 6.5, 4.5, 2.5]) The given array is clearly the result of a call to np.percentile(a, 50., axis=0). I thought the docs are auto-generated, and that the "array..." result of the docstring would be calculated by numpy while building the docs? Or am I misunderstanding something here? Cheers, Andreas.

2 1