## Numerical linear algebra: how to compute $b^TA^{-1}b$ efficiently

What is the most efficient way to compute $b^TA^{-1}b$ for a given $A$ and $b$?

Do we have to calculate $A^{-1}b$, or is this not necessary?

edit: I forgot to mention that A is symmetric and positive definite and sparse (so usually you'd use the conjugate gradient method).

What I have is a convex quadratic $x^TAx + b^Tx$. The minimum of this is at $2Ax+b=0$, and if you plug this minimum into the original form, then you get $x^T(-b/2)+b^Tx=b^Tx/2$ and this leads you to have to compute $-1/4\cdot b^TA^{-1}b$. So another way to pose the question is: can you find the height at the minimum faster than the location of the minimum?

-
Do we have to calculate $A^{−1}b$. Yes. Decompose $A$ with your favorite matrix decomposition, apply to $b$, and take dot product of the result with $b$. – J. M. Aug 17 2011 at 3:23
Thanks for your answer. I added extra information namely that A is symmetric and positive definite. What is the reason that this cannot be done more efficiently than first forming $A^{-1}b$? – Jules Aug 17 2011 at 3:30
Rough answer: this is "essentially" computing an element of the inverse, which is "essentially" computing a determinant (by the adjoint formula), and it is known that computing a determinant is as hard as a matrix product. Let me add that the complexity of matrix product is in practice $O(n^3)$ on floating-point numbers, and in theory $O(n^{2+\omega})$ for some $0<\omega<0.38$ whose exact value is still unknown. – Federico Poloni Aug 17 2011 at 7:33
Thanks! So computing the whole inverse is essentially not harder than computing just one entry of it. On the other hand this only lets you compute diagonal entries of a matrix with special structure, so perhaps there is a better method. – Jules Aug 17 2011 at 10:21
If you compute $b^T_1A^{-1}b_1$, $b^T_2A^{-1}b_2$, and $(b_1+b_2)^TA^{-1}(b_1+b_2)$, you can get by subtraction $b^T_1A^{-1}b_2+b_2^TA^{-1}b_1$, which will give you arbitrary entries of the inverse, not just diagonal ones. – Peter Shor Aug 17 2011 at 13:14
This is possibly an answer from a practical point of view: If you use the CG method for solving $x=A^{-1}b$ then $b^T A^{-1}b$ can be obtained along the way. However, it has been shown that computing $b^T A^{-1}b$ during the iteration can converge faster than first solving for $x$ and then multiplying $b^T x$. See "Z. Strakos and P. Tichy, On efficient numerical approximation of the bilinear form c*A-1b , SIAM Journal on Scientific Computing (SISC), 33, 2011, pp. 565-587" and the references therein for the positive definite case.