File:  [local] / rpl / lapack / lapack / dgbequ.f
Revision 1.7: download - view: text, annotated - select for diffs - revision graph
Tue Dec 21 13:53:24 2010 UTC (13 years, 4 months ago) by bertrand
Branches: MAIN
CVS tags: rpl-4_1_3, rpl-4_1_2, rpl-4_1_1, rpl-4_1_0, rpl-4_0_24, rpl-4_0_22, rpl-4_0_21, rpl-4_0_20, rpl-4_0, HEAD
Mise à jour de lapack vers la version 3.3.0.

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

CVSweb interface <joel.bertrand@systella.fr>