dots - multiplications of nd-arrays
A call will look like
[i j k l] = indices;
c = dots([i l],a,[i j k],b,[l j]);
meaning: for all i,l,
c(i,l) = sum_(j,k) a(i,j,k)*b(l,j));
The sum is done on the indices of a and b that are not indices of c

"dots" works in place without use of extra memory: no use of "repmat" or "Tony's trick".

"dots" is fairly fast, thanks to the use of "ddot" (BLAS) and a MEX file.

The functionality is very similar to Paolo de Vela's "multiprod". The ambition is work more economicaly in CPU and RAM terms. Futher, the call syntax is different, but this will obviously be a matter of taste.