Mailman 3 December 2006 - NumPy-Discussion

Automatic matrices
by Paulo Jose da Silva e Silva Dec. 21, 2006

Dec. 21, 2006

Hello, If a numpy user is specially concerned with numerical linear algebra (or more generally with Math), it may find unconvenient the use of the dot function instead of the * operator. This behavior may be specially unpleasant for someone migrating from Matlab. I believe that this is the may reason for the existence of the matrix class within numpy. However, after trying to use the matrix class I have came across a major roadblock: many numpy/scipy functions return an array by default and not matrices. Then, we need then to add many conversion calls to the 'mat' function in our code. This is also unconvenient. I have had the idea of trying to write some code to automatically call the mat function for me. I have a very simple and inefficient prototype now that can be downloaded at: http://www.ime.usp.br/~pjssilva/matrix_import.py The use is simple. Instead of importing a numerical module with import, use the special class MatrixfiedModule from the above file. Let me give an example: --- ipython session --- In [2]:import matrix_import as mi In [3]:num = mi.MatrixfiedModule('numpy') Importing numpy In [4]:la = mi.MatrixfiedModule('scipy.linalg') Importing scipy.linalg In [5]:A = num.random.rand(3,4) Importing numpy.random In [6]:Q, R = la.qr(A) In [7]:la.norm(Q*R - A) Out[7]:6.0555516793379748e-16 ----- End session ----- For now the solution is very inefficient: every function call to a MatrixfiedModule function is wrapped on the fly to search for array return values ad convert them to matrix. This can certainly be improved the wrapping all the functions in the original module first. I plan to add this possibility soon. It is also incomplete: The automatic conversion only happens for return values of module function. It doesn't try to deal with special objects like finfo(float).eps or mgrid[0:9.,0:6.]. I am not sure how to deal with this. I can donate the code to scipy if there is any interest. Any comments? Best, Paulo -- Paulo José da Silva e Silva Professor Assistente, Dep. de Ciência da Computação (Assistant Professor, Computer Science Dept.) Universidade de São Paulo - Brazil e-mail: pjssilva(a)ime.usp.br Web: http://www.ime.usp.br/~pjssilva Teoria é o que não entendemos o (Theory is something we don't) suficiente para chamar de prática. (understand well enough to call) (practice)

5 11

Assignment when there is nothing to assign
by Keith Goodman Dec. 21, 2006

Dec. 21, 2006

I have the following two lines in a for loop: if len(idx): x[idx,i] = x[idx,i-1] where idx is the output of a where statement. I need the 'if len(idx)' line to prevent an error when idx is empty. Would it make sense to allow x[idx,i] = x[idx,i-1] when idx is empty instead of raising an error? The error I get is ValueError: array is not broadcastable to correct shape

1 1

ANN: PyTables 1.4 released
by Francesc Altet Dec. 21, 2006

Dec. 21, 2006

