Python: conversion of a iterated assignment with an atomic assignment using numpy is not working when matrix height > 256

Question

I'm working using numpy 1.6.2 and python 2.7.
Given an N x M x D matrix A and a matrix I that contains a list of indices.
I have to fill a zeros matrix ACopy with the sum of element of A according to the indeces found in I (see code).

Here is my code:

ACopy = zeros(A.shape)
for j in xrange(0, size(A, 0)):
  i = I[j]
  ACopy[j, i, :] = A[j, i, :] + A[j, i + 1, :]

Indices matrix:

I = array([2, 0, 3, 2, 1])

A matrix:

A = array([[[ 0,  1,  2],
        [ 3,  4,  5],
        [ 6,  7,  8],
        [ 9, 10, 11],
        [12, 13, 14]],

       [[15, 16, 17],
        [18, 19, 20],
        [21, 22, 23],
        [24, 25, 26],
        [27, 28, 29]],

       [[30, 31, 32],
        [33, 34, 35],
        [36, 37, 38],
        [39, 40, 41],
        [42, 43, 44]],

       [[45, 46, 47],
        [48, 49, 50],
        [51, 52, 53],
        [54, 55, 56],
        [57, 58, 59]],

       [[60, 61, 62],
        [63, 64, 65],
        [66, 67, 68],
        [69, 70, 71],
        [72, 73, 74]]])

I try to improve my code avoiding the for loop in this way:

r = r_[0:len(I)]
ACopy[r, I, :] = A[r, I, :] + A[r, I + 1, :]

I noticed that the output matrices ACopy are different and I can't understand why. Any idea?

Thank you all!

EDIT: I'm computing a lot of matrices and I try with np.array_equals(ACopy1,ACopy2), where ACopy1 is the output of the first method and ACopy2 the output of the second method. Sometimes the matrices are the same, but not everytime. The two methods output should be the same or are there any bordeline case?

EDIT2: I noticed that this strange behaviour happens only when matrix height is bigger than 256. Here is my test suite:

from numpy import *

w = 5
h = 257

for i in xrange(1000):
  Z = random.rand(w, h, 5)

  I = (random.rand(w) * h - 1).astype(uint8)
  r = r_[0:w]
  ZCopy = zeros(Z.shape)
  ZCopy2 = zeros(Z.shape)

  for j in xrange(0, size(Z, 0)):
    i = I[j]
    ZCopy[j, i, :] = Z[j, i, :] + Z[j, i + 1, :]

  ZCopy2[r, I, :] = Z[r, I, :] + Z[r, I + 1, :]

  if (ZCopy - ZCopy2).any() > 0:
    print(ZCopy, ZCopy2, I)
    raise ValueError

Answer 1

I get the problem!
I cast the matrix I to uint8 and so matrix I elements are between 0 and 255.
I resolved using I = (random.rand(w) * h - 1).astype(uint32)