I have a $B$ matrix: $B = B_{ij}$

I need to find closed matrix formulation of:

$$\sum_i \sum_j \sum_m \sum_n B_{ij} B_{jm} B_{mn} B_{ni}$$

But I am so confused!

Edit by Henrik:

Originally, it was asked to express

$$\sum_i \sum_j \sum_m \sum_n B_{mi} B_{mj} B_{ni} B_{nj}$$

in terms of matrices.

actually, there is also a condition: i is not equal to j.


migrated from mathematica.stackexchange.com Sep 24 '18 at 14:04

This question came from our site for users of Wolfram Mathematica.

  • $\begingroup$ This site is about Mathematica, the software. Should we move your question over to Mathematics? $\endgroup$ – Johu Sep 23 '18 at 0:38
  • $\begingroup$ Could you also check, if you have not made any mistakes with the indexes. $\endgroup$ – Johu Sep 23 '18 at 0:44
  • $\begingroup$ @DavidG.Stork I am not sure whether your edit was meaningfull. By swapping some of the indices, you actually changed the question dramatically. For example, Louis' (imho correct) answer became invalid under this modification. $\endgroup$ – Henrik Schumacher Sep 23 '18 at 7:41
  • $\begingroup$ @LouisB I'd suggest to undelete your post. Imho, it is a perfect answer to the original question. $\endgroup$ – Henrik Schumacher Sep 23 '18 at 7:53

If this is indeed a Mathematica question, then first note that:

$$\sum _j B_{i,j} B_{j,k}\equiv (B.B)_{i,j}$$


$$\sum _i B_{i,i}\equiv \operatorname{Tr}[B]$$

So, the Mathematica equivalent of:

$$\sum _i \sum _j \sum _m \sum _n B_{i,j} B_{j,m} B_{m,n} B_{n,i} $$


Tr[B . B . B . B]


Tr[MatrixPower[B, 4]]

For the original form of the question, note that:

$$ B_{i,j}\equiv B^T{}_{j,i} $$

So, the Mathematica equivalent of:

$$\sum _i \sum _j \sum _m \sum _n B_{m,i} B_{m,j} B_{n,i} B_{n,j}$$




The OP added the requirement that terms where $i=j$ should not be included. Without explaining why, you can use the following to compute this version:

$$\sum _i \sum _j \sum _m \sum _n B_{m,i} B_{n,i} B_{m,j} B_{n,j} \ \left(1-\delta _{i,j}\right)\equiv \operatorname{Tr}\left[B^T.B.B^T.B\right]-\ \operatorname{Tr}\left[\left(B^T.B\right)^2\right]$$

For your example, $B=\left( \begin{array}{cc} 1 & 2 \\ 3 & 4 \\ \end{array} \right)$, we have:

B = {{1, 2}, {3, 4}};
Tr[Transpose[B] . B . Transpose[B] . B] - Tr[ (Transpose[B] . B)^2 ]


  • $\begingroup$ thank you so much. firstly, subscripts "k" and "l" should change to "m" and "n", secondly, I checked your answer to my original question with an example, and I didn't find it correct. $\endgroup$ – farhadpti Sep 24 '18 at 14:21
  • $\begingroup$ @farhadpti I fixed the typo. How did you check the answer? $\endgroup$ – Carl Woll Sep 24 '18 at 14:24
  • $\begingroup$ consider B = [1 2; 3 4]. the answer of sigma form and matrix form are 771 and 892 respectively. $\endgroup$ – farhadpti Sep 24 '18 at 14:29
  • $\begingroup$ @farhadpti I get 892 for both. $\endgroup$ – Carl Woll Sep 24 '18 at 14:33
  • $\begingroup$ sorry, Yes you are right! I made a simple mistake in coding!! thank you so much. I updated the question. actually there is a condition that "i" should not be equal to "j". so some elements of B.transpose(B).B.transpose(B) must be removed. is it possible not to calculate them at all? what is the mathematical correct way for that? $\endgroup$ – farhadpti Sep 24 '18 at 14:43

Carl's answer is perfect as it is. His use of MatrixPower inspired me to think a little about computational complexity: Matrix-matrix multiplication has complexity $O(n^3)$ (in the straight-forward implementations, not these theoretically fancy but practically irrelevant algorithms). If we can get rid of one or two of them, we can even afford a transposition to speed up the computation:

n = 5000;
B = RandomReal[{-1, 1}, {n, n}];
a = Tr[B.B.B.B]; // RepeatedTiming // First
b = Tr[MatrixPower[B, 4]]; // RepeatedTiming // First
c = With[{A = B.B}, Total[A Transpose[A], 2]]; // RepeatedTiming // First
a == b == c





For the original question, LouisB's answer (c below) seems to be both correct and efficient:

n = 5000;
B = RandomReal[{-1, 1}, {n, n}];
a = Tr[B.Transpose[B].B.Transpose[B]]; // RepeatedTiming // First
b = With[{u = Flatten[B.Transpose[B]]}, u.u]; // RepeatedTiming // First
c = Total[(Transpose[B].B)^2, 2]; // RepeatedTiming // First
a == b == c






This is a closed matrix formulation of $B^4$.

  • 4
    $\begingroup$ It's the trace of $B^4$, no? $\endgroup$ – LouisB Sep 23 '18 at 2:00

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.