=========================== Announcing PyTables 1.4 =========================== 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. It is based on the HDF5 library for doing the I/O and leverages the numarray/NumPy/Numeric packages so as to deliver the data to the user in convenient in-memory containers. This is a new major release of PyTables, and probably the last major one of the 1.x series (i.e. with numarray at the core). On it, we have implemented better code to deal with table buffers, enhanced the capability for reading native HDF5 files, enhanced support for 64-bit platforms (but not with Python 2.5: see ``Special Warning`` section below), better support for AIX, optional automatic parent creation and the traditional amount of bug fixes. Go to the PyTables web site for downloading the beast: http://www.pytables.org/ or keep reading for more info about the new features and bugs fixed. Changes more in depth ===================== Improvements: - Table buffers code refactored: now each Row read iterator has its own buffers, completely independent of their table (although write iterators still share a single buffer in the same table). This separation makes the logic of buffering much more clear and less prone to errors (in fact, some of them have been solved). Performance and memory consumption are more or less equal than before. - When flushing the complete file (i.e. when calling File.flush()), only the buffers of those nodes that are alive (i.e. referenced from user code) are actually flushed. This brings much better efficiency (and also stability) to situations where one has to flush (and hence, close) files with many nodes on it. - Better support for AIX by renaming the internal LONLONG_MAX C constant (it was used internally by the xlc compiler). Thanks to Brian Granger for the report. - Added optional automatic parent creation support during node creation, copying and moving operations. See the release notes for more information. - Improved support for Python2.4 and 64-bit platforms (but beware, there are still known issues when using Python2.5 in combination with 64-bit platforms). Thanks to Gerard Vermeulen for his patches for Win64 platforms. - Implemented a workaround for a leak present in numarray --> Numeric conversions when using the array protocol, as can be seen in: http://comments.gmane.org/gmane.comp.python.numeric.general/12563 The workaround can potentially be far slower than the array protocol (because a copy of the arrays is always made), but at least the new code doesn't leak anymore. Bug fixes: - Previously, when the size for memory compounds type was less than the size of the type on disk (for example, when one have padding or aligned fields), PyTables was unable to read info on them. This has been fixed. This allows reading general compound types in HDF5 files written with other tools than PyTables. - When many tables with indexed columns were created simultaneously, a bug make PyTables to crash. This has been fixed (for more info, see bug #26). - Fixed a typo in the code that prevented recognizing complex data in non-PyTables files. - Table.createIndex() now refuses to index complex columns. - Now, it is possible to index several nested columns that hangs from the same column parent. Fixes bug #24. - Fixed a typo in nctoh5 utility that prevented using filters properly. Thanks to Lou Wicker for reporting this. - When setting/appending an array in-memory to an Array (or descendant) object and they have mismatched byteorders, the array was set/appended without being byteswapped first. This has been fixed. Thanks to Elias Collas for the report. Deprecated features: - None Backward-incompatible changes: - Please, see ``RELEASE-NOTES.txt`` file. Special Warning for Python 2.5 and 64-bit platforms users ========================================================= Unfortunately, and due to problems with the combination numarray 1.5.2, Python2.5 and 64-bit platforms, PyTables cannot be safely used yet in such scenario. This will be solved either when numarray can address this issue (hopefully with numarray 1.5.3), or when PyTables 2.x series (with NumPy at its core) will be out. Important note for Windows users ================================ If you are willing to use PyTables with Python 2.4 or 2.5 in Windows platforms, you will need to get the HDF5 library compiled for MSVC 7.1, aka .NET 2003. It can be found at: ftp://ftp.ncsa.uiuc.edu/HDF/HDF5/current/bin/windows/5-165-win-net.ZIP Users of Python 2.3 on Windows will have to download the version of HDF5 compiled with MSVC 6.0 available in: ftp://ftp.ncsa.uiuc.edu/HDF/HDF5/current/bin/windows/5-165-win.ZIP Platforms ========= This version has been extensively checked on quite a few platforms, like Linux on Intel32 (Pentium), Win on Intel32 (Pentium), Linux on Intel64 (Itanium2), FreeBSD on AMD64 (Opteron), Linux on PowerPC (and PowerPC64) and MacOSX on PowerPC. For other platforms, chances are that the code can be easily compiled and run without further issues. Please, contact us in case you are experiencing problems. Resources ========= Go to the PyTables web site for more details: http://www.pytables.org About the HDF5 library: http://hdf.ncsa.uiuc.edu/HDF5/ About numarray: http://www.stsci.edu/resources/software_hardware/numarray About NumPy: http://numpy.scipy.org/ To know more about the company behind the PyTables development, see: http://www.carabos.com/ Acknowledgments =============== Thanks to various the 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. Many thanks also to SourceForge who have helped to make and distribute this package! And last but not least, a big thank you to Acusim (http://www.acusim.com/) for sponsoring many of the job done for releasing this version of PyTables. Share your experience ===================== Let us know of any bugs, suggestions, gripes, kudos, etc. you may have. ---- **Enjoy data!** -- The PyTables Team

1 0

PyArray_DIMS problem
by Gennan Chen Dec. 20, 2006

Dec. 20, 2006

Hi! I have problem with this function call under FC6 X86_64 for my own numpy extension printf("\n %d %d %d", PyArray_DIM(imgi,0),PyArray_DIM(imgi, 1),PyArray_DIM(imgi,2)) it gave me 166 256 256 if I tried: int *dim; dim = PyArray_DIMS(imgi) printf("\n %d %d %d", dim[0], dim[1], dim[2]); it gave me 166 0 256 Numpy version: In [2]: numpy.__version__ Out[2]: '1.0.2.dev3487' I did test it under OS X 10.4.8 on MacPro. Those two methods gave me the exact results. So, what happens here ?? Gen-Nan Chen

