how to do scatter/gather operations in numpy

Question

lets say i have arrays:

a = array((1,2,3,4,5))
indices = array((1,1,1,1))

and i perform operation:

a[indices] += 1

the result is

array([1, 3, 3, 4, 5])

in other words, the duplicates in indices are ignored

if I wanted the duplicates not to be ignored, resulting in:

array([1, 6, 3, 4, 5])

how would I go about this?

the example above is somewhat trivial, what follows is exactly what I am trying to do:

def inflate(self,pressure):
    faceforces = pressure * cross(self.verts[self.faces[:,1]]-self.verts[self.faces[:,0]], self.verts[self.faces[:,2]]-self.verts[self.faces[:,0]])
    self.verts[self.faces[:,0]] += faceforces
    self.verts[self.faces[:,1]] += faceforces
    self.verts[self.faces[:,2]] += faceforces

def constrain_lengths(self):
    vectors = self.verts[self.constraints[:,1]] - self.verts[self.constraints[:,0]]
    lengths = sqrt(sum(square(vectors), axis=1))
    correction = 0.5 * (vectors.T * (1 - (self.restlengths / lengths))).T
    self.verts[self.constraints[:,0]] += correction
    self.verts[self.constraints[:,1]] -= correction

def compute_normals(self):
    facenormals = cross(self.verts[self.faces[:,1]]-self.verts[self.faces[:,0]], self.verts[self.faces[:,2]]-self.verts[self.faces[:,0]])
    self.normals.fill(0)
    self.normals[self.faces[:,0]] += facenormals
    self.normals[self.faces[:,1]] += facenormals
    self.normals[self.faces[:,2]] += facenormals
    lengths = sqrt(sum(square(self.normals), axis=1))
    self.normals = (self.normals.T / lengths).T

Ive been getting some very buggy results as a result of duplicates being ignored in my indexed assignment operations.

Answer 1

I don't know of a way to do it that is any faster than:

for face in self.faces[:,0]:
    self.verts[face] += faceforces

You could also make self.faces into an array of 3 dictionaries where the keys correspond to the face and the value to the number of times it needs to be added. You'd then get code like:

for face in self.faces[0]:
    self.verts[face] += self.faces[0][face]*faceforces

which might be faster. I do hope that someone comes up with a better way because I wanted to do this when trying to help someone speed-up their code earlier today.