fortran-lang/stdlib

Rank-revealing QR decomposition

loiseaujc opened this issue · 0 comments

Motivation

The stdlib_linalg module includes the qr decomposition. As far as I know, the stdlib_linalg_lapack module exposes interfaces for xGEQP3 which is the lapack routine for computing the QR with pivot which can be use as a rank-revealing factorization. I don't think it'll be a lot of work to extend the qr interface to allow for switching between the standard qr and the pivoting qr. The updated call to qr could look something like

call qr(A, Q, R, pivots)

where pivots is an integer(ilp), optional, intent(out) dummy argument (possibly allocatable). If pivots is not passed, the driver falls back to the regular qr factorization.

Prior Art

  • Currently, stdlib_linalg exposes qr(A, Q, R, overrite_a, storage, err) which computes only the regular QR factorization.

  • scipy exposes scipy.linalg.qr(a, overwrite_a=False, lwork=None, mode="full", pivoting=False, check_finite=True) for the QR factorization of A. Using pivoting=True returns Q, R, p where p is the vector of pivots.

  • The LinearAlgebra module in Julia exposes LinearAlgebra.qr(A, pivot=NoPivot(); blocksize) -> F where one can switch between regular QR factorization (pivot=NoPivot()) and pivoted QR (pivot=ColumnNorm()).

Additional Information

Pinging @perazz, @jvdp1 and @jalvesz.