num::woodbury_solver

Defined in header <linear/factorization/woodbury.hpp>

template <factorization F>
class woodbury_solver;

Solves with \(A_{base} + PQ^{T}\) through the retained factor of \(A_{base}\).

The Woodbury identity gives \(x = y - WG^{-1}P^{T}y\) with \(y = A_{base}^{-T}b\), \(W = A_{base}^{-T}Q\) and \(G = I + P^{T}W\), and the mirror identity for forward solves. Construction costs one rank-p transpose solve. \(A_{base}^{-1}P\) and \(A_{base}^{-T}W\) are formed on first use, so a caller that only solves transposes never pays for them.

The base factorization is referenced, not copied, and must outlive this object.

Members

member	description
woodbury_solver(const F &base, low_rank_update update)
idx rank() const	The rank p of the update.
idx size() const	The order n of the corrected matrix.
const mat<real> &left() const	P.
const mat<real> &transpose_right() const	\(W = A_{base}^{-T}Q\).
const mat<real> &inverse_left() const	\(A_{base}^{-1}P\), formed on first use.
const mat<real> &transpose_right_squared() const	\(A_{base}^{-T}W\), formed on first use.
void right_solve(mat<real> &y, vec<real> &scratch) const	Replace y by \(yG^{-1}\), solving \(G^{T}y_i^{T} = y_i^{T}\) per row.
vec<real> inverse_diagonal(view<const real> base_diagonal) const	\(diag(A_{new}^{-1})\) from the base diagonal. \(diag(A_{new}^{-1}) = diag(A_{base}^{-1}) - diag(A_{base}^{-1}PG^{-1}Q^{T}A_{base}^{-1})\), which costs no solve beyond the two rank-p blocks already held.