cupyx.scipy.sparse.linalg.lsmr#

cupyx.scipy.sparse.linalg.lsmr(A, b, x0=None, damp=0.0, atol=1e-06, btol=1e-06, conlim=100000000.0, maxiter=None)[source]#

Iterative solver for least-squares problems.

lsmr solves the system of linear equations Ax = b. If the system is inconsistent, it solves the least-squares problem min ||b - Ax||_2. A is a rectangular matrix of dimension m-by-n, where all cases are allowed: m = n, m > n, or m < n. B is a vector of length m. The matrix A may be dense or sparse (usually sparse).

Parameters:
  • A (ndarray, spmatrix or LinearOperator) – The real or complex matrix of the linear system. A must be cupy.ndarray, cupyx.scipy.sparse.spmatrix or cupyx.scipy.sparse.linalg.LinearOperator.

  • b (cupy.ndarray) – Right hand side of the linear system with shape (m,) or (m, 1).

  • x0 (cupy.ndarray) – Starting guess for the solution. If None zeros are used.

  • damp (float) –

    Damping factor for regularized least-squares. lsmr solves the regularized least-squares problem

    min ||(b) - (  A   )x||
        ||(0)   (damp*I) ||_2
    

    where damp is a scalar. If damp is None or 0, the system is solved without regularization.

  • atol (float) – Stopping tolerances. lsmr continues iterations until a certain backward error estimate is smaller than some quantity depending on atol and btol.

  • btol (float) – Stopping tolerances. lsmr continues iterations until a certain backward error estimate is smaller than some quantity depending on atol and btol.

  • conlim (float) – lsmr terminates if an estimate of cond(A) i.e. condition number of matrix exceeds conlim. If conlim is None, the default value is 1e+8.

  • maxiter (int) – Maximum number of iterations.

Returns:

  • x (ndarray): Least-square solution returned.

  • istop (int): istop gives the reason for stopping:

    0 means x=0 is a solution.
    
    1 means x is an approximate solution to A*x = B,
    according to atol and btol.
    
    2 means x approximately solves the least-squares problem
    according to atol.
    
    3 means COND(A) seems to be greater than CONLIM.
    
    4 is the same as 1 with atol = btol = eps (machine
    precision)
    
    5 is the same as 2 with atol = eps.
    
    6 is the same as 3 with CONLIM = 1/eps.
    
    7 means ITN reached maxiter before the other stopping
    conditions were satisfied.
    
  • itn (int): Number of iterations used.

  • normr (float): norm(b-Ax)

  • normar (float): norm(A^T (b - Ax))

  • norma (float): norm(A)

  • conda (float): Condition number of A.

  • normx (float): norm(x)

Return type:

tuple

References

D. C.-L. Fong and M. A. Saunders, “LSMR: An iterative algorithm for sparse least-squares problems”, SIAM J. Sci. Comput., vol. 33, pp. 2950-2971, 2011.