# Iterative single variable solutions in large linear systems

I have a system where $A$ is a large $n\times n$ marix with fast MVMs. It may have many nonzero entries (albeit in a structured way so as to allow fast MVMs), and is not necessarily diagonally dominant.

However, $A$ is weakly positive definite.

$Y$ is a $n\times m$ dense matrix.

I understand iterative Krylov-subspace-based methods exist for finding $X=A^{-1}Y$, and they perform well. Are there any optimizations that can be implemented if I am only interested in finding: $$\textbf{x}=\text{diag }Y^\top A^{-1}Y$$ In other words, I only want to solve for the $i$-th entry of the solution $Y^\top\textbf{x}_i$ where $A\textbf{x}_i=\textbf{y}_i$ for $Y=[\textbf{y}_i]_i$.

MVM stands for matrix-vector multiplication. In my application this takes only ~$n$ runtime compared to $n^2$ when dense.

• You might want to spell out what "MVM" means. Feb 5 '17 at 7:31
• @WolfgangBangerth done
– VF1
Feb 5 '17 at 13:13
• As far as I know, no. Feb 6 '17 at 20:11
• I believe the general case is still a research question at this point. I can think of several applications in statistics, ML, and physics that will immediately benefit from an efficient solution. Mar 7 '17 at 18:15
• Do you know anything specific about the eigenvalues of the matrix $A$? If they are clustered, then we may use Krylov methods to construct a polynomial approximation $B = p(A)$ minimizing $\|BA - I \|_2$ that may be used as a cheap surrogate for the matrix inverse. Otherwise, Krylov methods may not be of much use.... Mar 7 '17 at 18:19

An approach by Papandreou and Yuille for diagonal $$A$$ relates variance estimation to the expectation of a quadratic form. The logic follows more generally: since $$A$$ is PD so is its inverse. Then $$Z\sim N(0, A^{-1})$$ is well-formed, and $$\textbf{x}=\mathbb{E}[(Y^\top Z)^2]$$.
If we can sample $$Z$$ then an approximation is possible (and converges). The linked paper has success when we can find the square root of $$A$$: solving $$AZ=A^{1/2}G$$ for samples $$G\sim N(0,I)$$ offers a way of sampling $$Z$$ since $$A^{-1/2}G\sim N(0,A^{-1})$$.
When an operator $$A^{1/2}$$ is available, this provides an estimation opportunity for variance.
• How do you plan to solve $AZ=A^{1/2}G$, when the entire procedure has to be cheaper than solving a linear system with $A$, otherwise there would be no point in using it? Feb 7 '17 at 7:22
• Right, I'm saying this is an approach only when you have cheap access to $A^{1/2}$ (say, $A$ is a large Kronecker product). If so, then you're only solving one linear system, $AZ=A^{1/2}G$, where $Z,A^{1/2}G$ are vectors. This is definitely better than the multiple linear systems required for $A^{-1}Y$, where $Y$ is a matrix. That's why I didn't accept my answer. Of course, this is only cheaper if the number of samples of $G$ that you make is less than the width of $Y$.