Return to dla_porpvgrw.f CVS log

Up to [local] / rpl / lapack / lapack

Patches pour OS/2

    1:       DOUBLE PRECISION FUNCTION DLA_PORPVGRW( UPLO, NCOLS, A, LDA, AF, 
    2:      $                                        LDAF, WORK )
    3: *
    4: *     -- LAPACK routine (version 3.2.2)                                 --
    5: *     -- Contributed by James Demmel, Deaglan Halligan, Yozo Hida and --
    6: *     -- Jason Riedy of Univ. of California Berkeley.                 --
    7: *     -- June 2010                                                    --
    8: *
    9: *     -- LAPACK is a software package provided by Univ. of Tennessee, --
   10: *     -- Univ. of California Berkeley and NAG Ltd.                    --
   11: *
   12:       IMPLICIT NONE
   13: *     ..
   14: *     .. Scalar Arguments ..
   15:       CHARACTER*1        UPLO
   16:       INTEGER            NCOLS, LDA, LDAF
   17: *     ..
   18: *     .. Array Arguments ..
   19:       DOUBLE PRECISION   A( LDA, * ), AF( LDAF, * ), WORK( * )
   20: *     ..
   21: *
   22: *  Purpose
   23: *  =======
   24: * 
   25: *  DLA_PORPVGRW computes the reciprocal pivot growth factor
   26: *  norm(A)/norm(U). The "max absolute element" norm is used. If this is
   27: *  much less than 1, the stability of the LU factorization of the
   28: *  (equilibrated) matrix A could be poor. This also means that the
   29: *  solution X, estimated condition numbers, and error bounds could be
   30: *  unreliable.
   31: *
   32: *  Arguments
   33: *  =========
   34: *
   35: *     UPLO    (input) CHARACTER*1
   36: *       = 'U':  Upper triangle of A is stored;
   37: *       = 'L':  Lower triangle of A is stored.
   38: *
   39: *     NCOLS   (input) INTEGER
   40: *     The number of columns of the matrix A. NCOLS >= 0.
   41: *
   42: *     A       (input) DOUBLE PRECISION array, dimension (LDA,N)
   43: *     On entry, the N-by-N matrix A.
   44: *
   45: *     LDA     (input) INTEGER
   46: *     The leading dimension of the array A.  LDA >= max(1,N).
   47: *
   48: *     AF      (input) DOUBLE PRECISION array, dimension (LDAF,N)
   49: *     The triangular factor U or L from the Cholesky factorization
   50: *     A = U**T*U or A = L*L**T, as computed by DPOTRF.
   51: *
   52: *     LDAF    (input) INTEGER
   53: *     The leading dimension of the array AF.  LDAF >= max(1,N).
   54: *
   55: *     WORK    (input) DOUBLE PRECISION array, dimension (2*N)
   56: *
   57: *  =====================================================================
   58: *
   59: *     .. Local Scalars ..
   60:       INTEGER            I, J
   61:       DOUBLE PRECISION   AMAX, UMAX, RPVGRW
   62:       LOGICAL            UPPER
   63: *     ..
   64: *     .. Intrinsic Functions ..
   65:       INTRINSIC          ABS, MAX, MIN
   66: *     ..
   67: *     .. External Functions ..
   68:       EXTERNAL           LSAME, DLASET
   69:       LOGICAL            LSAME
   70: *     ..
   71: *     .. Executable Statements ..
   72: *
   73:       UPPER = LSAME( 'Upper', UPLO )
   74: *
   75: *     DPOTRF will have factored only the NCOLSxNCOLS leading minor, so
   76: *     we restrict the growth search to that minor and use only the first
   77: *     2*NCOLS workspace entries.
   78: *
   79:       RPVGRW = 1.0D+0
   80:       DO I = 1, 2*NCOLS
   81:          WORK( I ) = 0.0D+0
   82:       END DO
   83: *
   84: *     Find the max magnitude entry of each column.
   85: *
   86:       IF ( UPPER ) THEN
   87:          DO J = 1, NCOLS
   88:             DO I = 1, J
   89:                WORK( NCOLS+J ) =
   90:      $              MAX( ABS( A( I, J ) ), WORK( NCOLS+J ) )
   91:             END DO
   92:          END DO
   93:       ELSE
   94:          DO J = 1, NCOLS
   95:             DO I = J, NCOLS
   96:                WORK( NCOLS+J ) =
   97:      $              MAX( ABS( A( I, J ) ), WORK( NCOLS+J ) )
   98:             END DO
   99:          END DO
  100:       END IF
  101: *
  102: *     Now find the max magnitude entry of each column of the factor in
  103: *     AF.  No pivoting, so no permutations.
  104: *
  105:       IF ( LSAME( 'Upper', UPLO ) ) THEN
  106:          DO J = 1, NCOLS
  107:             DO I = 1, J
  108:                WORK( J ) = MAX( ABS( AF( I, J ) ), WORK( J ) )
  109:             END DO
  110:          END DO
  111:       ELSE
  112:          DO J = 1, NCOLS
  113:             DO I = J, NCOLS
  114:                WORK( J ) = MAX( ABS( AF( I, J ) ), WORK( J ) )
  115:             END DO
  116:          END DO
  117:       END IF
  118: *
  119: *     Compute the *inverse* of the max element growth factor.  Dividing
  120: *     by zero would imply the largest entry of the factor's column is
  121: *     zero.  Than can happen when either the column of A is zero or
  122: *     massive pivots made the factor underflow to zero.  Neither counts
  123: *     as growth in itself, so simply ignore terms with zero
  124: *     denominators.
  125: *
  126:       IF ( LSAME( 'Upper', UPLO ) ) THEN
  127:          DO I = 1, NCOLS
  128:             UMAX = WORK( I )
  129:             AMAX = WORK( NCOLS+I )
  130:             IF ( UMAX /= 0.0D+0 ) THEN
  131:                RPVGRW = MIN( AMAX / UMAX, RPVGRW )
  132:             END IF
  133:          END DO
  134:       ELSE
  135:          DO I = 1, NCOLS
  136:             UMAX = WORK( I )
  137:             AMAX = WORK( NCOLS+I )
  138:             IF ( UMAX /= 0.0D+0 ) THEN
  139:                RPVGRW = MIN( AMAX / UMAX, RPVGRW )
  140:             END IF
  141:          END DO
  142:       END IF
  143: 
  144:       DLA_PORPVGRW = RPVGRW
  145:       END