Let $D$ be an arbitrary diagonal matrix and let $A$ be an orthogonal matrix ($A'A = AA' = I$). How to compute the following matrix inverse efficiently?
$$(D + ADA^T)^{-1}$$
Hints or references are also welcome. Thanks a lot!
$\begingroup$
$\endgroup$
3
-
3$\begingroup$ note that $(D+ADA^T)^{-1}=A(DA+AD)^{-1}$; I would be surprised if $DA+AD$ has any structure that would simplify finding an inverse. $\endgroup$– Carlo BeenakkerCommented Aug 6, 2018 at 6:15
-
2$\begingroup$ We get a lot of questions about "update the inverse of a matrix after adding a diagonal correction" at scicomp.se (see e.g. scicomp.stackexchange.com/a/10282/4405), and the answer is almost invariably that there's nothing better than recomputing the inverse from scratch. There is a little more structure here, but my first guess would be no (unless $D=I$). $\endgroup$– Federico PoloniCommented Aug 6, 2018 at 6:36
-
$\begingroup$ Are $A'$ and $A^T$ the same thing? $\endgroup$– Federico PoloniCommented Aug 8, 2018 at 6:56
Add a comment
|
3 Answers
$\begingroup$
$\endgroup$
To begin with, the matrix in question can well be degenerate, consider for example
\begin{equation*} D=\left( \begin{array}{cc} 1 & 0 \\ 0 & -1 \end{array}\right), A=\left( \begin{array}{cc} 0 & -1 \\ 1 & 0 \end{array}\right). \end{equation*}
I suspect that because of this you should not count for much for a general $D$.
Things are better when $D>0$, because in this case $$D+ADA^T=D^{1/2}(I+BB^T)D^{1/2},\,B=D^{-1/2}AD^{1/2}.$$ Computing $(I+C)^{-1}$ where $C$ is a positive definite matrix with a known inverse $C^{-1}$ is a standard problem (computing it numerically, that is). As for good references, I simply do not remember (I am through with numerical analysis since a decade ago), but the wiki page https://en.wikipedia.org/wiki/Numerical_linear_algebra may be a good place to start.
$\begingroup$
$\endgroup$
I would like to add to Alex's idea. It is good to first make note of the fact that the given matrix which we are trying to invert is a symmetric matrix. Additionally, like Alex suggested having D>0 helps a lot, as, it actually makes the whole matrix to be positive definite (not only C (from Alex's answer)) (obviously! because of the additional Identity matrix) (actually, the fact that it is positive definite is more clearly visible if we keep it in the given original form). So, if D>0, the given matrix as whole is Symmetric Positive Definite (SPD) and, having it to be SPD makes things very comfortable (Both, to choose the right iterative method and also to prove the convergence of the method). At this point you have to ask if the matrix is Sparse or Dense (as you are asking about efficiency, I presume you are dealing with big matrix sizes). So, if it is Sparse then, Conjugate Gradient is the best to choose for SPD matrices (https://en.wikipedia.org/wiki/Conjugate_gradient_method) else, if the matrix is Dense, Cholesky decomposition is the way to go for SPD matrices and also for present structure in general as Federico's comment suggests (more efficient than LU decomposition) (https://en.wikipedia.org/wiki/Cholesky_decomposition). Additionally if the matrix is sparse, to make things more efficient in terms of storage while implementation, it is better to use sparse storage data structures (see https://en.wikipedia.org/wiki/Sparse_matrix).
(PS: the suggestions are from the well established methods, you can additionally look for some preconditioning which helps a lot for little effort! (see https://en.wikipedia.org/wiki/Preconditioner), Multigrid could a nice choice I believe for an SPD matrix for preconditioning.)
$\begingroup$
$\endgroup$
An answer to your question heavily depends on whether your matrix is sparse or dense. If it is dense (or even sparse, but not larger than a few thousands rows/columns), than any direct method, like based on LU-factorisation, will do. If it is sparse, the problem becomes a lot more complex. I suggest you visit the following blog page and ask any questions there: http://comecau.blogspot.com/2018_09_05_archive.html
