I want to make a table in which each element is a matrix product:

<< Developer`
n = 2;
m = 6;

matrix = Map[ToPackedArray[#] &, Table[Sin[1.0*i*j], {i, 1, m}, {j, 1, n}, {k, 1, n}]];

vector = {IdentityMatrix[n][[1]]}\[Transpose];

Do[Flatten[Table[ #1\[ConjugateTranspose].(matrix[[#2]].#1) &[vector, i], {i, 1, m}]], {ite, 1, 2*10^5}] // AbsoluteTiming
(*{7.430010, Null}*)

In the code above there are two important parts:

First part = Choosing one matrix from the lists of matrices matrix and multiply it by vector.

Second part = Multiplying vector\[ConjugateTranspose] by the result of the first part.

I thought if I convert matrix[[#2]].#1, the first part, into one matrix multiplication, it would be more efficient as instead of many multiplications I calculate all of them in one go:

matrix2 = ToPackedArray[Flatten[Table[matrix[[i]]\[Transpose], {i, 1, m}], 1]];

matrix2 is the desired matrix to do the next step. If you use matrix//MatrixForm, then both the matrices look like each other.

Do[x = matrix2.vector; Flatten[Table[ #1\[ConjugateTranspose].(x[[#2 ;; #2 + n - 1,All]]) &[vector, i], {i, 1, n*m, n}]], {ite, 1, 2*10^5}] // AbsoluteTiming
(*{9.080013, Null}*)

I had to still choose different parts of (matrix2.#1),first part, to multiply that part with #1\[ConjugateTranspose],second part.

As you see although I calculate all the matrix multiplications in the first part at once, it is slower. I thought it might be because of this part of code [[#2 ;; #2 + n - 1, All]] in which I had to choose different parts of the result. I tried:

Do[x = matrix2.vector; Flatten[Table[ (x[[#1 ;; #1 + n - 1, All]]) &[i], {i, 1, n*m, n}]], {ite, 1, 2*10^5}] // AbsoluteTiming
(*{5.290007, Null}*)

Do[Flatten[Table[ (matrix[[#2]].#1) &[vector, i], {i, 1, m}]], {ite, 1, 2*10^5}] // AbsoluteTiming
(*{4.630007, Null}*)

At each iteration for matrix2 I just do the calculation once but for that of matrix at each iteration there are 6 multiplications. In the below I remove the table for matrix2 and replace it with Do for matrix:

Do[x = matrix2.vector, {ite, 1, 2*10^5}] // AbsoluteTiming
(*{0.240000, Null}*)

Do[ Do[(matrix[[#2]].#1) &[vector, i], {i, 1, m}], {ite, 1, 2*10^5}] // AbsoluteTiming
(*{2.620004, Null}*)

Now the calculations for matrix2 are faster. Thus, I believe that choosing parts of matrix2.vector causes the calculation gets slower.

How can I overcome this bottleneck? How can I choose parts of a matrix faster?


1 Answer 1


You can increase the performance by transposing your matrix and using Dot without explicit [[...]]

<< Developer`
n = 2;
m = 6;

matrix = Table[Sin[1.0*i*j], {i, m}, {j, n}, {k, n}];
matrix3 = ToPackedArray@Transpose@matrix;

vector = N@Transpose@{IdentityMatrix[n]};
vectorHC = ConjugateTranspose@vector;

Flatten[vectorHC.matrix3.vector, {{1, 3, 2, 4, 5}}] == 
 Flatten@Table[#1\[ConjugateTranspose].(matrix[[#2]].#1) &[vector, i], {i, 1, m}]
(* True *)

Do[Flatten[vectorHC.matrix3.vector, {{1, 3, 2, 4, 5}}], {ite, 10^5}] // AbsoluteTiming
(* {0.894130, Null} *)

Do[Flatten[Table[#1\[ConjugateTranspose].(matrix[[#2]].#1) &[vector, i], {i, m}]], 
  {ite, 10^5}] // AbsoluteTiming
(* {8.979821, Null} *)

Converting vector to a floating-point array (N[...]) also increase the performance a bit.


Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

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