Pau Gargallo wrote:
I think that the operation Angus wants is the following (at least I would like that one ;-)
[snip]
But sometimes it would be useful (at least for me) to have: a.shape = I+(n,k) b.shape = I+(k,m) and to get only: dot2( a,b ).shape == I+(n,m) dot2( a,b )[i] == dot2( a[i], b[i] )
This would be a natural extension of the scalar product (a*b)[i] == a[i]*b[i] If dot2 was a kind of ufunc, this will be the expected behaviour, while the current dot's behaviour will be obtained by dot2.outer(a,b).
Does this make any sense?
Thank-you Pau, this looks exactly like what I want. For further explanation, here, I believe, is an implementation of the desired routine using a loop. It would, however, be great to do this using quicker (ufunc?) machinery. Pau, can you confirm that this is the same as the routine you're interested in? def dot2(a,b): '''Returns dot product of last two dimensions of two 3-D arrays, threaded over first dimension.''' try: assert a.shape[1] == b.shape[2] assert a.shape[0] == b.shape[0] except AssertionError: print "incorrect input shapes" res = zeros( (a.shape[0], a.shape[1], a.shape[1]), dtype=float ) for i in range(a.shape[0]): res[i,...] = dot( a[i,...], b[i,...] ) return res I think the 'arrays of 2-D matrices' comment (which I've snipped out, oh well) captures the idea well. Angus. -- Angus McMorland email a.mcmorland@auckland.ac.nz mobile +64-21-155-4906 PhD Student, Neurophysiology / Multiphoton & Confocal Imaging Physiology, University of Auckland phone +64-9-3737-599 x89707 Armourer, Auckland University Fencing Secretary, Fencing North Inc. -- Angus McMorland email a.mcmorland@auckland.ac.nz mobile +64-21-155-4906 PhD Student, Neurophysiology / Multiphoton & Confocal Imaging Physiology, University of Auckland phone +64-9-3737-599 x89707 Armourer, Auckland University Fencing Secretary, Fencing North Inc.