Purpose
=======
LA_GELSY computes the minimum-norm least squares solution to one
or more real or complex linear systems A*x = b using a complete
orthogonal factorization of A. Matrix A is rectangular and may be
rankdeficient. The vectors b and corresponding solution vectors x are
the columns of matrices denoted B and X, respectively.
The routine computes a QR factorization of A with column pivoting:
A * P = Q * [ R11 R12 ]
[ 0 R22 ]
where R11 is the largest leading submatrix whose estimated condition
number is less than 1/RCOND. The order of R11, RANK, is the effective
rank of A. R22 is considered to be negligible, and R12 is annihilated
by orthogonal (unitary) transformations from the right, yielding the
complete orthogonal (unitary) factorization
A * P = Q * [ T11 0 ] * Z
[ 0 0 ]
The minimum-norm least squares solution is then
x = P * Z^H [ T11^-1 * Q1^H * b ]
[ 0 ]
where Q1 consists of the first RANK columns of Q.
=========
SUBROUTINE LA_GELSY( A, B, RANK=rank, &
JPVT= jpvt, RCOND= rcond, INFO= info )
(), INTENT(INOUT) :: A(:,:),
INTEGER, INTENT(OUT), OPTIONAL :: RANK
INTEGER, INTENT(INOUT), OPTIONAL :: JPVT(:)
REAL(), INTENT(IN), OPTIONAL :: RCOND
INTEGER, INTENT(OUT), OPTIONAL :: INFO
where
::= REAL | COMPLEX
::= KIND(1.0) | KIND(1.0D0)
::= B(:,:) | B(:)
Arguments
=========
A (input/output) REAL or COMPLEX array, shape (:,:).
On entry, the matrix A.
On exit, A has been overwritten by details of its complete
orthogonal factorization.
B (input/output) REAL or COMPLEX array, shape (:,:) with
size(B,1) = max(size(A,1),size(A,2)) or shape (:) with
size(B) = max(size(A,1), size(A,2)).
On entry, the matrix B.
On exit, rows 1 to size(A,2) contain the solution matrix X .
If size(A,1) >= size(A,2) and RANK = size(A,2), the residual
sum-of-squares for the solution vector in a column of B is
given by the sum of squares of elements in rows size(A,2)+1 :
size(A,1) of that column.
RANK Optional (output) INTEGER.
The effective rank of A, i.e., the order of the submatrix R11.
This is the same as the order of the submatrix T11 in the
complete orthogonal factorization of A.
JPVT Optional (input/output) INTEGER array, shape (:) with
size(JPVT) = size(A,2).
On entry, if JPVT(i) /= 0, the i-th column of A is an initial
column, otherwise it is a free column.
Before the QR factorization of A, all initial columns are
permuted to the leading positions; only the remaining free
columns are moved as a result of column pivoting during the
factorization.
On exit, if JPVT(i) = k, then the i-th column of the matrix
product A*P was the k-th column of A.
RCOND Optional (input) REAL.
RCOND is used to determine the effective rank of A. This is
defined as the order of the largest leading triangular
submatrix R11 in the QR factorization of A, with pivoting,
whose estimated condition number < 1/RCOND.
Default value: 10*max(size(A,1),size(A,2))*BEPSILON(1.0_),
where is the working precision.
INFO Optional (output) INTEGER.
= 0: successful exit
< 0: if INFO = -i, the i-th argument had an illegal value
If INFO is not present and an error occurs, then the program
is terminated with an error message.