Previous Next Contents Generated Index Doc Set Home



Solution to a Least Squares Problem for a Full-Rank General Matrix (Simple Driver)

The subroutines described in this section solve an overdetermined or underdetermined linear system AX = B or ATX = B by using a QR or LQ factorization of a general matrix A. The matrix A is assumed to have full rank. Note that the simple driver xGELSS and the expert driver xGELSX are also available for problems in which A is rank-deficient.

Calling Sequence

CALL DGELS 
(TRANSA, M, N, NRHS, DA, LDA, DB, LDB, DWORK, LDWORK, 
INFO)
CALL SGELS 
(TRANSA, M, N, NRHS, SA, LDA, SB, LDB, SWORK, LDWORK, 
INFO)
CALL ZGELS 
(TRANSA, M, N, NRHS, ZA, LDA, ZB, LDB, ZWORK, LDWORK, 
INFO)
CALL CGELS 
(TRANSA, M, N, NRHS, CA, LDA, CB, LDB, CWORK, LDWORK, 
INFO)






void dgels 
(char trans, int m, int n, int nrhs, double *da, int 
lda, double *db, int ldb, int *info)
void sgels 
(char trans, int m, int n, int nrhs, float *sa, int 
lda, float *sb, int ldb, int *info)
void zgels 
(char trans, int m, int n, int nrhs, doublecomplex *za, 
int lda, doublecomplex *zb, int ldb, int *info)
void cgels 
(char trans, int m, int n, int nrhs, complex *ca, int 
lda, complex *cb, int ldb, int *info)

Arguments

TRANSA

Determines the problem to solve. Legal values for TRANSA are shown below. Any values not shown below are illegal.

'N' or 'n'

Find a minimum norm solution to AX = B for M < N or minimize ||B-AX|| for M N.

'T' or 't'

Find a minimum norm solution to ATX = B for M N or minimize ||B-ATX|| for M < N.

M

Number of rows of the matrix A. M 0.

N

Number of columns of the matrix A. N 0.

NRHS

Number of right-hand sides, equal to the number of columns of the matrices B and X. NRHS 0.

xA

On entry, the matrix A.

On exit, if M N then A is overwritten by details of its QR factorization as computed by xGEQRF.

On exit, if M < N then A is overwritten by details of its LQ factorization as computed by xGELQF.

LDA

Leading dimension of the array A as specified in a dimension or type statement. LDA max(1, M).

xB

On entry, the matrix B of right-hand side vectors, stored columnwise. B is M×NRHS if TRANSA = 'N' or 'n'; B is N×NRHS if TRANSA =
'T' or 't'.
On exit, B is overwritten by the solution vectors, stored columnwise as shown below.

M, N

TRANSA

Matrix B

M N

'N' or 'n'

Rows 1 to N of B contain the least squares solution vectors. Residual sum of squares for the solution in each column is given by the sum of squares of elements N+1 to M in that column.

M < N

'N' or 'n'

Rows 1 to N of B contain the minimum norm solution vectors.

M N

'T' or 't'

Rows 1 to M of B contain the minimum norm solution vectors.

M < N

'T' or 't'

Rows 1 to M of B contain the least squares solution vectors. Residual sum of squares for the solution in each column is given by the sum of squares of elements M+1 to N in that column.

LDB

Leading dimension of the array B as specified in a dimension or type statement. LDB max(1, M, N).

xWORK

Scratch array with a dimension of LDWORK.

LDWORK

Leading dimension of the array WORK as specified in a dimension or type statement. LDWORK min(M, N) + max(1, M, N, NRHS).

INFO

On exit:

INFO = 0

Subroutine completed normally.

INFO < 0

The ith argument, where i = |INFO|, had an illegal value.

Sample Program




      PROGRAM TEST
      IMPLICIT NONE
C
      INTEGER           LDA, LDB, LDWORK, M, N, NRHS
      PARAMETER        (M = 4)
      PARAMETER        (N = 3)
      PARAMETER        (LDA = M)
      PARAMETER        (LDB = M)
      PARAMETER        (LDWORK = M + N)
      PARAMETER        (NRHS = 1)
C
      DOUBLE PRECISION  A(LDA,N), B(LDB,NRHS), WORK(LDWORK)
      INTEGER           ICOL, INFO, IROW
C
      EXTERNAL          DGELS
      INTRINSIC         ABS
C
C     Initialize the array A to store the coefficient matrix A
C     shown below.  Initialize the array B to store the right
C     hand side vector b shown below.
C
C         1   6   2         96
C     A = 1  -2  -8    b = 192
C         1  -2   4        192
C         1   6  14        -96
C
      DATA A / 1.0D0,  1.0D0,  1.0D0, 1.0D0,
     $         6.0D0, -2.0D0, -2.0D0, 6.0D0,
     $         2.0D0, -8.0D0,  4.0D0, 1.4D1 /
      DATA B / 9.6D1, 1.92D2, 1.92D2, -9.6D1 /
C
C     Print the initial value of the arrays.
C
      PRINT 1000
      PRINT 1010, ((A(IROW,ICOL), ICOL = 1, N), IROW = 1, M)
      PRINT 1020
      PRINT 1030, B
C
C     Compute and print the solution.
C
      CALL DGELS ('NO TRANSPOSE A', M, N, NRHS, A, LDA, B, LDB,
     $            WORK, LDWORK, INFO)
      IF (INFO .NE. 0) THEN
        PRINT 1040, ABS(INFO)
        STOP 1
      END IF
      PRINT 1050
      PRINT 1030, B
C
 1000 FORMAT (1X, 'A:')
 1010 FORMAT (3(3X, F8.4))
 1020 FORMAT (/1X, 'b:')
 1030 FORMAT (1X, F8.4)
 1040 FORMAT (1X, 'Illegal value for argument #', I2,
     $        ' in DGELS.')
 1050 FORMAT (/1X, 'min(b - Ax):')
C
      END
 

Sample Output

 
 A:
     1.0000     6.0000     2.0000
     1.0000    -2.0000    -8.0000
     1.0000    -2.0000     4.0000
     1.0000     6.0000    14.0000



 b:
  96.0000
 192.0000
 192.0000
 -96.0000



 min(b - Ax):
 148.0000
 -14.0000
  -8.0000
 -96.0000






Previous Next Contents Generated Index Doc Set Home