Efficient diagonal update of matrix inverse

Question

I am computing $(kI + A)^{-1}$ in an iterative algorithm where $k$ changes in each iteration. $I$ is an $n$-by-$n$ identity matrix, $A$ is an $n$-by-$n$ precomputed symmetric positive-definite matrix. Since $A$ is precomputed I may invert, factor, decompose, or do anything to $A$ before the algorithm starts. $k$ will converge (not monotonically) to the sought output.

Now, my question is if there is an efficient way to compute the inverse that does not involve computing the inverse of a full $n$-by-$n$ matrix?

Answer 1

Turns out this was very simple. I'll just write it here if someone else has need for it: $$(kI+A)^{-1}=(kI+PDP^{-1})^{-1}=(P(D + kI)P^{-1})^{-1}=P(D + kI)^{-1}P^{-1}$$ where $A=PDP^{-1}$ is the eigenvalue decomposition. And the inverse of a diagonal matrix is quickly computed as the matrix with diagonal elements the reciprocal of the diagonal elements of the original matrix. This is roughly $50$ times faster than the original solution for my data on my computer.

-- Tommy L

Answer 2

EDIT. 1. The Tommy L 's method is not better than the naive method.

Indeed, the complexity of the calculation of $(kI+A)^{-1}$ is $\approx n^3$ blocks (addition-multiplication).

About the complexity of $P(D+kI)^{-1}P^{-1}=QP^{-1}$ (when we know $P,D,P^{-1}$); the complexity of the calculation of $Q$ is $O(n^2)$ and the one of $QP^{-1}$ is $\approx n^3$ blocks as above. Note that one works with a fixed number of digits.

2. When the real $k$ can take many (more than $n$) values, one idea is to do the following

The problem is equivalent to calculate the resolvent of $A$: $R(x)=(xI-A)^{-1}=\dfrac{Adjoint(xI-A)}{\det(xI-A)}=$

$\dfrac{P_0+\cdots+P_{n-1}x^{n-1}}{a_0+\cdots+a_{n-1}x^{n-1}+x^n}$.

STEP 1. Using the Leverrier iteration, we can calculate the matrices $P_i$ and the scalars $a_j$

$P_{n-1}:=I:a_{n-1}:=-Trace(A):$

for $k$ from $n-2$ by $-1$ to $0$ do

$P_k:=P_{k+1}A+a_{k+1}I:a_k:=-\dfrac{1}{n-k}Trace(P_kA):$

od:

During this step, we must calculate the exact values of the $P_i,a_j$ (assuming that $A$ is exactly known), which requires a very large number of digits. Then the complexity of the calculation is $O(n^4)$ (at least) but it's done only one time.

STEP 2. We put $x:=-k_1,-k_2,\cdots$. Unfortunately, the time of calculation -with a fixed number of significative digits- of $R(-k)$, is larger than the time of calculation of $(-kI-A)^{-1}$ !!!