[Numpy-discussion] Fwd: Advanced selection, duplicate indices, and augmented assignment
Timothy Hochberg
tim.hochberg at ieee.org
Mon Jan 8 12:49:46 EST 2007
On 1/6/07, Travis Oliphant <oliphant at ee.byu.edu> wrote:
>
> Timothy Hochberg wrote:
> >
> >
> >
> > On 1/6/07, *Robert Kern* < robert.kern at gmail.com
> > <mailto:robert.kern at gmail.com>> wrote:
> >
> > Sean R. Lynch wrote:
> >
> > >>>> x = zeros((3,))
> > >>>> x[array([0, 1, 1])] += array([1, 1, 1])
> > >>>> x
> > > array([ 1., 1., 0.])
> > >
> > > If this worked the way I was hoping, the output would be [1 2 0]
> > because
> > > it would add to element 1 twice due to its duplication in the
> > advanced
> > > selection array.
> > >
> > > Is the current behavior intentional or is it an accident of
> > > implementation?
> >
> > It is an unavoidable consequence of the way Python interprets that
> > code and the
> > way numpy arrays are fundamentally implemented. See Travis
> > Oliphant's, Tim
> > Hochberg's and my posts in the thread "Histograms via indirect
> > index arrays" for
> > more details.
> >
> >
> http://projects.scipy.org/pipermail/numpy-discussion/2006-March/thread.html#6877
> >
> >
> > Do we have to revisit that thread? I seem to recall it getting kind of
> > cranky. To avoid reliving that, I will attempt dredge up the relevant
> > issue:
> >
> > "a[indx]+=b" should be the same as "a[indx]=a[indx]+b". All else
> > follow from that. If staring at that for a while doesn't enlighten
> > you, then you will have to read that thread.
> >
> > [ I suspect that in theory we could make the += form behave as you
> > expect, but that would break the identity above, which would confuse a
> > bunch of people]
> >
>
> I don't think we could make it work as he expects (and not radically
> change NumPy itself so that other operations are very different) because
> of the limitations of Python.
>
> The point was that it is too complicated to do any differently (and is
> probably impossible).
Oh, you really shouldn't use words like "impossible" unless you want someone
to figure out how to do it. Some people view that as a challenge ;-)
It appears like it is possible, although it requires a certain amount
dubious hackery. Basically, you need a hidden tag on arrays that keeps track
of whether that array was the result of an inplace op and what kind. This
gets combined with a little sys.getrefcount magic to determine if the object
in question is a temporary. Below is a sketch of how it probably could be
accomplished. Not that I necessarily think it's a good idea, particularly
since the implementation is so dicey.
Oh, it's quite possible there are some corner and not so corner cases that
leak through too.
import numpy as np
import sys
class BecauseImContrary(np.ndarray):
"""Super hacktastic way to get repeated values to work
>>> a = np.arange(5).view(BecauseImContrary)
>>> indices = [0, 1, 0, 2]
# The basic test case, does repeated adds work.
>>> a1 = a.copy()
>>> a1[indices] += 1
>>> print a1
[2 2 3 3 4]
# Now some tests to make sure it doesn't add things when it
shouldn't
>>> a1 = a.copy()
>>> a1[indices] = a1[indices] + 1
>>> print a1
[1 2 3 3 4]
>>> def nested(x):
... y = np.ones_like(x)
... y += 998
... return y
>>> hasattr(nested(a1), '_iop')
False
>>> a[indices] = a[indices]
>>> print a
[0 1 2 3 4]
>>> a1 = a.copy()
>>> a1 += 1
>>> print a1
[1 2 3 4 5]
"""
def __iadd__(self, other):
result = np.ndarray.__iadd__(self, other)
count = sys.getrefcount(self)
assert count >= 6
if count == 6:
result._iop = 'add'
return result
def __setitem__(self, key, value):
count = sys.getrefcount(value)
assert count >= 5
# Only works for 1D indices, cause I'm lazy.
if count == 5 and hasattr(value, '_iop') and isinstance(key, (
np.ndarray, list)):
key = np.asarray(key)
for i in key:
self[i] = 0
value = np.ones_like(key) * value
for i, x in zip(key, value):
self[i] += x
else:
np.ndarray.__setitem__(self, key, value)
if __name__ == '__main__':
import doctest, scratch
doctest.testmod(scratch)
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://mail.python.org/pipermail/numpy-discussion/attachments/20070108/aa1b71dc/attachment.html>
More information about the NumPy-Discussion
mailing list