Suppose I have a square matrix $A_{mn}$ where $m=(i,j)$ and $n=(k,l)$. For example, if $i,j,k,l=1..N$ the matrix is $N^2\times N^2$. For simplicity, let me generate it by

n = 3;
A = RandomReal[1, {n, n, n, n}];

I want to find an inverse of this matrix: $A^{-1}_{(i,j)(k,l)}$. I understand how to it by introducing another square matrix $B=A_{mn}$ and inverting it:

B = ConstantArray[0, {n^2, n^2}];
 B[[i + n ( j - 1)]][[k + n ( l - 1)]] = A[[i]][[j]][[k]][[l]];
 , {i, 1, n}, {j, 1, n}, {k, 1, n}, {l, 1, n}
Bi = Inverse[B];
Ai = Table[
   Bi[[i + n ( j - 1)]][[k + n ( l - 1)]], {i, 1, n}, {j, 1, n}, {k, 
    1, n}, {l, 1, n}];

I believe there should be a better way of doing it(maybe using some tensor commands)?

B = ArrayFlatten[Transpose[A, {3, 1, 4, 2}]];
Ai = Transpose[ArrayReshape[Inverse[B], {n, n, n, n}], {2, 1, 4, 3}];
  • $\begingroup$ Thanks, it works! Is there an easy way to understand why you transpose these layers? $\endgroup$ Apr 22 at 9:24
  • 3
    $\begingroup$ Well, I just do it to get the ordering right. Please don't ask to figure out the right permutation! I just tried a couple of permutations until it worked! ;) $\endgroup$ Apr 22 at 9:26

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.