4 4

A question about argmax and argsort
by asaf david Dec. 20, 2006

Dec. 20, 2006

Hello Let's say i have an N sized array, and i want to get the positions of the K largest items. for K = 1 this is simply argmax. is there any way to generalize it for k !=1? currently I use argsort and take only K items from it, but I'm paying an additional ~lg(N)... Thanks in advance, asaf

4 3

sum of two arrays with different shape?
by zhang yunfeng Dec. 20, 2006

Dec. 20, 2006

Hi, I'm newbie to Numpy. When reading tutorials at http://www.scipy.org/Tentative_NumPy_Tutorial, I found a snippet about addition of two arrays with different shape, Does it make sense? If array shapes are not same, why it doesn't throw out an error? see the code below (taken from the above webpage) array a.shape is (4,) and y.shape is (3,4) and a+y ? ------------------------------------------- >>> y = arange(12) >>> y array([ 0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11]) >>> y.shape = 3,4 # does not modify the total number of elements >>> y array([[ 0, 1, 2, 3], [ 4, 5, 6, 7], [ 8, 9, 10, 11]]) It is possible to operate with arrays of diferent dimensions as long as they fit well. >>> 3*a # multiply each element of a by 3 array([ 30, 60, 90, 120]) >>> a+y # sum a to each row of y array([[10, 21, 32, 43], [14, 25, 36, 47], [18, 29, 40, 51]]) -------------------------------------------- -- http://my.opera.com/zhangyunfeng

3 4

Profiling numpy ? (parts written in C)
by David Cournapeau Dec. 20, 2006

Dec. 20, 2006

Hi, Following the discussion on clip and other functions which *may* be slow in numpy, I would like to know if there is a way to easily profile numpy, ie functions which are written in C. For example, I am not sure to understand why a function like take(a, b) with a a double 256x4 array and b a 8000x256 int array takes almost 200 ms on a fairly fast CPU; in the source code, I can see that numpy uses memmove, and I know memmove to be slower than memcpy. Is there an easy way to check that this is coming from memmove (case in which nothing much can be done to improve the situation I guess), and not from something else ? cheers, David

7 17

(no subject)
by Bandler, Derek Dec. 19, 2006

Dec. 19, 2006

Hi, I would like to get information on the software licenses for numpy & numeric. On the sourceforge home for the packages, the listed license is OSI-Approved Open Source </softwaremap/trove_list.php?form_cat=14> . Is it possible to get more information on this? A copy of the document would be useful. Thank you. Best regards, Derek Bandler

3 2

Subclasses - use of __finalize__
by Colin J. Williams Dec. 18, 2006

Dec. 18, 2006

2 1

Python 2.5 on win and Gnuplot and alter_code1
by hermann.ruckerbauer＠qimonda.com Dec. 18, 2006

Dec. 18, 2006

Hello *, Sorry for this maybe stupid question, I might still miss some understanding of pyton or other stuff ... I have problems getting Gnuplot running under python 2.5 on Windows. My understanding is, that Gnuplot requires Numeric. Unfortionally there is only a numeric for python 2.4 available on the net. This binary does not install under python 2.5. Now I checked numpy, and it seems to able to act as replacement for the numeric package. I'm not really sure if this is the old.numeric module (just replace each numeric with numpy.old.numeric) or if this requires some more changes on the source files (e. g. exchang numeric with numpy and change some of the function calls). So I tried to convert my Gnuplot .py files with the alter_code1. But somehow I'm too stupid to get it running . First absolutel nothing happende when doing e.g. import numpy.oldnumeric.alter_code1 as noa noa.converttree(<top_level>) No error message, but also no file was touched, altough the command took a while to be executed. As I had no idea how to give the top level I tried a little bit around and copied the alter_code1.py into the gnuplot directory. And this file was touched ... It got the comment: At the beginning of alter_code1.py, but non of the other files have been touched. After some copy and paste trials i added on all .py files in the Gnuplot folder the comment from the beginning of the alter_code1.py. And now magically all files in this folders have been touched by running the converttree. BUT: only the comment "## Automatically adapted for numpy.oldnumeric Dec 15, 2006 by " Was added at the beginning of each file, none of the commands have been touched ... Does anybody have some hints for me to get gnuplot running under python 2.5 for windows ? Thanks and regards Hermann

