Solving Generalization of Saddle point problem

Question

I am interested in knowing if there is a generalization of the Uzawa iteration for the linear problems of the form

$$\left[ \begin{array}{ccc}A& B^T&0\\ B&0&C^T\\ 0&C&0 \end{array} \right]\left[\begin{array}{c}x\\ y\\ z\end{array} \right] = \left[\begin{array}{c} f\\ g\\ h \end{array} \right]$$

Here $A$, is symmetric and positive definite, and in my case it is diagonal so trivial to invert. $A\in \mathbb{R}^{n\times n}$, $B\in \mathbb{R}^{m\times n}$, and $C\in \mathbb{R}^{l\times m}$.

With $n \approx 3m $ and $l \ll m$.

For some background, I need to solve this system for every step in a time simulation. I hope to use the previous time step's solution as a good guess for the new solution, perhaps after extrapolating.

I have a way to solve this using a triple nested iteration, but I am interested in something like the Uzawa iteration because I feel I have a good starting guess.

Any help would be greatly appreciated. Do these matrices have a name?

Answer 1

Looking at your equations, you see that $x=A^{-1}(f-B^Ty)$. Since $A$ is diagonal and invertible, you can easily transform your linear system by a block-Gauss elimination:

$$\left[ \begin{array}{cc}-BA^{-1}B^T& C^T\\ C&0 \end{array} \right]\left[\begin{array}{c}y\\ z\end{array} \right] = \left[\begin{array}{c} g - A^{-1}f\\ h \end{array} \right].$$

This is now only a problem in $y,z$, but $x$ can always be recovered via $x=A^{-1}(f-B^Ty)$ which is again cheap to compute. Now the problem is in saddle-point form, which allows the use of Uzawa-type techniques. Let us denote $S_1=BA^{-1}B^T$. Since $B$ has full column rank, the matrix $S_1$ is symmetric and positive definite. Using that $y = S_1^{-1}(C^Tz+A^{-1}f-g)$, we can now form a second Schur complement as

$$\left[ \begin{array}{cc}-S_1& C^T\\ 0 & S_2 \end{array} \right]\left[\begin{array}{c}y\\ z\end{array} \right] = \left[\begin{array}{c} g - A^{-1}f\\ h + S_1^{-1}(g-A^{-1}f) \end{array} \right],$$

where $S_2 = CS_1^{-1}C^T$ which is as well symmetric and positive definite, since $C$ has full column rank. Now you can use your favorite solver to solve for $z$, then recover $y$ and finally $x$. Note that while $S_1$ can still be explicitly computed for practical purposes, $S_1^{-1}$ and therefore $S_2$ will probably be dense. Therefore it might be good to use a nested Krylov subspace method (e.g. conjugate gradients). In this way, both, the inner and outer iterations only require multiplications with the respective Schur-complements, thereby avoid their direct assembly. For instance an application of $S_2$ would be decomposed into applying $C^T$ then an (inner) solve with $S_1^{-1}$ and another application of $C$. The same can be done for $S_1$ if you don't want to explicitly compute it.

Let me finally remark that in a practical implementation of Uzawa-type iterations it is possible to efficiently update $y$ alongside $z$ and thereby obtain an efficient algorithm. Also inexact Uzawa approaches and suitable preconditioners for your subproblems can possibly be considered.