# Fastest way to do vector / matrix multiplication with constant matrix

Suppose we have a fixed n x n matrix, call it $B$. What is fastest way to evaluate $v^t B\;v$ over many different vectors $v$?

Since $B$ is a constant matrix, does this allow the number of operations needed to evaluate $v^t B\;v$ for an arbitrary $v$ to be reduced?

-
You might want to learn the way Mathematica interfaces with CUDA. In particular, read about the dot product and loading to memory. If you load the matrix to CUDA memory, it will be the fastest you can go. –  Hector Oct 1 '13 at 11:06
Floating point matrix multiplication is done using the MKL in Mathematica, and it is already very fast. It's much faster than a naive C implementation, especially for large matrices. What you can do to speed things up is try reduce the overhead of Mathematica before it calls the MKL, e.g. make sure that B is a packed array. Otherwise it's going to be as fast as it can get, unless the matrices are tiny. Whether CUDA helps depends on your hardware. –  Szabolcs Oct 1 '13 at 13:38
Loading the vectors into the GPU memory might be slower than just doing the math on the CPU. –  s0rce Oct 1 '13 at 13:49
You can't get away without at least $N^2$ multiplications. There could be some simplifications if your $B$ is symmetric or if $v$ is a random vector. Perhaps if you described what you were trying to achieve, then we could help you with that. –  rm -rf Oct 1 '13 at 14:13
Assuming that v is a vector of vectors, v == {v1, v2, ...}, the direct way is this: Table[x.b.x, {x, v}]. This is quite fast. We can improve a bit by doing part of the multiplication in a single go: r = v.b; Table[v[[i]].r[[i]], {i, Length[v]}]. This is a bit faster. I don't know how to do the Table part in one go while avoiding unpacking the arrays. –  Szabolcs Oct 1 '13 at 14:23

If the vectors v come one at a time and must be handled that way, as opposed to being collected for batch analysis, then I don't see how you can beat simply v.B.v. However, if you can collect them as the rows of a many x n matrix V then Total[V.B*V,{2}] is very fast. (If the vectors come "naturally" as the columns of an n x many matrix U then use Total[U*B.U].) However, there may be version and/or platform differences, so YMMV.
On my machine with Mathematica V9 this is indeed nearly one order of magnitude faster than the Table approach. –  sebhofer Oct 2 '13 at 9:47