This set of 3 files is aimed at efficiently extending matrix multiplication to higher dimensional arrays. It can also be used to compute the outer product of two arrays, and perform a trace over any two dimensions of an array.
In computing the matrix product of two arrays, first an outer product is taken, and then the resulting array is traced over the appropriate two dimensions. In my tests this has been 25-30 times faster than using "for" loops alone.
Example: to calculate the product of A and B, summing on the 2nd dimension of A and the 5th dimension of B, one would write C = gmdmp(A, 2, B, 5). Dimensions 2 of A and 5 of B need not be of identical length.
Errors: Due to use of the "squeeze" function, problems may be encountered with "diagsum" if there are singleton dimensions which are to be retained. This will usually just result in the output array having those singleton dimensions removed. However, an error may be encountered if the input array is of rank >2 and its last dimension is a singleton. |
