Return to dgeequ.f CVS log

Up to [local] / rpl / lapack / lapack

Patches pour OS/2

    1:       SUBROUTINE DGEEQU( M, N, A, LDA, R, C, ROWCND, COLCND, AMAX,
    2:      $                   INFO )
    3: *
    4: *  -- LAPACK routine (version 3.2) --
    5: *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
    6: *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
    7: *     November 2006
    8: *
    9: *     .. Scalar Arguments ..
   10:       INTEGER            INFO, LDA, M, N
   11:       DOUBLE PRECISION   AMAX, COLCND, ROWCND
   12: *     ..
   13: *     .. Array Arguments ..
   14:       DOUBLE PRECISION   A( LDA, * ), C( * ), R( * )
   15: *     ..
   16: *
   17: *  Purpose
   18: *  =======
   19: *
   20: *  DGEEQU computes row and column scalings intended to equilibrate an
   21: *  M-by-N matrix A and reduce its condition number.  R returns the row
   22: *  scale factors and C the column scale factors, chosen to try to make
   23: *  the largest element in each row and column of the matrix B with
   24: *  elements B(i,j)=R(i)*A(i,j)*C(j) have absolute value 1.
   25: *
   26: *  R(i) and C(j) are restricted to be between SMLNUM = smallest safe
   27: *  number and BIGNUM = largest safe number.  Use of these scaling
   28: *  factors is not guaranteed to reduce the condition number of A but
   29: *  works well in practice.
   30: *
   31: *  Arguments
   32: *  =========
   33: *
   34: *  M       (input) INTEGER
   35: *          The number of rows of the matrix A.  M >= 0.
   36: *
   37: *  N       (input) INTEGER
   38: *          The number of columns of the matrix A.  N >= 0.
   39: *
   40: *  A       (input) DOUBLE PRECISION array, dimension (LDA,N)
   41: *          The M-by-N matrix whose equilibration factors are
   42: *          to be computed.
   43: *
   44: *  LDA     (input) INTEGER
   45: *          The leading dimension of the array A.  LDA >= max(1,M).
   46: *
   47: *  R       (output) DOUBLE PRECISION array, dimension (M)
   48: *          If INFO = 0 or INFO > M, R contains the row scale factors
   49: *          for A.
   50: *
   51: *  C       (output) DOUBLE PRECISION array, dimension (N)
   52: *          If INFO = 0,  C contains the column scale factors for A.
   53: *
   54: *  ROWCND  (output) DOUBLE PRECISION
   55: *          If INFO = 0 or INFO > M, ROWCND contains the ratio of the
   56: *          smallest R(i) to the largest R(i).  If ROWCND >= 0.1 and
   57: *          AMAX is neither too large nor too small, it is not worth
   58: *          scaling by R.
   59: *
   60: *  COLCND  (output) DOUBLE PRECISION
   61: *          If INFO = 0, COLCND contains the ratio of the smallest
   62: *          C(i) to the largest C(i).  If COLCND >= 0.1, it is not
   63: *          worth scaling by C.
   64: *
   65: *  AMAX    (output) DOUBLE PRECISION
   66: *          Absolute value of largest matrix element.  If AMAX is very
   67: *          close to overflow or very close to underflow, the matrix
   68: *          should be scaled.
   69: *
   70: *  INFO    (output) INTEGER
   71: *          = 0:  successful exit
   72: *          < 0:  if INFO = -i, the i-th argument had an illegal value
   73: *          > 0:  if INFO = i,  and i is
   74: *                <= M:  the i-th row of A is exactly zero
   75: *                >  M:  the (i-M)-th column of A is exactly zero
   76: *
   77: *  =====================================================================
   78: *
   79: *     .. Parameters ..
   80:       DOUBLE PRECISION   ONE, ZERO
   81:       PARAMETER          ( ONE = 1.0D+0, ZERO = 0.0D+0 )
   82: *     ..
   83: *     .. Local Scalars ..
   84:       INTEGER            I, J
   85:       DOUBLE PRECISION   BIGNUM, RCMAX, RCMIN, SMLNUM
   86: *     ..
   87: *     .. External Functions ..
   88:       DOUBLE PRECISION   DLAMCH
   89:       EXTERNAL           DLAMCH
   90: *     ..
   91: *     .. External Subroutines ..
   92:       EXTERNAL           XERBLA
   93: *     ..
   94: *     .. Intrinsic Functions ..
   95:       INTRINSIC          ABS, MAX, MIN
   96: *     ..
   97: *     .. Executable Statements ..
   98: *
   99: *     Test the input parameters.
  100: *
  101:       INFO = 0
  102:       IF( M.LT.0 ) THEN
  103:          INFO = -1
  104:       ELSE IF( N.LT.0 ) THEN
  105:          INFO = -2
  106:       ELSE IF( LDA.LT.MAX( 1, M ) ) THEN
  107:          INFO = -4
  108:       END IF
  109:       IF( INFO.NE.0 ) THEN
  110:          CALL XERBLA( 'DGEEQU', -INFO )
  111:          RETURN
  112:       END IF
  113: *
  114: *     Quick return if possible
  115: *
  116:       IF( M.EQ.0 .OR. N.EQ.0 ) THEN
  117:          ROWCND = ONE
  118:          COLCND = ONE
  119:          AMAX = ZERO
  120:          RETURN
  121:       END IF
  122: *
  123: *     Get machine constants.
  124: *
  125:       SMLNUM = DLAMCH( 'S' )
  126:       BIGNUM = ONE / SMLNUM
  127: *
  128: *     Compute row scale factors.
  129: *
  130:       DO 10 I = 1, M
  131:          R( I ) = ZERO
  132:    10 CONTINUE
  133: *
  134: *     Find the maximum element in each row.
  135: *
  136:       DO 30 J = 1, N
  137:          DO 20 I = 1, M
  138:             R( I ) = MAX( R( I ), ABS( A( I, J ) ) )
  139:    20    CONTINUE
  140:    30 CONTINUE
  141: *
  142: *     Find the maximum and minimum scale factors.
  143: *
  144:       RCMIN = BIGNUM
  145:       RCMAX = ZERO
  146:       DO 40 I = 1, M
  147:          RCMAX = MAX( RCMAX, R( I ) )
  148:          RCMIN = MIN( RCMIN, R( I ) )
  149:    40 CONTINUE
  150:       AMAX = RCMAX
  151: *
  152:       IF( RCMIN.EQ.ZERO ) THEN
  153: *
  154: *        Find the first zero scale factor and return an error code.
  155: *
  156:          DO 50 I = 1, M
  157:             IF( R( I ).EQ.ZERO ) THEN
  158:                INFO = I
  159:                RETURN
  160:             END IF
  161:    50    CONTINUE
  162:       ELSE
  163: *
  164: *        Invert the scale factors.
  165: *
  166:          DO 60 I = 1, M
  167:             R( I ) = ONE / MIN( MAX( R( I ), SMLNUM ), BIGNUM )
  168:    60    CONTINUE
  169: *
  170: *        Compute ROWCND = min(R(I)) / max(R(I))
  171: *
  172:          ROWCND = MAX( RCMIN, SMLNUM ) / MIN( RCMAX, BIGNUM )
  173:       END IF
  174: *
  175: *     Compute column scale factors
  176: *
  177:       DO 70 J = 1, N
  178:          C( J ) = ZERO
  179:    70 CONTINUE
  180: *
  181: *     Find the maximum element in each column,
  182: *     assuming the row scaling computed above.
  183: *
  184:       DO 90 J = 1, N
  185:          DO 80 I = 1, M
  186:             C( J ) = MAX( C( J ), ABS( A( I, J ) )*R( I ) )
  187:    80    CONTINUE
  188:    90 CONTINUE
  189: *
  190: *     Find the maximum and minimum scale factors.
  191: *
  192:       RCMIN = BIGNUM
  193:       RCMAX = ZERO
  194:       DO 100 J = 1, N
  195:          RCMIN = MIN( RCMIN, C( J ) )
  196:          RCMAX = MAX( RCMAX, C( J ) )
  197:   100 CONTINUE
  198: *
  199:       IF( RCMIN.EQ.ZERO ) THEN
  200: *
  201: *        Find the first zero scale factor and return an error code.
  202: *
  203:          DO 110 J = 1, N
  204:             IF( C( J ).EQ.ZERO ) THEN
  205:                INFO = M + J
  206:                RETURN
  207:             END IF
  208:   110    CONTINUE
  209:       ELSE
  210: *
  211: *        Invert the scale factors.
  212: *
  213:          DO 120 J = 1, N
  214:             C( J ) = ONE / MIN( MAX( C( J ), SMLNUM ), BIGNUM )
  215:   120    CONTINUE
  216: *
  217: *        Compute COLCND = min(C(J)) / max(C(J))
  218: *
  219:          COLCND = MAX( RCMIN, SMLNUM ) / MIN( RCMAX, BIGNUM )
  220:       END IF
  221: *
  222:       RETURN
  223: *
  224: *     End of DGEEQU
  225: *
  226:       END