Mailman 3 March 2006 - NumPy-Discussion

Re: [SciPy-user] multiplying sparse matrices in scipy-0.4.8?
by Travis Oliphant 28 Mar '06

28 Mar '06

Ed Schofield wrote: > > > Yes, this is a bug. > > Travis, could you please take a look at this? The FORTRAN functions > dcscmucsc and dcscmucsr both seem to be returning incorrect values in > indptr. In fact, could you please explain the Python code in matmat() > that calls these functions? I'd like to understand what the "while 1" > loop is for, and in particular why we have > c, rowc, ptrc, irow, kcol, ierr = func(*args) > when c, rowc, etc are part of *args anyway. The arguments are input and output arguments because that is my understanding of f2py's preference if something is both input and ouput rather than using an in-place argument. The while True: loop is so that resizing can occur if the guess on how many non-zero elements to reserve for the output matrix is wrong. There were some problems with the re-entry that I've now fixed. The tests you added to SciPy now pass. -Travis

1 0

Ransom Proposals
by Tim Hochberg 28 Mar '06

28 Mar '06

Now that the fortran keyword is dead, I'll take time off from rejoicing, to briefly go over other proposals that were involved in that thread lest they get lost in the shuffle. * reshape -- The origination of that whole, enormous thread was with reshape and what its return behaviour should be. Currently it returns a view if possible, otherwise it returns a copy. That's evil behaviour. I believe that it should either always return a copy or return a view if possible and otherwise raise an exception. I think I actually prefer the later, although it's probably a more disruptive change, since its more flexible when combined with asarray and friends. One possible compromise would be to have numpy.reshape always copy, while array.reshape always returns a view. * ascontiguous -- In private email Chris Barker mentioned that the name ascontiguous was confusing, or at least not to his taste and suggested "something like" ascontiguous_array. I don't like that one, but it might worth considering something that matches asarray and asanyarray. ascontigarray looks okay to me, but it's quite possible that I've staring at this too long and that's just cryptic. * ndmin -- I still think ndmin should be spun off to a separate function. It's trivial to implement a fuction, call it paddims for lack of a better name (asatleastndimensionarray seems kind of long and cryptic!). It should have minimal performance since no data gets copied, and if it makes a difference I would be willing to implement it in C in need be so that that performance impact was minimized. def paddims(a, n): "return a view of 'a' with at least 'n' dimensions" dims = a.shape b = a.view() b.shape = (1,)*(n - len(dims)) + dims return b * copy -- Yes, I understand that the copy flag is probably not going to change for backward compatibility reasons if nothing else, but there was one last point I wanted to make about the copy flag. One of the warts that results from the copy flag, and I think that this is a common problem for functions that take parameters that switch their mode of operation, is that some combinations of input become nonsensical. Consider these five possibilities: array(anarray, dtype=matchingtype, order=matchingorder, copy=True) # OK; copied array(anarray, dtype=matchingtype, order=matchingorder, copy=False) # OK; not copied array(anarray, dtype=nonmatchingtype, order=nonmatchingorder, copy=True) # OK; copied array(anarray, dtype=nonmatchingtype, order=nonmatchingorder, copy=False) # Ugly; copied array(nonarray, dtype=whatever, order=whatever, copy=False) # Ugly; copied [Note that I've folded nonmatchingtype and nonmatchingorder together since they have the same effect] Of these five possibilities, two have results where the arguments and the action taken become uncoupled. One way to address this would be to change the name of the copy flag to something that matches reality: force_copy. However, that seems kind of pointless, since it still introduces as the underlying problem that some of the modes the array function can operate in are kind of bogus. Compare this to the case where the two primitives are array and asarray: array(anarray, dtype=matchingtype, order=matchingorder) # copied array(anarray, dtype=nonmatchingtype, order=nonmatchingorder) # copied asarray(anarray, dtype=matchingtype, order=matchingorder) # not copied asarray(anarray, dtype=nonmatchingtype, order=nonmatchingorder) # copied asarray(nonarray, dtype=whatever, order=whatever, copy=False) # copied There's still five cases, so the interface hasn't narrowed any[*], but all the possible argument combinations make sense (or raise a straightforward error). And think how much easier this behaviour is to explain! Anyway that's it for now and hopefully for a while. Regards, -tim [*] In reality it does narrow the interface because we already have asarray, but now you really need as array whereas before it was really just shorthand for array(copy=False).

13 47

scipybooksample | dtype
by Jelle Feringa / EZCT Architecture & Design Research 28 Mar '06

28 Mar '06

Hi all, I'm familiarizing myself with scipy, and working my way through the http://numeric.scipy.org/scipybooksample.pdf tutorial. the tutorial states: ------------------------------------- >>> a = array([[1,2,3],[4,5,6]]) >>> a.shape (2, 3) >>> a.dtype <type 'int32 arrtype'> ------------------------------------- when i run this code on my machine, having numpy 0.9.6 installed, i get the following: >>> a = numpy.array([[1,2,3],[4,5,6]]) >>> a.dtype dtype('<i4') >>> numpy.version.version '0.9.6' >>> >>> seems that the way dtypes are being printed changed since the scipybooksample.pdf was released, would anyone of you be able to point me out where I can read up to what the updated dtypes signify? it's a bit cryptic to the uninitiated. i do understand that at the rate scipy/numpy is progressing that it can be hard to get the documentation consist with the current version. //btw: the new ndimage module is a great piece of work, thanks so much for all your efforts! -jelle