3 3

Re: [Numpy-discussion] Unexpected output usingnumpy.ndarray and__radd__
by Mark Hoffmann Dec. 18, 2006

Dec. 18, 2006

Excellent, thank you - it solved the problem! /Mark -----Original Message----- From: numpy-discussion-bounces(a)scipy.org [mailto:numpy-discussion-bounces@scipy.org] On Behalf Of Tim Hochberg Sent: 18. december 2006 14:06 To: Discussion of Numerical Python Subject: Re: [Numpy-discussion] Unexpected output usingnumpy.ndarray and__radd__ Mark Hoffmann wrote: > I appreciate the answer and the solution suggestion. I see that it is possible to make a work around by subclassing from ndarray. Still, in the "print a+c" statement, I don't understand why a.__add__(c) doesn't return NotImplemented (because ndarray shouldn't recognize the Cyclehist class) and directly call c.__radd__(a) implemented in my Cyclehist class. I tried the exactly same programme using Python 2.4.1 and Scipy 0.3.2 (based on numeric/numarray) and the result of the "print a+c" didn't get nested as I expect. > > Regards, > Mark > I'm not sure what this is doing -- it looks kind of bizzare -- however, you can fix this case without resorting to subclassing to ndarray. Just toss an '__array_priority__ = 10' up at the top of the class definition and it will use your __methods__ in preference to the ndarrays. I don't have time to look into this further right now unfortunately. -tim > -----Original Message----- > From: numpy-discussion-bounces(a)scipy.org > [mailto:numpy-discussion-bounces@scipy.org] On Behalf Of Stefan van > der Walt > Sent: 18. december 2006 10:36 > To: Discussion of Numerical Python > Subject: Re: [Numpy-discussion] Unexpected output using numpy.ndarray > and__radd__ > > Hi Mark > > On Mon, Dec 18, 2006 at 08:30:20AM +0100, Mark Hoffmann wrote: > >> The following issue has puzzled me for a while. I want to add a >> numpy.ndarray and an instance of my own class. I define this >> operation by implementing the methods __add__ and __radd__. My programme (including output) looks like: >> >> #!/usr/local/bin/python >> >> import numpy >> >> class Cyclehist: >> def __init__(self,vals): >> self.valuearray = numpy.array(vals) >> >> def __str__(self): >> return 'Cyclehist object: valuearray = '+str(self.valuearray) >> >> def __add__(self,other): >> print "__add__ : ",self,other >> return self.valuearray + other >> >> def __radd__(self,other): >> print "__radd__ : ",self,other >> return other + self.valuearray >> >> c = Cyclehist([1.0,-21.2,3.2]) >> a = numpy.array([-1.0,2.2,-2.2]) >> print c + a >> print a + c >> > > In the first instance, c.__add__(a) is called, which works fine. In the second, a.__add__(c) is executed, which is your problem, since you rather want c.__radd__(a) to be executed. A documentation snippets: > > """For instance, to evaluate the expression x-y, where y is an instance of a class that has an __rsub__() method, y.__rsub__(x) is called if x.__sub__(y) returns NotImplemented. > > Note: If the right operand's type is a subclass of the left operand's type and that subclass provides the reflected method for the operation, this method will be called before the left operand's non-reflected method. This behavior allows subclasses to override their ancestors' operations.""" > > Since a.__add__ does not return NotImplemented, c.__radd__ is not called where you expect it to be. I am not sure why broadcasting takes place here, maybe someone else on the list can elaborate. > > To solve your problem, you may want to look into subclassing ndarrays, as described at http://www.scipy.org/Subclasses. > > Cheers > Stéfan > _______________________________________________ > Numpy-discussion mailing list > Numpy-discussion(a)scipy.org > http://projects.scipy.org/mailman/listinfo/numpy-discussion > _______________________________________________ > Numpy-discussion mailing list > Numpy-discussion(a)scipy.org > http://projects.scipy.org/mailman/listinfo/numpy-discussion > > > _______________________________________________ Numpy-discussion mailing list Numpy-discussion(a)scipy.org http://projects.scipy.org/mailman/listinfo/numpy-discussion

1 0

Re: [Numpy-discussion] Unexpected output using numpy.ndarray and__radd__
by Mark Hoffmann Dec. 18, 2006

Dec. 18, 2006

I appreciate the answer and the solution suggestion. I see that it is possible to make a work around by subclassing from ndarray. Still, in the "print a+c" statement, I don't understand why a.__add__(c) doesn't return NotImplemented (because ndarray shouldn't recognize the Cyclehist class) and directly call c.__radd__(a) implemented in my Cyclehist class. I tried the exactly same programme using Python 2.4.1 and Scipy 0.3.2 (based on numeric/numarray) and the result of the "print a+c" didn't get nested as I expect. Regards, Mark -----Original Message----- From: numpy-discussion-bounces(a)scipy.org [mailto:numpy-discussion-bounces@scipy.org] On Behalf Of Stefan van der Walt Sent: 18. december 2006 10:36 To: Discussion of Numerical Python Subject: Re: [Numpy-discussion] Unexpected output using numpy.ndarray and__radd__ Hi Mark On Mon, Dec 18, 2006 at 08:30:20AM +0100, Mark Hoffmann wrote: > The following issue has puzzled me for a while. I want to add a > numpy.ndarray and an instance of my own class. I define this operation > by implementing the methods __add__ and __radd__. My programme (including output) looks like: > > #!/usr/local/bin/python > > import numpy > > class Cyclehist: > def __init__(self,vals): > self.valuearray = numpy.array(vals) > > def __str__(self): > return 'Cyclehist object: valuearray = '+str(self.valuearray) > > def __add__(self,other): > print "__add__ : ",self,other > return self.valuearray + other > > def __radd__(self,other): > print "__radd__ : ",self,other > return other + self.valuearray > > c = Cyclehist([1.0,-21.2,3.2]) > a = numpy.array([-1.0,2.2,-2.2]) > print c + a > print a + c In the first instance, c.__add__(a) is called, which works fine. In the second, a.__add__(c) is executed, which is your problem, since you rather want c.__radd__(a) to be executed. A documentation snippets: """For instance, to evaluate the expression x-y, where y is an instance of a class that has an __rsub__() method, y.__rsub__(x) is called if x.__sub__(y) returns NotImplemented. Note: If the right operand's type is a subclass of the left operand's type and that subclass provides the reflected method for the operation, this method will be called before the left operand's non-reflected method. This behavior allows subclasses to override their ancestors' operations.""" Since a.__add__ does not return NotImplemented, c.__radd__ is not called where you expect it to be. I am not sure why broadcasting takes place here, maybe someone else on the list can elaborate. To solve your problem, you may want to look into subclassing ndarrays, as described at http://www.scipy.org/Subclasses. Cheers Stéfan _______________________________________________ Numpy-discussion mailing list Numpy-discussion(a)scipy.org http://projects.scipy.org/mailman/listinfo/numpy-discussion

2 1

Unexpected output using numpy.ndarray and __radd__
by Mark Hoffmann Dec. 18, 2006

Dec. 18, 2006

> Hi, > > The following issue has puzzled me for a while. I want to add a > numpy.ndarray and an instance of my own class. I define this operation > by implementing the methods __add__ and __radd__. My programme > (including output) looks like: > > #!/usr/local/bin/python > > import numpy > > class Cyclehist: > def __init__(self,vals): > self.valuearray = numpy.array(vals) > > def __str__(self): > return 'Cyclehist object: valuearray = '+str(self.valuearray) > > def __add__(self,other): > print "__add__ : ",self,other > return self.valuearray + other > > def __radd__(self,other): > print "__radd__ : ",self,other > return other + self.valuearray > > c = Cyclehist([1.0,-21.2,3.2]) > a = numpy.array([-1.0,2.2,-2.2]) > print c + a > print a + c > > # ---------- OUTPUT ---------- > # > # addprob $ addprob.py > # __add__ : Cyclehist object: valuearray = [ 1. -21.2 3.2] [-1. > 2.2 -2.2] > # [ 0. -19. 1.] > # __radd__ : Cyclehist object: valuearray = [ 1. -21.2 3.2] -1.0 > # __radd__ : Cyclehist object: valuearray = [ 1. -21.2 3.2] 2.2 > # __radd__ : Cyclehist object: valuearray = [ 1. -21.2 3.2] -2.2 > # [[ 0. -22.2 2.2] [ 3.2 -19. 5.4] [ -1.2 -23.4 1. ]] > # addprob $ > # > # ---------------------------- > > > I expected the output of "c+a" and "a+c" to be identical, however, the > output of "a+c" gets nested in an elementwise fashion. Can anybody > explain this? Is it a bug or a feature? I'm using Python 2.4.4c1 and > numpy 1.0. I tried the programme using an older version of Python and > numpy and there the result of "c+a" and "a+c" are identical. > > > Regards, > > Mark Hoffmann >

2 1

Anyone have a "little" shooting-method function to share
by David L Goldsmith Dec. 17, 2006

Dec. 17, 2006

Hi! I tried to send this earlier: it made it into my sent mail folder, but does not appear to have made it to the list. I need to numerically solve: (1-t)x" + x' - x = f(t), x(0) = x0, x(1) = x1 I've been trying to use (because it's the approach I inherited) an elementary finite-difference discretization, but unit tests have shown that that approach isn't working. After a little review, I believe I understand the problem: near t=0, the thing is like an advection-diffusion equation with diffusion as strong as advection (a case I'm sure is treated somewhere in the literature, but not as easily findable as:) near t=1, the thing is like the advection-diffusion equation I did easily find treated, namely one where advection dominates. Based on that treatment, I understand why an FD approach (with a uniform grid) would fail around t=0, and from there it is easy to accept that the fact that the thing changes its "nature" over the course of its "life," implies that (a uniform grid) FD is probably not a very good approach for this equation anywhere on its domain. I could try (and maybe will have to) a variable grid FD approach, but I'd first like to try a "shooting" method (for some reason my intuition is telling me that in this case this might be more efficient). This is where you all come in: I understand the algorithm and could (may have to) code it myself, but if anyone out there already has code for this (in Python using numpy, and preferably already under test), might you be willing to share? Thanks in advance, David Goldsmith PS: In the interim, I realized that if something like this "pre-exists" in a module, that module might be scipy. Sure enough, scipy has integrate.ode and integrate.odeint; Google-ing scipy: integrate.odeint help led me to Travis's 2004 "SciPy Tutorial," where I find "There are many optional inputs and outputs available when using odeint which can help tune the solver. These additional inputs and outputs are not needed much of the time, however..." Well, is one of those inputs an algorithm specifier? Alternatively, is there a "shooting" algorithm implemented elsewhere in scipy? ------------------------------------------------------------------------- Using Tomcat but need to do more? Need to support web services, security? Get stuff done quickly with pre-integrated technology to make your job easier Download IBM WebSphere Application Server v.1.0.1 based on Apache Geronimo http://sel.as-us.falkag.net/sel?cmd=lnk&kid=120709&bid=263057&dat=121642

3 4

Re: [Numpy-discussion] Pyrex and numpy
by Evan Lapisky Dec. 15, 2006

Dec. 15, 2006

> I was just going to try pyrex out with numpy to see how it compares with > weave (which is cool but quirky). My first attempt ended in failure: I > tried to compile the demo in in numpy/doc/pyrex and got this error: > > c_numpy.pxd:99:22: Array element cannot be a Python object > > > Does anyone who uses pyrex see this? Does anyone know what it's from? > Not that I deleted numpyx.c, since otherwise pyrex isn't invoked at all? > I had the same problem. I don't know why it didn't work, but the pyrex example from http://scipy.org/PerformancePython worked just fine. -Evan

2 1

Unexpected output using numpy.ndarray and __radd__
by Mark Hoffmann Dec. 15, 2006

Dec. 15, 2006

Hi, The following issue has puzzled me for a while. I want to add a numpy.ndarray and an instance of my own class. I define this operation by implementing the methods __add__ and __radd__. My programme (including output) looks like: #!/usr/local/bin/python import numpy class Cyclehist: def __init__(self,vals): self.valuearray = numpy.array(vals) def __str__(self): return 'Cyclehist object: valuearray = '+str(self.valuearray) def __add__(self,other): print "__add__ : ",self,other return self.valuearray + other def __radd__(self,other): print "__radd__ : ",self,other return other + self.valuearray c = Cyclehist([1.0,-21.2,3.2]) a = numpy.array([-1.0,2.2,-2.2]) print c + a print a + c # ---------- OUTPUT ---------- # # addprob $ addprob.py # __add__ : Cyclehist object: valuearray = [ 1. -21.2 3.2] [-1. 2.2 -2.2] # [ 0. -19. 1.] # __radd__ : Cyclehist object: valuearray = [ 1. -21.2 3.2] -1.0 # __radd__ : Cyclehist object: valuearray = [ 1. -21.2 3.2] 2.2 # __radd__ : Cyclehist object: valuearray = [ 1. -21.2 3.2] -2.2 # [[ 0. -22.2 2.2] [ 3.2 -19. 5.4] [ -1.2 -23.4 1. ]] # addprob $ # # ---------------------------- I expected the output of "c+a" and "a+c" to be identical, however, the output of "a+c" gets nested in an elementwise fashion. Can anybody explain this? Is it a bug or a feature? I'm using Python 2.4.4c1 and numpy 1.0. I tried the programme using an older version of Python and numpy and there the result of "c+a" and "a+c" are identical. Regards, Mark Hoffmann

1 0

migrating to numpy
by Christian Meesters Dec. 15, 2006

Dec. 15, 2006

Hi I was looking for a guide how to migrate old numeric code to numpy and couldn't find anything on the web. Any pointers for me? TIA and sorry for bothering, since this surely was discussed already, Christian

2 1

Histograms of extremely large data sets
by Cameron Walsh Dec. 15, 2006

Dec. 15, 2006

Hi all, I'm trying to generate histograms of extremely large datasets. I've tried a few methods, listed below, all with their own shortcomings. Mailing-list archive and google searches have not revealed any solutions. Method 1: import numpy import matplotlib data=numpy.empty((489,1000,1000),dtype="uint8") # Replace this line with actual data samples, but the size and types are correct. histogram = pylab.hist(data, bins=range(0,256)) pylab.xlim(0,256) pylab.show() The problem with this method is it appears to never finish. It is however, extremely fast for smaller data sets, like 5x1000x1000 (1-2 seconds) instead of 500x1000x1000. Method 2: import numpy import matplotlib data=numpy.empty((489,1000,1000),dtype="uint8") # Replace this line with actual data samples, but the size and types are correct. bins=numpy.zeros((256),dtype="uint32") for val in data.flat: bins[val]+=1 barchart = pylab.bar(xrange(256),bins,align="center") pylab.xlim(0,256) pylab.show() The problem with this method is it is incredibly slow, taking up to 30 seconds for a 1x1000x1000 sample, I have neither the patience nor the inclination to time a 500x1000x1000 sample. Method 3: import numpy data=numpy.empty((489,1000,1000),dtype="uint8") # Replace this line with actual data samples, but the size and types are correct. a=numpy.histogram(data,256) The problem with this one is: Traceback (most recent call last): File "<stdin>", line 1, in <module> File "/usr/local/lib/python2.5/site-packages/numpy/lib/function_base.py", line 96, in histogram n = sort(a).searchsorted(bins) ValueError: dimensions too large. It seems that iterating over the entire array and doing it manually is the slowest possible method, but that the rest are not much better. Is there a faster method available, or do I have to implement method 2 in C and submit the change as a patch? Thanks and best regards, Cameron.

8 21

Can we change how fortran compiler version strings are handled?!
by Brian Granger Dec. 14, 2006

Dec. 14, 2006

Hi, I have been doing quite a bit of numpy evangelism here at my work and slowly people are starting to use it. One of the main things people are interested in is f2py. But, I am finding that there is one persistent problem that keeps coming up when people try to install numpy on various systems: In three cases I have found that numpy failed to find and use a fortran compiler because the version string didn't match what was hardcoded into numpy.distutils. The reality is that version strings are in no way "standardized". In the most recent cases, we had a version of the lahey compiler that had the extra word "Express" and in another case, the xlf version string on a supercomputer was completely different. What is crazy to me is that this simple mismatch prevents numpy from even trying the compiler. Can we please change how Numpy handles the version string of fortran compilers? My suggestion would be to simply print the version string, but to attempt to use the compiler no matter what the version string is. That way, the success or failure of using the fortran compiler will be determined by the actual compiler, not its version string. There could be some other smart way of handling this, but I think it should be dealt with to make the installation process easier. I am willing to work up a patch if there is agreement on what should be done. Oh, the other difficult thing is that in the current arrangement, numpy.distutils doesn't print an error message that is easy to debug. It just silently does find the compiler rather than saying why. Thanks Brian

2 1

Pyrex and numpy
by Tim Hochberg Dec. 14, 2006

Dec. 14, 2006

I was just going to try pyrex out with numpy to see how it compares with weave (which is cool but quirky). My first attempt ended in failure: I tried to compile the demo in in numpy/doc/pyrex and got this error: c_numpy.pxd:99:22: Array element cannot be a Python object Does anyone who uses pyrex see this? Does anyone know what it's from? Not that I deleted numpyx.c, since otherwise pyrex isn't invoked at all? -tim

2 1

linalg.lstsq for complex
by Bill Baxter Dec. 13, 2006

Dec. 13, 2006

Is this code from linalg.lstsq for the complex case correct? lapack_routine = lapack_lite.zgelsd lwork = 1 rwork = zeros((lwork,), real_t) work = zeros((lwork,),t) results = lapack_routine(m, n, n_rhs, a, m, bstar, ldb, s, rcond, 0, work, -1, rwork, iwork, 0) lwork = int(abs(work[0])) rwork = zeros((lwork,),real_t) a_real = zeros((m,n),real_t) bstar_real = zeros((ldb,n_rhs,),real_t) results = lapack_lite.dgelsd(m, n, n_rhs, a_real, m, bstar_real, ldb, s, rcond, 0, rwork, -1, iwork, 0) lrwork = int(rwork[0]) work = zeros((lwork,), t) rwork = zeros((lrwork,), real_t) results = lapack_routine(m, n, n_rhs, a, m, bstar, ldb, s, rcond, The middle call to dgelsd looks unnecessary to me. At the very least, allocating astar_real and bstar_real shouldn't be necessary since they aren't referenced anywhere else in the lstsq function. The lapack documentation for zgelsd also doesn't mention any need to call dgelsd to compute the size of the work array. --bb

1 0

Re: [Numpy-discussion] rollaxis
by Pierre GM Dec. 13, 2006

Dec. 13, 2006

> Generate a column-permutation tuple and use fancy indexing: Works like a charm, thanks a lot !

1 0

rollaxis
by Pierre GM Dec. 13, 2006

Dec. 13, 2006

All, I have a ND array whose axes I want to reorganize, so that axis "i" is at the end while the others stay in their relative position. What's the easiest ?

2 2

fromfile and tofile access with a tempfile.TemporaryFile()
by Tim Hirzel Dec. 13, 2006

Dec. 13, 2006

Hi, Does anyone know how to get fromfile and tofile to work from a tempfile.TemporaryFile? Or if its not possible? I am getting this: >>> import tempfile >>> f = tempfile.TemporaryFile() >>> f <open file '<fdopen>', mode 'w+b' at 0x01EE1728> >>> a = numpy.arange(10) >>> a.tofile(f) Traceback (most recent call last): File "<input>", line 1, in ? IOError: first argument must be a string or open file thanks! tim

5 10

.byteswap() and copy/view dilemma
by Francesc Altet Dec. 13, 2006

Dec. 13, 2006

Hi, I'm a bit confused about which cases the .byteswap() method in NumPy would return a copy or a view. From the docstrings: """ a.byteswap(False) -> View or copy. Swap the bytes in the array. Swap the bytes in the array. Return the byteswapped array. If the first argument is TRUE, byteswap in-place and return a reference to self. """ >From the above description, it is not clear to me whether the .byteswap(False) would return a copy or a view. However, the official NumPy book seems to explain this clearer: """ Byteswap the elements of the array and return the byteswapped array. If the argument is True, then byteswap in-place and return a reference to self. Otherwise, return a copy of the array with the elements byteswapped. The data-type descriptor is not changed so the array will have changed numbers. """ My experiments also show that .byteswap(False) does do a copy: In [154]:a=numpy.array([1,2,3]) In [155]:b=a.byteswap() In [156]:b Out[156]:array([16777216, 33554432, 50331648]) In [157]:a[1]=0 In [158]:b Out[158]:array([16777216, 33554432, 50331648]) I'm wondering if I'm the only one that finds this docstring confusing. Regards, -- Francesc Altet | Be careful about using the following code -- Carabos Coop. V. | I've only proven that it works, www.carabos.com | I haven't tested it. -- Donald Knuth

1 0