sgeco(l) — SunSoft Performance Library

NAME

sgeco - compute the LU factorization and estimate the condition number of a general matrix A. If the condition number is not needed then xGEFA is slightly faster. It is typical to follow a call to xGECO with a call to xGESL to solve Ax = b or to xGEDI to compute the determinant and inverse of A.

SYNOPSIS

CALL DGECO (DA, LDA, N, IPIVOT, DRCOND, DWORK)

CALL SGECO (SA, LDA, N, IPIVOT, SRCOND, SWORK)

CALL ZGECO (ZA, LDA, N, IPIVOT, DRCOND, ZWORK)

CALL CGECO (CA, LDA, N, IPIVOT, SRCOND, CWORK)

ARGUMENTS

xAOn entry, the matrix A.
On exit, an LU factorization of A.

LDALeading dimension of the array A as specified in a dimension or
type statement. LDA >= max(1,N).

NOrder of the matrix A. N >= 0.

IPIVOTOn exit, a vector of pivot indices.

xRCONDOn exit, an estimate of the reciprocal condition number of A.
0.0 <= RCOND <= 1.0. As the value of RCOND gets smaller, operations with A such as solving Ax = b may become less stable. If RCOND satisfies RCOND + 1.0 = 1.0 then A may be singular to working precision.

xWORKScratch array with a dimension of N.

SAMPLE PROGRAM

      PROGRAM TEST
      IMPLICIT NONE
C
      INTEGER           IAXEQB, LDA, LDB, N
      PARAMETER        (IAXEQB = 0)
      PARAMETER        (N = 3)
      PARAMETER        (LDA = N)
      PARAMETER        (LDB = LDA)
      DOUBLE PRECISION A(LDA,N), B(LDB), RCOND, WORK(N)
      INTEGER           ICOL, IPIVOT(N), IROW, JOB
C
      EXTERNAL          DGECO, DGESL
C
C     Initialize the array A to store the matrix A shown below.
C     Initialize the array B to store the vector b shown below.
C
C         1 2 2        15
C     A = 2 1 2    b = 15
C         2 2 1        15
C
      DATA A / 1.0D0, 3∗2.0D0, 1.0D0, 3∗2.0D0, 1.0D0 /
      DATA B / 3∗1.5D1 /
      PRINT 1000
      PRINT 1010, ((A(IROW,ICOL), ICOL = 1, N), IROW = 1, N)
      PRINT 1020
      PRINT 1030, B
      CALL DGECO (A, LDA, N, IPIVOT, RCOND, WORK)
      PRINT 1040, RCOND
      IF ((RCOND + 1.0D0) .EQ. 1.0D0) THEN
        PRINT 1060
      END IF
      JOB = IAXEQB
      CALL DGESL (A, LDA, N, IPIVOT, B, JOB)
      PRINT 1050
      PRINT 1030, B
C
1000 FORMAT (1X, ’A:’)
1010 FORMAT (3(3X, F4.1))
1020 FORMAT (/1X, ’b:’)
1030 FORMAT (3X, F4.1)
1040 FORMAT (/1X, ’Estimated reciprocal condition number:’, F7.4)
1050 FORMAT (/1X, ’A∗∗(-1) ∗ b:’)
1060 FORMAT (1X, ’A may be singular to working precision.’)
C
      END

SAMPLE OUTPUT

A:
    1.0    2.0    2.0
    2.0    1.0    2.0
    2.0    2.0    1.0

b:
   15.0
   15.0
   15.0

Estimated reciprocal condition number: 0.1429

A∗∗(-1) ∗ b:
    3.0
    3.0
    3.0