Solving new linear system that comes from an $p$ enrichment

Question

Let's say I am solving a simple Poisson problem using a Mixed (DG) finite element method. If we use orthogonal polynomials as basis functions we can write the finite-dimensional linear system as

$$ \left\lbrack\begin{array}{cc}I_1& S_{11} \\ -S_{11}^\mathrm{T} \end{array} \right\rbrack \left\lbrack\begin{array}{c} q_1\\u_1 \end{array} \right\rbrack = \left\lbrack \begin{array}{c} f_1\\g_1 \end{array} \right\rbrack, $$ for appropriately sized matrices $S_{11},I_1$ and right hand sides $f_1$ and $g_1$, Suppose I solve this linear system and obtain solutions $\hat{q}_1$ and $\hat{u}_1$.

Now consider the same underlying PDE with an enriched basis, polynomials of one degree higher, say. This will result in a new linear system with the following structure:

$$ \left\lbrack \begin{array}{cc} I_1& 0& S_{11}& S_{12}\\ 0&I_2& S_{21}&S_{22}\\ -S_{11}^\mathrm{T} & -S_{21}^{\mathrm{T}} & 0 &0\\ -S_{12}^\mathrm{T} & -S_{22}^\mathrm{T}& 0 & 0\\ \end{array} \right\rbrack \left\lbrack\begin{array}{cc}q_1\\q_2\\u_1\\u_2 \end{array} \right\rbrack= \left\lbrack\begin{array}{cc}f_1\\f_2\\g_1\\g_2 \end{array} \right\rbrack, $$

Where the vectors $q_2$ and $u_2$ contain all the new degrees of freedom.

My question are as follows,

1. How can we use the solution to the smaller problem to construct the solution of the larger problem?
2. Recursion: Suppose we have a fast way to compute the solution of the smaller problem (for all right hand sides) how can we use this method to find solution of the larger problem. I have a vague idea of some sort of preconditioner along the lines of multigrid.

Notes:

We can rewrite the large linear system as $$ \left\lbrack \begin{array}{cc} I_1& 0& S_{11}& S_{12}\\ 0&I_2& S_{21}&S_{22}\\ -S_{11}^\mathrm{T} & -S_{21}^{\mathrm{T}} & 0 &0\\ -S_{12}^\mathrm{T} & -S_{22}^\mathrm{T}& 0 & 0\\ \end{array} \right\rbrack \left\lbrack\begin{array}{cc}q_1-\hat{q}_1\\ q_2\\ u_1-\hat{u}_1\\ u_2\\ \end{array} \right\rbrack = \left\lbrack \begin{array}{c}0\\f_2\\0\\g_2 \end{array} \right\rbrack, $$ and after shuffling the order of the equations and unknowns get

$$ \left\lbrack \begin{array}{cccc} I_1& S_{11}&0& S_{12}\\ -S_{11}^\mathrm{T}& 0 &-S_{21}^\mathrm{T}& 0\\ 0& S_{21} &I_2&S_{22}\\ -S_{12}^\mathrm{T}&0 & -S_{22}^\mathrm{T} &0 \end{array} \right\rbrack \left\lbrack \begin{array}{c} q_1-\hat{q}_1\\ u_1-\hat{u}_1\\ q_2\\ u_2\end{array} \right\rbrack = \left\lbrack\begin{array}{c}0\\ 0\\ f_2\\g_2 \end{array} \right\rbrack. $$

The thing to notice here is that the original linear system appears in the upper left hand corner and the (presumably smaller) related linear system appears in the bottom right hand corner.

The $I$ matrices are really mass matrices, but since we chose orthogonal polynomials as the basis functions they are diagonal and therefore trivial to invert

I am especially interested in links to the literature or a name of a technique like this that I can use to search the literature. For the first question I know that we could use previous solution as a guess for an iterative solver but I am looking for something more clever than this.

Answer 1

What you are thinking of is something that uses the structure of the augmented matrix to make solution of the system simpler. For example, one could be tempted to think of forming the Schur complement with regard to the bottom right $2\times 2$ block of the rewritten system.

But I don't think anything like this exists. If it would, then you would have discovered a purely algebraic way to efficiently solve the linear systems that arise from high order discretizations. We don't have such algorithms, though: the best algorithms we have use multigrid, but these algorithms make very specific use of the properties and entries of the matrices we are trying to invert, as opposed to just their algebraic structure. In other words, I don't think that by just reordering terms you will be able to come up with an efficient solver.