[Numpy-discussion] Fwd: Advanced selection, duplicate indices, and augmented assignment

[A. M. Archibald]

On 10/01/07, Charles R Harris <charlesr.harris at gmail.com> wrote:

> > You're the second or third person to recently come up with a useful
> application for using inplace ops with repeated indices only to be stymied
> by that fact that they don't actually work that way. I've been wondering if
> we should try to come up with a way to support this usage. Although it might
> be feasible to make += work, the implementation is horrible, and there are
> all sorts of inefficiencies built into that method of doing it.
> >
> > The natural approach would seem to be to hang another method off the
> ufuncs. That way you get support for all of the various operations. I was
> thinking something like:
> >
> >
> > unfunc.inplace (dest, indices, source)
>
> Yeah, something like that. A general indirect addressing primitive of some
> sort would be useful for these situations. The indices are sort of pointers
> into the destination, the problem is how to bind this to operators. It
> should almost be some sort lvalue analogous to **p = whatever, or in this
> case dest[index[j]] {assignment} src[j]. So maybe indirect(dest, indices,
> source, op). Then there is all the madness of fancy indexing (which I think
> I could have done without). Because the current applications are limited and
> could be met with the simplest 1D version, it would probably be best to
> start with the most basic functionality and then extend it if the need
> arose. Or maybe even some simple c extension that did the whole basic +=
> thing for one dimension.

I disagree that it should be hung off ufuncs. That doesn't allow, for
example, averaging of normal vectors (without some highly clever and
hard-to-follow axis specification).

I implemented one possible version in pure python, just to think about
the UI (and to see if it could be made fast using sorting; I don't
think it will actually help much, and I also note that whether it even
has any hope depends on whether there are more objects to accumulate
than destination slots, or fewer). The basic design is

def accumulate(dest, indices, values, combine=N.add);

indices should be a one-dimensional list of indices, though since dest
can be an n-dimensional array of r-dimensional values, it's actually a
two-dimensional array; values is a list of r-dimensional values, and
combine should take two r-dimensional values and yield another.

The reason for the r-dimensional values is to allow, for example,
averaging of normals, or combination in some other way of colours, or
summing of vectors (my original application).

This implementation does not allow multidimensional arrays of indices,
partly because it would be cumbersome to iterate over them, partly
because a reshape (and possible copy) is not terribly offensive, and
partly because it would make the calling convention (more) confusing.

I am not seriously suggesting this implementation for use, but for
consideration; if people find the primitive does what they want, one
could reasonably simply implement it in C.

A. M. Archibald
-------------- next part --------------
A non-text attachment was scrubbed...
Name: accumulate.py
Type: text/x-python
Size: 1397 bytes
Desc: not available
URL: <http://mail.python.org/pipermail/numpy-discussion/attachments/20070110/f3ec2e36/attachment.py>