HSL_MI27: Projected preconditioned conjugate gradient method for saddle-point systems

This package uses the projected preconditioned conjugate gradient method to solve \((n+m)\times (n+m)\) saddle-point systems of the form \[\left[ \begin{array}{cc} \mathbf{A} & \mathbf{B}^T \\ \mathbf{B} & -\mathbf{C} \end{array} \right] \left[ \begin{array}{c} \mathbf{x}\\\mathbf{y} \end{array} \right] = \left[ \begin{array}{c} \mathbf{C} \\ \mathbf{d} \end{array} \right],\] where \(\mathbf{A}\) is an \(n\times n\) real and symmetric matrix, \(\mathbf{C}\) is an \(m\times m\) real, symmetric and positive semi-definite (possibly zero) matrix, and \(m\leq n.\) A preconditioner of the form \[\mathbf{P} = \left[ \begin{array}{cc} \mathbf{G} & \mathbf{B}^T \\ \mathbf{B} & -\mathbf{C} \end{array} \right]\] must be available where \(\mathbf{G}\) is a real and symmetric matrix. The following assumptions are assumed to hold:

• if \(\mathbf{C}\) is positive definite, both \(\mathbf{A + B}^T\mathbf{C}^{-1} \mathbf{B}\) and \(\mathbf{G + B}^T\mathbf{C}^{-1} \mathbf{B}\) are positive definite;

• if \(\mathbf{C} =0\) and the columns of the \(n\times (n-m)\) matrix \(\mathbf{Z}\) span the nullspace of \(\mathbf{B},\) both \(\mathbf{Z}^T\mathbf{ A Z}\) and \(\mathbf{Z}^T\mathbf{ G Z}\) are positive definite;

• if \(\mathbf{C}=\mathbf{E D E}^{\mathit{T}},\) where \(\mathbf{d}\) is a \(p\times p\) nonsingular matrix with \(0<p<m,\) the columns of the \(m\times (m-p)\) matrix \(\mathbf{F}\) span the nullspace of \(\mathbf{C}\) and the columns of the \(n\times (n-m+p)\) matrix \(\mathbf{Z}\) span the nullspace of \(\mathbf{F}^{\mathit{T}} \mathbf{B},\) then both \(\mathbf{Z}^{\mathit{T}}\left( \mathbf{A + B}^{\mathit{T}} \mathbf{E D}^{-1} \mathbf{E}^{\mathit{T}} \mathbf{B} \right)\mathbf{Z}\) and \(\mathbf{Z}^{\mathit{T}}\left( \mathbf{G + B}^{\mathit{T}} \mathbf{E D}^{-1} \mathbf{E}^{\mathit{T}} \mathbf{B} \right)\mathbf{Z}\) are positive definite.

If these assumptions do not hold, then negative curvature may occur and, consequently, the method terminates with an error.

The projected preconditioned conjugate gradient method iteratively finds the vector \(\mathbf{x}\) and then, once \(\mathbf{x}\) has been computed to a high enough level of accuracy, the vector \(\mathbf{y}\) is computed by performing one additional solve with the preconditioner \(\mathbf{P}.\) Reverse communication is used for preconditioning and matrix-vector products of the form \(\mathbf{A s},\) \(\mathbf{B s},\) \(\mathbf{B^{\mathit{T}} s}\) and \(\mathbf{Cs}.\) HSL_MI13 may be used to efficiently form suitable preconditioners and carry out the required preconditioning solves; HSL_MC65 may be used to form the required matrix-vector products. HSL_MI13 and HSL_MC65 are both available as part of the current version of HSL.