Re: tensor multiplication

[Paul van Delst <pvandelst(at)ncep.noaa.gov>]

Daniel Luebbert wrote:
> 
> Hi,
> 
> does anybody out there know an efficient and elegant way (i.e., without
> for-loops)
> to implement a tensor multiplication in IDL?
> 
> What I mean is this:
> IDL can do a matrix multiplication, e.g. if I do
>         c = indgen(3,4)
>         d = indgen(4)
> then for
>         help, c#d
> I get
>         LONG ARRAY[3],
> and that's what I expect.
> 
> But now, when I take one more dimension, like
>         c = indgen(2,3,4)
>         d = indgen(4)
> then
>         help, c#d
> gives an error! (incompatible matrix dimensions...).
> What a would like to get is obviously an
>         ARRAY[2,3]
> 
> Does anybody know how?

Great question! I have always only thought about 2D matrices, but why
should nD be any different? Maybe when IDL says "matrix" is really means
2-D array?

how about

sz = size(c)
e = INTARR( sz(1), sz(2) )
FOR i = 0, sz(1) - 1 DO BEGIN
  tmp_e  = REFORM( c[i,*,*] ) # d
  e[i,*] = TEMPORARY( tmp_e )
ENDFOR

I know it is not a great solution (I haven't tested it, just typed it)
but something like this should work. If you encapsulate it in it's own
function, you would have a general tensor mult. method.

I'm sure the Gumley's, Fanning's, Markwardt's, and JD Smith's of the
world will have more elegant answers.

paulv

p.s. Can someone explain to me the utility/need for having both the #
*and* ## operator? I understand their operation but why both?
Convenience? Performance? 
-- 
Paul van Delst           Ph:  (301) 763-8000 x7274 
CIMSS @ NOAA/NCEP        Fax: (301) 763-8545
Rm.202, 5200 Auth Rd.    Email: pvandelst@ncep.noaa.gov
Camp Springs MD 20746