# Schur complement formulation of linear system

Consider a system of the following form:

$$(A+K)x=b$$

where $$A$$ is symmetric, positive definite and block diagonal (in fact, a block diagonal matrix made of stiffness matrices arising from FEM discretizations at multiple timesteps, grouped together). Multiplying our linear system by $$P^t=K^tA^{-1}+I$$, we get $$(K+K^t + K^tA^{-1}K+A)x = P^tb$$

We assume that $$K$$ is non-symmetric, has full rank and that $$K+K^t$$ is positive definite, so that this is a symmetric and positive definite system of equations. The only inconvenient here is the presence of $$A^{-1}$$, and the fact that one needs to solve a (block diagonal) system to obtain $$P^tb$$.

Yet another way is to further rewrite everything recognizing a Schur complement:

$$\left\lbrack \begin{array}{cc} A& -K^t\\ -K^t& -K-K^t-A \end{array} \right\rbrack \left\lbrack \begin{array}{c} p\\ x\end{array} \right\rbrack= \left\lbrack \begin{array}{c}-b\\ -b\end{array}\right\rbrack$$

Which one of the three systems would you rather work with, and with which solver would you work? I am interested in wall-clock time and parallelizability

We assume our blocks to be large, e.g. of size $$1e8-1e10$$.

To connect the first two formulations, an idea might be to solve formulation one with pre-conditioned GMRES, and use $$P^t$$ as preconditioner. Applying $$P^t$$ can be done with the (preconditioned) conjugate gradient, and this is highly parallelizable since $$A$$ has a diagonal block structure.

But in this way, since we are solving essentially formulation two, we'd still need to come up with a good preconditioner for it.

Edit. Also because formulation two seems much worse conditioned than formulation one, as numerical experiments show.

• Could you elaborate and write out the structure of Pt? Commented Mar 16, 2023 at 13:53
• @rchilton1980of course, I edited the question Commented Mar 16, 2023 at 14:05
• Guessing K is not symmetric? Commented Mar 16, 2023 at 14:07
• No, $K$ has a block structure like the time matrix from the implicit Euler method. And the blocks are not constant, and not symmetric. @rchilton1980 Commented Mar 16, 2023 at 14:52