zgeequ.f

Go to the documentation of this file.

      SUBROUTINE zgeequ( M, N, A, LDA, R, C, ROWCND, COLCND, AMAX,
     $                   info )
*
*  -- LAPACK routine (version 3.2) --
*     Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd..
*     November 2006
*
*     .. Scalar Arguments ..
      INTEGER            info, lda, m, n
      DOUBLE PRECISION   amax, colcnd, rowcnd
*     ..
*     .. Array Arguments ..
      DOUBLE PRECISION   c( * ), r( * )
      COMPLEX*16         a( lda, * )
*     ..
*
*  Purpose
*  =======
*
*  ZGEEQU computes row and column scalings intended to equilibrate an
*  M-by-N matrix A and reduce its condition number.  R returns the row
*  scale factors and C the column scale factors, chosen to try to make
*  the largest element in each row and column of the matrix B with
*  elements B(i,j)=R(i)*A(i,j)*C(j) have absolute value 1.
*
*  R(i) and C(j) are restricted to be between SMLNUM = smallest safe
*  number and BIGNUM = largest safe number.  Use of these scaling
*  factors is not guaranteed to reduce the condition number of A but
*  works well in practice.
*
*  Arguments
*  =========
*
*  M       (input) INTEGER
*          The number of rows of the matrix A.  M >= 0.
*
*  N       (input) INTEGER
*          The number of columns of the matrix A.  N >= 0.
*
*  A       (input) COMPLEX*16 array, dimension (LDA,N)
*          The M-by-N matrix whose equilibration factors are
*          to be computed.
*
*  LDA     (input) INTEGER
*          The leading dimension of the array A.  LDA >= max(1,M).
*
*  R       (output) DOUBLE PRECISION array, dimension (M)
*          If INFO = 0 or INFO > M, R contains the row scale factors
*          for A.
*
*  C       (output) DOUBLE PRECISION array, dimension (N)
*          If INFO = 0,  C contains the column scale factors for A.
*
*  ROWCND  (output) DOUBLE PRECISION
*          If INFO = 0 or INFO > M, ROWCND contains the ratio of the
*          smallest R(i) to the largest R(i).  If ROWCND >= 0.1 and
*          AMAX is neither too large nor too small, it is not worth
*          scaling by R.
*
*  COLCND  (output) DOUBLE PRECISION
*          If INFO = 0, COLCND contains the ratio of the smallest
*          C(i) to the largest C(i).  If COLCND >= 0.1, it is not
*          worth scaling by C.
*
*  AMAX    (output) DOUBLE PRECISION
*          Absolute value of largest matrix element.  If AMAX is very
*          close to overflow or very close to underflow, the matrix
*          should be scaled.
*
*  INFO    (output) INTEGER
*          = 0:  successful exit
*          < 0:  if INFO = -i, the i-th argument had an illegal value
*          > 0:  if INFO = i,  and i is
*                <= M:  the i-th row of A is exactly zero
*                >  M:  the (i-M)-th column of A is exactly zero
*
*  =====================================================================
*
*     .. Parameters ..
      DOUBLE PRECISION   one, zero
      parameter( one = 1.0d+0, zero = 0.0d+0 )
*     ..
*     .. Local Scalars ..
      INTEGER            i, j
      DOUBLE PRECISION   bignum, rcmax, rcmin, smlnum
      COMPLEX*16         zdum
*     ..
*     .. External Functions ..
      DOUBLE PRECISION   dlamch
      EXTERNAL           dlamch
*     ..
*     .. External Subroutines ..
      EXTERNAL           xerbla
*     ..
*     .. Intrinsic Functions ..
      INTRINSIC          abs, dble, dimag, max, min
*     ..
*     .. Statement Functions ..
      DOUBLE PRECISION   cabs1
*     ..
*     .. Statement Function definitions ..
      cabs1( zdum ) = abs( dble( zdum ) ) + abs( dimag( zdum ) )
*     ..
*     .. Executable Statements ..
*
*     Test the input parameters.
*
      info = 0
      IF( m.LT.0 ) THEN
         info = -1
      ELSE IF( n.LT.0 ) THEN
         info = -2
      ELSE IF( lda.LT.max( 1, m ) ) THEN
         info = -4
      END IF
      IF( info.NE.0 ) THEN
         CALL xerbla( 'ZGEEQU', -info )
         return
      END IF
*
*     Quick return if possible
*
      IF( m.EQ.0 .OR. n.EQ.0 ) THEN
         rowcnd = one
         colcnd = one
         amax = zero
         return
      END IF
*
*     Get machine constants.
*
      smlnum = dlamch( 'S' )
      bignum = one / smlnum
*
*     Compute row scale factors.
*
      DO 10 i = 1, m
         r( i ) = zero
   10 continue
*
*     Find the maximum element in each row.
*
      DO 30 j = 1, n
         DO 20 i = 1, m
            r( i ) = max( r( i ), cabs1( a( i, j ) ) )
   20    continue
   30 continue
*
*     Find the maximum and minimum scale factors.
*
      rcmin = bignum
      rcmax = zero
      DO 40 i = 1, m
         rcmax = max( rcmax, r( i ) )
         rcmin = min( rcmin, r( i ) )
   40 continue
      amax = rcmax
*
      IF( rcmin.EQ.zero ) THEN
*
*        Find the first zero scale factor and return an error code.
*
         DO 50 i = 1, m
            IF( r( i ).EQ.zero ) THEN
               info = i
               return
            END IF
   50    continue
      ELSE
*
*        Invert the scale factors.
*
         DO 60 i = 1, m
            r( i ) = one / min( max( r( i ), smlnum ), bignum )
   60    continue
*
*        Compute ROWCND = min(R(I)) / max(R(I))
*
         rowcnd = max( rcmin, smlnum ) / min( rcmax, bignum )
      END IF
*
*     Compute column scale factors
*
      DO 70 j = 1, n
         c( j ) = zero
   70 continue
*
*     Find the maximum element in each column,
*     assuming the row scaling computed above.
*
      DO 90 j = 1, n
         DO 80 i = 1, m
            c( j ) = max( c( j ), cabs1( a( i, j ) )*r( i ) )
   80    continue
   90 continue
*
*     Find the maximum and minimum scale factors.
*
      rcmin = bignum
      rcmax = zero
      DO 100 j = 1, n
         rcmin = min( rcmin, c( j ) )
         rcmax = max( rcmax, c( j ) )
  100 continue
*
      IF( rcmin.EQ.zero ) THEN
*
*        Find the first zero scale factor and return an error code.
*
         DO 110 j = 1, n
            IF( c( j ).EQ.zero ) THEN
               info = m + j
               return
            END IF
  110    continue
      ELSE
*
*        Invert the scale factors.
*
         DO 120 j = 1, n
            c( j ) = one / min( max( c( j ), smlnum ), bignum )
  120    continue
*
*        Compute COLCND = min(C(J)) / max(C(J))
*
         colcnd = max( rcmin, smlnum ) / min( rcmax, bignum )
      END IF
*
      return
*
*     End of ZGEEQU
*
      END