1 0

inplace add
by Tim Hochberg 28 Mar '06

28 Mar '06

Travis recently brought up the odd behaviour of adding an array to a list using inplace add, notably: >>> l = range(4) >>> a = arange(4) >>> l += a >>> l array([0, 2, 4, 6]) It was agreed that this was a deficiency in PyNumber_InplaceAdd and this is slated to be fixed. Hopefully for 2.5. Thanks Travis, Robert Kern and Armin Rigo for taking care of this. I had planned to write something here about, how, although this fixes the problem with lists, there is a similar problem with tuples, and since tuples don't have __iadd__, we could still expect to see this problem there. However, I now think that's not the case. The proposed resolution order will work something like: def iadd(a, b): if hasattr(a, '__iadd__'): tmp = a.__iadd__(b) if tmp is not NotImplemented: return tmp if hasattr(a, '__add__'): tmp = a.__add__(b) if tmp is not NotImplemented: return tmp if hasattr(b, '__radd__'): temp = b.__radd__ if tmp is not NotImplemented: return tmp raise TypeError("can't add %r and %r" % (type(a), type(b))) You'll notice if you try this, that "iadd(atuple, anarray)" throws an error because atuple refuses to be added to anything other than another tuple.That raises an error at "__add__" and prevents any oddities. So, we should be safe from both tuples and lists, and really they are the only types that are prevalent enough to worry about in this context. Still, I have a niggling feeling that iadd should really not be calling radd at all, but nothing concrete. -tim

1 0

arrays are column-wise or row-wise?
by Brian Blais 28 Mar '06

28 Mar '06

Hello, I thought I had read somewhere that the values in a numpy contiguous array are stored column-wise (moving down from the upper left to the lower left first) like C as opposed to row-wise (moving right from the upper left to the upper right first). However, when I write a simply Pyrex function to access the data pointer and print it out, it seems to go across the row. Is this correct? I am using numpy 0.9.3 and Python 2.3 on Linux. This also occurs for Numeric. I post the code below, with the sample output. thanks, bb -- ----------------- bblais(a)bryant.edu http://web.bryant.edu/~bblais #test_numpy.pyx: cdef extern from "numpy/arrayobject.h": struct PyArray_Descr: int type_num, elsize char type ctypedef class numpy.ArrayType [object PyArrayObject]: cdef char *data cdef int nd cdef int *dimensions, *strides cdef object base cdef PyArray_Descr *descr cdef int flags def test_numpy(ArrayType x): cdef int i cdef double *data cdef int length data=<double *>x.data length=x.dimensions[0]*x.dimensions[1] for i from 0 <= i < length: print data[i] #------------------------------------------------------ #testit.py: import numpy from test_numpy import * a=numpy.array([[1,2,3],[4,5,6],[7,8,9],[10,11,12]],numpy.Float64) print a test_numpy(a) #------------------------------------------------------ #output from program: >>> execfile('testit.py') [[ 1. 2. 3.] [ 4. 5. 6.] [ 7. 8. 9.] [ 10. 11. 12.]] 1.0 2.0 3.0 4.0 5.0 6.0 7.0 8.0 9.0 10.0 11.0 12.0

3 3

Re: [Numpy-discussion] Ransom Proposals
by Tim Hochberg 28 Mar '06

28 Mar '06

