outer idiom (was Re: [Numpy-discussion] Getting the indexes of the myarray.min())

[Álvaro Tejero Cantero]

Hello,

(for background, I am trying to get an array of interparticle distances,
in which calculation of supradiagonal elements is unneeded because of
symmetry. I am not doing anything about this now, though).

>>> import numarray.random_array as rdn
>>> N, D = 1000, 3 #number of particles and dimensions of space
>>> r = rnd.random([N,D])  # r[i] gives the D coordinates of particle i
                           # r[:,0] gives the all the x coordinates
What I was doing:

>>> r_rel = [[r[i]-r[j] for i in arange(N)] for j in arange N]

now Tim says:

> Try this:
> 
>  >>> import numarray as na
>  >>> r = na.arange(5)
>  >>> na.subtract.outer(r, r)
> array([[ 0, -1, -2, -3, -4],
>        [ 1,  0, -1, -2, -3],
>        [ 2,  1,  0, -1, -2],
>        [ 3,  2,  1,  0, -1],
>        [ 4,  3,  2,  1,  0]])
> 

but this gives
>>> subtract.outer(r,r).shape
(10, 3, 10, 3)

that is, subtracts y coordinates to x coordinates which is not intended.
AFAIK the outer solution is MUCH faster than the nested for loops, so
what I do now is

>>> r_rel = transpose(array([subtract.outer(r[:,0],r[:,0]),
                             subtract.outer(r[:,1],r[:,1]), 
                             subtract.outer(r[:,2],r[:,2])]))
>>> r_rel.shape #as with the double loop 
(10,10,3) 


My question is then if there is any more elegant way to do this,
especially giving as a result independence of the number of dimensions).

Maybe an "axis"  (=0 in this case?) keyword for the outer function would
be useful in this context?



Thanks for the helpful welcome to the list!,
 á.

PS. the problem with the min() function regarded a matrix of collision
times between particles, which is symmetrical. The elements are reals,
so I only expect to find one minimum.

-- 
Álvaro Tejero Cantero

http://alqua.org    --	documentos libres
			free documents