Speaker
Description
Consider saddle-point systems of the following form
$$ \hspace{5.5cm}
\mathcal{A}
\left( \begin{array}{c} u \\ v \end{array} \right) =
\left(
\begin{array}{cc}
A & B^T \\ B & 0 \\
\end{array}
\right) \left( \begin{array}{c}
u \\ v \end{array} \right) = \left(
\begin{array}{c} f \\ g \end{array} \right)
\hspace{5.5cm} (1)
$$
where $A \in \mathbb{R}^{n \times n}$ is large, sparse and symmetric positive definite and $B \in \mathbb{R}^{k \times n} (n > k)$ has full rank. Such systems arise in various scientific applications, including the finite element approximation of PDEs in which a constraint is enforced on a subset of unknowns via a global scalar Lagrange multiplier, numerical continuation methods for large nonlinear systems of equations, constrained optimization as well as in linear least squares problems [2].
Null-space methods represent one of possible solution approaches to solve such systems. Despite these solution approaches that compute a null-space $Z \in \mathbb{R}^{n \times (n-k)}$ basis of the undetermined block $B$ are traditional in various engineering communities and offer advantages in some situations, general techniques to obtain such bases useful for large and sparse $A$ have not been studied frequently in the past. There are more reasons for this. One of them may be that the resulting transformed system with the matrix $Z^TAZ$ is often dense, or ill-conditioned, or both. Another reason may be that there is a lack of efficient preconditioners for the transformed systems even when the projected matrix $Z^T AZ$ is reasonably sparse. The core of the talk is devoted to discussing various algorithms to find the null-space basis $Z$. Although in several of the approaches treatment of $B$ that does not have full-rank is not difficult, we do not treat such cases at length.
Initial motivation for studying the null-space approaches in this talk is slightly different from the symmetric saddle-point problem given above. Consider the linear least-squares (LS) problem
$$
\min_x \| Ax - b \|_2,
$$
where the system matrix $A \in \mathbb{R}^{m \times n} (m \geq n)$ and the right-hand side $b \in \mathbb{R}^m$ are given. Assume that a part of the rows of $A$ is dense and the rows have been initially permuted such that the dense rows are the last rows of $A$. If we assume a conformal partitioning of the vector $b$ we have
$$
A = \left( \begin{array}{c} A_s \\ A_d \end{array} \right),
\
A_s \in \mathbb{R}^{m_s \times n}, \
A_d \in \mathbb{R}^{m_d \times n}, \
b = \left( \begin{array}{c} b_s \\ b_d \end{array} \right), \
b_s \in \mathbb{R}^{m_s}, \
b_d \in \mathbb{R}^{m_d},
$$
where $m_s$ denotes the number of sparse rows of $A$, and $m_d$ is the number of dense rows, with $m = m_s + m_d$, $m_s \geq n > m_d \geq 1$ and $m_d \ll m_s$.
The normal equations that describe the solution [3] (not representing always the best way to get the solution) become
$$\hspace{5cm}
Cx = (C_s + A^T_d A_d)x = c, \
c = A^T_s b_s + A^T_d b_d,
\hspace{5cm} (2)
$$
where $C_s = A^T_s A_s$ is the reduced normal matrix. The solution of (2) can be obtained from the equivalent $(n + md) \times (n + md)$ system
$$ \hspace{6.3cm}
\left(
\begin{array}{cc}
C_s & A_d^T \\ A_d & -I \\
\end{array}
\right) \left( \begin{array}{c}
x \\ A_d x
\end{array} \right) = \left( \begin{array}{c}
c \\ 0
\end{array}\right)
\hspace{6.3cm} (3)
$$
This is the system nearly of the form (1) with $k = m_d, H = C_s, B = A_d$. But the $(2,2)$ block is nonzero and equal to the negative unit matrix of the appropriate dimension. Provided $A_s$ has full column rank, $C_s$ is symmetric positive definite. But the null-space based approach can be extended even if the $(2, 2)$ block is nonzero and $C_s$ only positive semidefinite. This positive definiteness is quite common in the LS problems even when the whole $A$ is of full column rank. The approach to solve (3) also in this case for the full rank LS problem (3) is motivated by the following simple lemma.
\textbf{Lemma}: Consider $A = A_s A_d,\ As \in \mathbb{R}^{m_s \times n}, A_d \in \mathbb{R}^{m_d \times n}$. If $A$ is of full column rank, then $C_s = A^T_s A_s$ is positive definite on the null space of $A_d$.
Success of any null-space approach then depends on constructing a suitable null-space basis that keeps the matrix $Z^T AZ$ sufficiently sparse and well-conditioned. Of course, it cannot be an orthogonal null-space basis computed, for example, by pivoted $QR$ factorization of $B^T$. We will
discuss other approaches to find the null-space basis heavily based on the fact that $B$ is wide.
That is, it has far fewer rows than columns. The approaches include computing an oblique basis $Z$ computed from more separate QR factorizations and an approach that adds rows of $B$ one by one and extends the basis incrementally. Another approach for constructing $Z$ is the right oblique
conjugation. Applied to $B \in \mathbb{R}^{k \times n}$ it yields $V \in \mathbb{R}^{n \times n}$ and lower trapezoidal $L \in \mathbb{R}^{k \times n}$ satisfying
$$ \hspace{8cm}
BV = L.
\hspace{8cm} (4)
$$The null-space basis is then formed by the last $n-k$ columns of $V$ . This approach is directly related to Thesis of Michele Benzi [1]. Techniques implied by this Thesis have motivated the right oblique conjugation, and they may still offer further computational possibilities. Not only for constructing $Z$ but possibly also for preconditioning of the linear systems with the matrix $Z^T AZ$.
Our experimental results show that the null-space bases obtained from the mentioned approaches are of high quality and can push forward further research in optimization and solving constrained least squares. Linear least squares problems that contain a small number of dense rows arising from practical applications are used to illustrate our ideas and to explore their potential for solving large-scale systems. Our work was originally motivated by the paper by Howell [4]. A significant part of the obtained results has been published in [5].
References
[1] M. Benzi, A Direct Row-Projection Method For Sparse Linear Systems, Ph.D. Thesis, NCSU, 1993.
[2] M. Benzi, G.H. Golub, and J. Liesen, Numerical solution of saddle point problems, Acta Numerica,
14:1–137, 2005.
[3] ̊A. Bj ̈orck, Numerical Methods for Least Squares Problems, SIAM, Philadelphia, 1996.
[4] J. S. Howell, Prestructuring sparse matrices with dense rows and columns via null space methods,
Numerical Linear Algebra with Applications, 25:1–30, 2018.
[5] J. A. Scott and M. T ̊uma, A null-space approach for large-scale symmetric saddle point systems with
a small and non zero (2,2) block, Numerical Algorithms, published online, 2022.