Charles R Harris wrote: > > > On 3/27/06, *Tim Hochberg* <tim.hochberg(a)cox.net > <mailto:tim.hochberg@cox.net>> wrote: > > Charles R Harris wrote: > [CHOP] > > > > > How about functions always return a copy unless they are carefully > > documented interface helper functions like asarray. That's how I > > generally think of them anyway. > > I don't really care as long as they (a) stop being evil and (b) > there's > some way to do stuff efficiently. Functions that always return copies > seem sort of useless to me, but I'd be happy to move to methods if > that > was the case. However, I suspect you won't make much headway with this > since my impression is that many of the people who do care about > functions also care passionately about efficiency as well. > > > Of course, the main virtue of asarray is that it helps write > functions > > that do precisely what Tim is arguing against: mixing list and array > > types and returning copies in the first case, views in the > second. So > > if we throw out functions we should throw out asarray also and > make a > > rule that numpy functions don't take lists. > > This isn't right. I think it's perfectly fine, desirable in fact, for > functions to operate on both list and array types. The problem only > arises when you return ether a view or copy of one of the original > inputs depending on what it was. In my experience, this is relatively > hard to do in all but the most trivial of functions since most > operations return a new array. However, it is true that I think people > should be strongly discouraged from doing this. For example: > > def add(a, b): > a = asarray(a) > b = asarray(b) > return a + b > > is fine, but: > > def iadd(a, b): > a = asarray(a) > b = asarray(b) > a += b > return a > > is not and should be rewritten. Possibly like: > > def iadd(a, b): > if not isinstance(a, ndarray): raise TypeError("'a' must be an > array") > b = asarray(b) > a += b > return a > > since += is a bit squirley. (Try somelist += somearray, to see > what I mean). > > > It's a syntax error, which I believe is the right thing. Really? If so, this is a version thing. This is what I get running on Python 2.4.1 and current numpy SVN (0.9.7.2286): >>> l = list(a) >>> l [999, 1, 2, 3, 4, 5, 6, 7, 8] >>> a array([999, 1, 2, 3, 4, 5, 6, 7, 8]) >>> l += a >>> l array([1998, 2, 4, 6, 8, 10, 12, 14, 16]) >>> a array([999, 1, 2, 3, 4, 5, 6, 7, 8]) > > > And then there are scalars... > > I'm curious what problems you forsee with scalars. > > > There was/is the whole problem of when to return scalars and if they > should be python scalars. The problem arises from trying to be > compatible with normal python stuff. Matlab got around this by making > everything 1) an array, and 2) double precision. On the other hand, > Matlab is not nearly is easy to use for nonstandard types. \ I think most of this problem will go away in Python with the introduction of the __index__ protocol, so for the moment at least I'm not worried about it. Regards, -tim

5 9

problem with ndarray attributes
by Christopher Hanley 28 Mar '06

28 Mar '06

Hi Travis, I have run across a feature in numpy which I am not sure how to handle. One of the new features of the ndarray is that named fields can be accessed as attributes. So, if a user created a recarray with fields named "address", "phone", and "names" they would be able to access all of the names in the recarray objects with "recarrayObject.names". I have found that this can cause problems when a user creates a table with a field name that happen to coincide with an existing class attribute or method. In particular, I am having a problem in pyfits where a user has a table with a column named "field". This column name is causing an attribute "field" to be created that overrides the existing class field method. It would seem to me that existing class methods and attributes should take precedence over the attributes created to correspond to column names. Is there a way to implement this behavior in ndarray? Thank you for your time and help, Chris

4 4

More about booleans as indices
by Pierre GM 28 Mar '06

28 Mar '06

Folks, I was playing with the compress method and extended slicing techniques, when I ran into this problem: I want to select the records of rec.array (`arec`) for which a given tag (`tag`) has some given values, as listed in `specialtaglist`. A basic arec[arec.tag in specialtaglist] raises a ValueError exception (ValueError: The truth value of an array with more than one element is ambiguous. Use a.any() or a.all()). A workaround I found was: arec[ [arec.tag[k] in specialtaglist for k in range(len(arec))] ] I was wondering whether there was a better/faster/simpler method ? Thx in advance

1 0

order flag again
by Tim Hochberg 27 Mar '06

27 Mar '06

I notice that some methods of the array object, notably reshape, take the order flag. I think this is a mistake. At this point, x.reshape() should always be returning a view. That means it's a cheap operation. Similarly, setting the order flag is also a cheap operation[*]. For that reason there's not a pressing efficiency reason to combine to the two operations and I feel that they should be broken apart. Thus we'd gain an operation 'setorder' or similar and lose the order flag everywhere but array and the various asXXXarray[**] functions since they need it to efficiently construct a correctly ordered array. Thus: b = a.reshape(newshape, order=neworder) would become: b = a.reshape(newshape).setorder(neworder) That doesn't look substantially better in that form, but compare: b = a.reshape(weeble, order=wobble) with: b = a.reshape(weeble).setorder(wobble) And the superiority of the second form begins to becomes apparent. The requirement that order be keyword parameter helps, but not enough IMO. The second form is orthogonal and self documenting and just all around more swell ;) Regards, -tim [*] I know this even though I haven't thought through what reshape plus the order flag is doing because reshape isn't allowed to return a copy. Thus setting the order flag in this context must amount to just rejiggering shape and strides. [**] And I suppose, maybe, ravel. However, ravel is yet another evil function, or rather method in this case, that sometimes returns a view and sometimes a copy. Ravels is really unnecessary, and I don't have the energy to wage a campaign against its evil ways (much to everyones relief I imagine), so I plan to just ignore it and recomend that no one use it . Ever.

3 5

using boolean arrays as indices
by Michael Sorich 27 Mar '06

27 Mar '06

If I have a 2D array and I want to use a boolean array to select the rows. If I use a[boolarray] it works, but if I use a[boolarray, :] an exception is raised. Is there a reason for this? eg >>>N.__version__ '0.9.6' >>>a = N.arange(10).reshape((2,5)) >>>a array([[0, 1, 2, 3, 4], [5, 6, 7, 8, 9]]) >>>a[N.array([True, False])] array([[0, 1, 2, 3, 4]]) >>>a[N.array([True, False]),:] Traceback (most recent call last): File "<input>", line 1, in ? IndexError: arrays used as indices must be of integer type

2 1