TRIQS / nda

C++ library for multi-dimensional arrays
https://triqs.github.io/nda
Other
13 stars 11 forks source link

Matrix multiplication between arrays #35

Closed jasonkaye closed 1 year ago

jasonkaye commented 1 year ago

The following type of operation is common: I have a matrix A, and a rank-3 array B, and I want to compute C_ikl = sum_j A_ij * B_jkl.

To do this with NDA, currently I am (1) taking a reshaped_view of B, to turn it into a rank-2 array with the proper leading dimension, (2) taking a matrix_view of the result, (3) performing the matrix-matrix product with the * operation, and (4) reshaping the result to get my rank-3 array C.

It would seem to be cleaner to abandon the matrix type, and simply define a special operator for matrix multiplication between arrays. Let's call this operator **. Now, A is a rank-2 array instead of a matrix, B is a rank 3 array, and I simply write

C = A ** B.

It is first checked that the last dimension of A matches the first dimension of B, and then it does the right thing. One could even allow A to be a rank-3 array, and to compute C_ijlm = sum_k A_ijk B_klm by C = A ** B; more generally, A and B could be arrays of any rank, as long as the last dimension of A matches the first dimension of B.