# Compute a double dot product between two tensors of rank 3 and 2

I would need help to calculate a double dot product between a rank 3 tensor A and a rank 2 tensor B (A:B) using mathematica.

Does someone know how to do that?

The double dot product is also known as the Frobenius inner product--in other words, it is the result of flattening the matrices and treating them as vectors.

So, here is another way to write it:

A = {
{{1, 2, 3}, {4, 5, 6}, {7, 8, 9}},
{{2, 0, 0}, {0, 3, 0}, {0, 0, 1}}
};
B = {{2, 1, 4}, {0, 3, 0}, {0, 0, 1}};

Flatten[A, {{1}, {2, 3}}].Flatten[B] (* -> {40, 14} *)


As far as I'm aware, for double ranked tensors, the double dot product is equal to:

$$A:B = \operatorname{Trace}( A \cdot B^T )$$

For this reason, a "hacked" solution to your problem would be

A = {{{1, 2, 3}, {4, 5, 6}, {7, 8, 9}}, {{2, 0, 0}, {0, 3, 0}, {0, 0,
1}}} ;
B = {{2, 1, 4}, {0, 3, 0}, {0, 0, 1}};
M = {, };
For[i = 1, i <= 2, i++, M[[i]] = Tr[A[[i]].Transpose[B]]]


Here I'm looping on the first index of A and then applying the previous formula. This works but it is not very elegant. You could probably automate the dimension of M instead of hardcoding it like I did.

The result will give M={40, 14}

• If this is the correct interpretation of the double dot product, a much more succinct version of your code is Tr[#.Transpose[B]] & /@ A. – march Sep 25 '15 at 5:51

Here is a straight-forward solution using TensorContract/TensorProduct :

A = {{{1,2,3},{4,5,6},{7,8,9}},{{2,0,0},{0,3,0},{0,0,1}}};
B = {{2,1,4},{0,3,0},{0,0,1}};

TensorContract[TensorProduct[A,B],{{2,4},{3,5}}]


{40, 14}

For large PackedArrays, the TensorProduct/TensorContract combo -- while being very expressive -- should be avoided as it easily eats up all memory (and it is slow, too). I found out that Tr[# B] & /@ A is a bit faster than Tr[#.Transpose[B]] & /@ A. However, it's even faster to use the following function:

cf = Compile[{{A, _Real, 2}, {B, _Real, 2}},
Total[Total[A B]],
RuntimeAttributes -> {Listable},
Parallelization -> True
]

A = RandomReal[{-1, 1}, {300, 400, 500}];
B = RandomReal[{-1, 1}, {400, 500}];
a = Tr[#.Transpose[B]] & /@ A; // AbsoluteTiming// First
b = Tr /@ (A.B\[Transpose]); // AbsoluteTiming// First
c = Total[Times[#, B], 2] & /@ A; // AbsoluteTiming// First
d = cf[A, B]; // AbsoluteTiming// First

(* 0.537068 *)
(* 0.524494 *)
(* 0.380716 *)
(* 0.079156 *)


Checking correctness:

Max[Abs[a - b]]
Max[Abs[a - c]]
Max[Abs[a - d]]

(* 1.13687*10^-13 *)
(* 6.82121*10^-13 *)
(* 5.68434*10^-13 *)


The functions Contract, multiDot from Exterior Differential Calculus and Symbolic Matrix Algebra perform contractions on nested lists.

It is convenient to think of an nth-level nested list as an nth-rank tensor. Contraction then produces lower rank tensors.

For use in the examples we define the following rank-3 and rank-4 tensors in three dimensions:

T3 = Array[t3, {3, 3, 3}];
T4 = Array[t4, {3, 3, 3, 3}];


The function Contract[x_,{i1,j1},{i2,j2},...] takes as arguments an nth-rank tensor x, and one or more lists of pairs of integers, indicating the positions of the indices to be contracted. It returns a tensor of rank n-2k, where k is the number of lists:

Contract[T3, {1, 2}]


(* {t3[1, 1, 1] + t3[2, 2, 1] + t3[3, 3, 1], t3[1, 1, 2] + t3[2, 2, 2] + t3[3, 3, 2], t3[1, 1, 3] + t3[2, 2, 3] + t3[3, 3, 3]} *)

The function multiDot[x_,y_,{i1,j1},{i2,j2},...] is similar to Contract, and generalizes the built-in function Dot. It takes as arguments two tensors -- x, y -- and one or more lists of pairs of integers. For each such pair {i, j}, it contracts the ith index of the first tensor (x) with the jth index of the second tensor (y), returning a tensor of rank m+n-2k, where m, n are the ranks of x, y, respectively, and k is the number of index pairs. (When the components of x and y are differential forms, multiDot multiplies them using Wedge).

 multiDot[T3, T3, {1, 1}, {2, 2}]//Dimensions


(* {3, 3} *)