LSMR¶
-
template<Number T, DataLocation L = Host>
class LSMR¶ Iterative solver for least-square problems using the LSMR algorithm [1]. It solves the system of linear equations \( \mathbf{A}\vec{x} = \vec{b} \). If the system is inconsistent, it solves the least-square problem
\[ \min_{\vec{x}} \| \mathbf{A}\vec{x} - \vec{b} \|_2 \]A
is a real rectangular matrix of dimension m-by-n andb
is a real vector of length m. The matrixA
may be dense or sparse (usually sparse).To improve the convergence rate of the solver, you can specify a right preconditioner \( \mathbf{P}^{-1} \), such that now LSMR solves the following least square problem
\[ \min_{\vec{x}} \| (\mathbf{AP})^\dagger \; \mathbf{P} \vec{x} - \vec{b} \|_2 \]More specifically, the Preconditioner object is a function that computes \( \vec{y} = \mathbf{P}^{-1} \mathbf{P}^{-T} \vec{x}\) for an input vector
x
and output vectory
.The code was adapted from the SciPy package (https://github.com/scipy/scipy.git).
References:
[1] D. C.-L. Fong and M. Saunders, ‘LSMR: An Iterative Algorithm for Sparse Least-Squares Problems’, SIAM J. Sci. Comput., vol. 33, no. 5, pp. 2950–2971, Jan. 2011, doi: 10.1137/10079687X.
- Template Parameters:
T – Data type used by the solver. Must be either
float
ordouble
.L – Location of the data.
Public Functions
-
LSMR(LSMROptions options = LSMROptions())¶
Creates a new LSMR solver.
- Parameters:
options – [in] solver parameters. See LSMROptions for more information.
-
virtual ~LSMR() = default¶
Default destructor.
-
inline LSMRMatResult operator()(MatrixMap<T, L> A, VectorArgs b, std::optional<VectorArgs> x0 = std::nullopt, std::optional<Preconditioner<T, L>> precond = std::nullopt)¶
Solves the least square problem
Ax = b
.- Parameters:
A – [in] coefficient matrix
b – [in] right-hand side of the equation
x0 – [in] initial guess for the solution (optional)
precond – [in] a Preconditioner to the system (optional). Suppose that \( \mathbf{L}^{-1}\) is the preconditioner of system, then
precond
must apply \( \mathbf{L}^{-1} \mathbf{L}^{-T}\) over a vectorx
.
- Returns:
An abstract object representing the solution of the linear system.
-
LSMRReport report()¶
- Returns:
Some relevant information about the last run. See LSMRReport for more information.
-
void set_params(LSMROptions options)¶
Set the parameters of the solver.
- Parameters:
options – [in] solver parameters. See LSMROptions for more information.
-
struct LSMROptions¶
Parameters of the LSMR solver.
Public Members
-
double damping_factor = 0¶
Damping factor \( \lambda \) for regularized least-squares problems. If \( \lambda \) is not
0
, LSMR solves the regularized least-squares problem::\[\begin{split} \min_x \| \begin{bmatrix} \vec{b} \\ 0 \end{bmatrix} - \begin{bmatrix} \mathbf{A}\vec{x} \\ \lambda \mathbf{I} \end{bmatrix} \|_2 \end{split}\]
-
double atol = 1E-6¶
Stopping criteria. If the system
Ax = b
seems to be consistent, LSMR terminates when \( \| \vec{r}\|_2 \leq atol * \| \mathbf{A}\vec{x} \|_2 + btol * \| \vec{b}\|\). Otherwise, it terminates when \( \left \| \mathbf{A}^H \vec{r}\|_2 \leq atol * \| \mathbf{A}\vec{r} \right \|_2 \).
-
index_t max_iter = kAuto¶
LSMR terminates when the number of iterations reaches
max_iter
. The default ismax_iter = min(m, n)
.
-
int verbose = 0¶
Set the verbosity level.
0
- No report.1
- Report only the residual\[ \| \mathbf{A}\vec{x} - \vec{b} \|_2 \]2
- Full iteration log.
-
double damping_factor = 0¶