Add a function to compute trace of matrix product efficiently

[jachymb]

I figured something like this could be useful (although not strictly necessary) pre-requisite for #3149

Consider rectangular matrices $A\in {\mathbb{R}}^{m\times n},B\in {\mathbb{R}}^{n\times m}$. We can observe the following:

$\mathrm{Trace}(A\cdot B)=\sum _i{A}_i\cdot (B^T{)}_i=\sum _i\sum _j{A}_{i,j}{B}_{j,i}$

The left-hand-side expression is elegant and concise way to write it. From programmers perspective, it's convenient especially when we already have matrices on input.

The middle expression seems like a good way to compute the value in a vectorized way. It's asymptotically faster than computing the whole product and then taking the trace.

The right-hand-side expression is just expanded form and may appear in the wild. It suggests computing it as sum of elementwise product of $A$ and $B^T$ which is maybe better computationally, but kinda obscures the linear-algebraic structure.

I thus suggest we add a function like real trace_dot(matrix a, matrix b). Note that we already have functions like trace_quad_form that do something related.

Not sure if the property $\mathrm{Trace}(A\cdot B)=\mathrm{Trace}(B\cdot A)$ can be useful here, but it allows us to choose (in the middle formula) to do either more dot products of shorter vectors or fewer dot-products of longer vectors. It should amount to the same number of multiplications but perhaps for hardware acceleration one of those is preferred.

For autodiff, I think there is

$\frac{\mathrm{d}}{\mathrm{d}A}\mathrm{Trace}(A\cdot B)=B^T$