Speeding up sums involving 16x16 matrices and 16x16x16x16 antisymmetric tensor

Question

I need to perform the following contraction involving four 16x16 matrices and a tensor $\theta$ that is 16x16x16x16.

T4 = Table[Simplify[Sum[VL[[m, m1]] VL[[n, n1]] VL[[p, p1]] VL[[q, q1]] θ4[[m1, n1,p1, q1]], {m1, 16}, {n1, 16}, {p1, 16}, {q1, 
 16}]], {m, 16}, {n,16}, {p, 16}, {q, 16}]

where \theta and VL are given as follows.

θ4=SparseArray[{{15, 3, 4, 5} -> g1, {15, 3, 5, 4} -> -g1, {15,4, 3, 5} -> -g1, {15, 4, 5, 3} -> g1, {15, 5, 3, 4} -> g1, {15, 5, 4, 3} -> -g1, {16, 6, 7, 8} -> g2, {16, 6, 8, 7} -> -g2, {16, 7, 6, 8} -> -g2, {16, 7, 8, 6} ->g2, {16, 8, 6, 7} -> g2, {16, 8, 7, 6} -> -g2, {15, 9, 10, 11} -> gt1, {15, 9, 11, 10} -> -gt1, {15, 10, 9, 11} -> -gt1, {15, 10, 11, 9} -> gt1, {15, 11, 9, 10} -> gt1, {15, 11, 10, 9} -> -gt1, {16, 12, 13, 14} -> gt2, {16, 12, 14, 13} -> -gt2, {16, 13, 12, 14} -> -gt2, {16, 13, 14, 12} -> gt2, {16, 14, 12, 13} -> gt2, {16, 14, 13, 12} -> -gt2},{16,16,16,16}];

VL={{Cosh[τ], 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, -Sinh[τ], 0}, {0, Cosh[ϕ], 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, -Sinh[ϕ]}, {0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}, {0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}, {0, 0, 0, 0,  1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}, {0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}, {0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0}, {0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0}, {0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0}, {0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0}, {0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0,  0}, {0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0}, {0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0}, {0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0}, {-Sinh[τ], 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, Cosh[τ], 0}, {0, -Sinh[ϕ], 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, Cosh[ϕ]}};

This should be a simple operation, but it takes ages (tens of minutes) and I have to abort the operation every single time. I don't know whether the difficulty is due to the large size of the matrices (but this shouldn't be since the matrices are sparse in entries). Thanks in advance for any possible help !

Answer 1

You can use TensorContract instead of Sum:

r = Activate @ TensorContract[
    Inactive[TensorProduct][VL, VL, VL, VL, θ4],
    {{2, 9}, {4, 10}, {6, 11}, {8, 12}}
]; //AbsoluteTiming

{0.419766, Null}

Using Activate with Inactive is a tip from @jose.

Answer 2

Using Dot[] appears to be about a factor of 3 faster than Carl Wolls solution (although I doubt that the timing for such short times is precise)

   vlT = Transpose[VL]
   r1 = Flatten[
   VL.Flatten[
     VL.\[Theta]4.vlT, {{2}, {1}, {4}, {3}}].vlT, {{2}, {1}, {4}, {3}}]; 

The timings on my computer are 0.000207 for the Dot[] solution and 0.0006595 for Carl's. Both results agree. (I converted VL to a SparseArray[] as well.)