3
$\begingroup$

Given \begin{equation} A_i=B+C_i \end{equation} where $A_i$,$B$ and $C_i$, $i=1,\dotsc,N$ are large square matrices, $B$ is symmetric, $C_i$ are zero matrices aside for a square block on the diagonal.

I would like to factor the $A_i$, preferably but not necessarily by $LU$ decomposition, but I don't want to do this repeatedly $N$ times, as $N$ may be large. Is there a method of doing this?

$\endgroup$
6
  • 1
    $\begingroup$ Are the nonzero blocks in $C_{i}$ much smaller than $B$ and always in the same position in the matrix? $\endgroup$ Oct 26, 2018 at 21:25
  • $\begingroup$ Is your intent to solve $A_i x = f_i$ for multiple $A_i$ and (possibly) multiple $f_i$? $\endgroup$
    – user7440
    Oct 26, 2018 at 22:48
  • $\begingroup$ The blocks are always much smaller than $B$; @user7440 yes, as part of time stepping $\endgroup$
    – DanielRch
    Oct 27, 2018 at 10:14
  • 1
    $\begingroup$ If you can guarantee $A_i$ and $B$ to be positive definite somehow, then you can use the product form Cholesky update of Goldfarb and Scheinberg. It's similar to the SMW formula, but is far more numerically stable. It is implemented in interior-point solvers like SeDuMi for example $\endgroup$ Oct 31, 2018 at 17:19
  • $\begingroup$ Conversely, if you can extend the Goldfarb and Scheinberg technique to indefinite symmetric matrices while guaranteeing stability, then that would be a nice contribution worth of publication in Math. Prog. Series A :-) $\endgroup$ Oct 31, 2018 at 17:20

2 Answers 2

4
$\begingroup$

The structure of $C_i$ gives a hope that $C_i$ can be described as a low-rank update. If that is the case (size of the square block is much smaller than the size of the matrix $B$), one can apply classical methods for low-rank updates based on Sherman-Morrisson-Woodbury identity.

While updating the LU factorization is not the easiest task (the required pivoting might be sensitive to the performed modification) it is certainly possible. To begin with, you might want to read the classic low-rank update section of QR in Golub, Van Loan (Section 12.6 in 2nd edition). Then, proceed with the intricacies of updating the $PA=LU$ factorization.

Note, since $B$ is symmetric, one can take advantage of Cholesky factorization whenever $C_i$ is also symmetric.

However, if $C_i$ is not (enough) low-rank, I don't think there is much hope on any savings here.

$\endgroup$
2
  • $\begingroup$ The SMW identity often leads to instabilities, do I have reason to believe it would be stable in this case? $\endgroup$
    – DanielRch
    Oct 29, 2018 at 19:05
  • 1
    $\begingroup$ @DanielRch you do not, in general. There are just not a lot of other options. $\endgroup$
    – Anton Menshov
    Oct 29, 2018 at 21:09
2
$\begingroup$
  • If you plan on solving $A_i x = f_i$, one alternative option is to use a Krylov iterative solver (like GMRES) with $B$ as preconditioner and $B^{-1} f_i$ as starting vector. Then you can prove that all the residuals will belong to the range of $C_i$. The iterations will be performed in the range of $C_i$ (which reduces the dimension of vectors from $N$ down to the rank of $C_i$).
  • For a reference with GMRES, you could check Section 4.2 of Heikkola, Rossi, Toivanen, "A parallel fictitious domain method for three-dimensional Helmholtz equation", SISC, vol. 24, pp. 1567-1588 (2003).
  • Finally, you could recycle the Krylov subspace at every step $i$, as long as $C_{i+1}$ and $f_{i+1}$ are relatively "close" to $C_{i}$ and $f_i$.
$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.