Return to zla_porpvgrw.f CVS log

Up to [local] / rpl / lapack / lapack

Mise à jour de lapack vers la version 3.3.0.

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