LAPACK  3.4.1
LAPACK: Linear Algebra PACKage
zgbequb.f
Go to the documentation of this file.
00001 *> \brief \b ZGBEQUB
00002 *
00003 *  =========== DOCUMENTATION ===========
00004 *
00005 * Online html documentation available at 
00006 *            http://www.netlib.org/lapack/explore-html/ 
00007 *
00008 *> \htmlonly
00009 *> Download ZGBEQUB + dependencies 
00010 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/zgbequb.f"> 
00011 *> [TGZ]</a> 
00012 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/zgbequb.f"> 
00013 *> [ZIP]</a> 
00014 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/zgbequb.f"> 
00015 *> [TXT]</a>
00016 *> \endhtmlonly 
00017 *
00018 *  Definition:
00019 *  ===========
00020 *
00021 *       SUBROUTINE ZGBEQUB( M, N, KL, KU, AB, LDAB, R, C, ROWCND, COLCND,
00022 *                           AMAX, INFO )
00023 * 
00024 *       .. Scalar Arguments ..
00025 *       INTEGER            INFO, KL, KU, LDAB, M, N
00026 *       DOUBLE PRECISION   AMAX, COLCND, ROWCND
00027 *       ..
00028 *       .. Array Arguments ..
00029 *       DOUBLE PRECISION   C( * ), R( * )
00030 *       COMPLEX*16         AB( LDAB, * )
00031 *       ..
00032 *  
00033 *
00034 *> \par Purpose:
00035 *  =============
00036 *>
00037 *> \verbatim
00038 *>
00039 *> ZGBEQUB computes row and column scalings intended to equilibrate an
00040 *> M-by-N matrix A and reduce its condition number.  R returns the row
00041 *> scale factors and C the column scale factors, chosen to try to make
00042 *> the largest element in each row and column of the matrix B with
00043 *> elements B(i,j)=R(i)*A(i,j)*C(j) have an absolute value of at most
00044 *> the radix.
00045 *>
00046 *> R(i) and C(j) are restricted to be a power of the radix between
00047 *> SMLNUM = smallest safe number and BIGNUM = largest safe number.  Use
00048 *> of these scaling factors is not guaranteed to reduce the condition
00049 *> number of A but works well in practice.
00050 *>
00051 *> This routine differs from ZGEEQU by restricting the scaling factors
00052 *> to a power of the radix.  Baring over- and underflow, scaling by
00053 *> these factors introduces no additional rounding errors.  However, the
00054 *> scaled entries' magnitured are no longer approximately 1 but lie
00055 *> between sqrt(radix) and 1/sqrt(radix).
00056 *> \endverbatim
00057 *
00058 *  Arguments:
00059 *  ==========
00060 *
00061 *> \param[in] M
00062 *> \verbatim
00063 *>          M is INTEGER
00064 *>          The number of rows of the matrix A.  M >= 0.
00065 *> \endverbatim
00066 *>
00067 *> \param[in] N
00068 *> \verbatim
00069 *>          N is INTEGER
00070 *>          The number of columns of the matrix A.  N >= 0.
00071 *> \endverbatim
00072 *>
00073 *> \param[in] KL
00074 *> \verbatim
00075 *>          KL is INTEGER
00076 *>          The number of subdiagonals within the band of A.  KL >= 0.
00077 *> \endverbatim
00078 *>
00079 *> \param[in] KU
00080 *> \verbatim
00081 *>          KU is INTEGER
00082 *>          The number of superdiagonals within the band of A.  KU >= 0.
00083 *> \endverbatim
00084 *>
00085 *> \param[in] AB
00086 *> \verbatim
00087 *>          AB is DOUBLE PRECISION array, dimension (LDAB,N)
00088 *>          On entry, the matrix A in band storage, in rows 1 to KL+KU+1.
00089 *>          The j-th column of A is stored in the j-th column of the
00090 *>          array AB as follows:
00091 *>          AB(KU+1+i-j,j) = A(i,j) for max(1,j-KU)<=i<=min(N,j+kl)
00092 *> \endverbatim
00093 *>
00094 *> \param[in] LDAB
00095 *> \verbatim
00096 *>          LDAB is INTEGER
00097 *>          The leading dimension of the array A.  LDAB >= max(1,M).
00098 *> \endverbatim
00099 *>
00100 *> \param[out] R
00101 *> \verbatim
00102 *>          R is DOUBLE PRECISION array, dimension (M)
00103 *>          If INFO = 0 or INFO > M, R contains the row scale factors
00104 *>          for A.
00105 *> \endverbatim
00106 *>
00107 *> \param[out] C
00108 *> \verbatim
00109 *>          C is DOUBLE PRECISION array, dimension (N)
00110 *>          If INFO = 0,  C contains the column scale factors for A.
00111 *> \endverbatim
00112 *>
00113 *> \param[out] ROWCND
00114 *> \verbatim
00115 *>          ROWCND is DOUBLE PRECISION
00116 *>          If INFO = 0 or INFO > M, ROWCND contains the ratio of the
00117 *>          smallest R(i) to the largest R(i).  If ROWCND >= 0.1 and
00118 *>          AMAX is neither too large nor too small, it is not worth
00119 *>          scaling by R.
00120 *> \endverbatim
00121 *>
00122 *> \param[out] COLCND
00123 *> \verbatim
00124 *>          COLCND is DOUBLE PRECISION
00125 *>          If INFO = 0, COLCND contains the ratio of the smallest
00126 *>          C(i) to the largest C(i).  If COLCND >= 0.1, it is not
00127 *>          worth scaling by C.
00128 *> \endverbatim
00129 *>
00130 *> \param[out] AMAX
00131 *> \verbatim
00132 *>          AMAX is DOUBLE PRECISION
00133 *>          Absolute value of largest matrix element.  If AMAX is very
00134 *>          close to overflow or very close to underflow, the matrix
00135 *>          should be scaled.
00136 *> \endverbatim
00137 *>
00138 *> \param[out] INFO
00139 *> \verbatim
00140 *>          INFO is INTEGER
00141 *>          = 0:  successful exit
00142 *>          < 0:  if INFO = -i, the i-th argument had an illegal value
00143 *>          > 0:  if INFO = i,  and i is
00144 *>                <= M:  the i-th row of A is exactly zero
00145 *>                >  M:  the (i-M)-th column of A is exactly zero
00146 *> \endverbatim
00147 *
00148 *  Authors:
00149 *  ========
00150 *
00151 *> \author Univ. of Tennessee 
00152 *> \author Univ. of California Berkeley 
00153 *> \author Univ. of Colorado Denver 
00154 *> \author NAG Ltd. 
00155 *
00156 *> \date November 2011
00157 *
00158 *> \ingroup complex16GBcomputational
00159 *
00160 *  =====================================================================
00161       SUBROUTINE ZGBEQUB( M, N, KL, KU, AB, LDAB, R, C, ROWCND, COLCND,
00162      $                    AMAX, INFO )
00163 *
00164 *  -- LAPACK computational routine (version 3.4.0) --
00165 *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
00166 *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
00167 *     November 2011
00168 *
00169 *     .. Scalar Arguments ..
00170       INTEGER            INFO, KL, KU, LDAB, M, N
00171       DOUBLE PRECISION   AMAX, COLCND, ROWCND
00172 *     ..
00173 *     .. Array Arguments ..
00174       DOUBLE PRECISION   C( * ), R( * )
00175       COMPLEX*16         AB( LDAB, * )
00176 *     ..
00177 *
00178 *  =====================================================================
00179 *
00180 *     .. Parameters ..
00181       DOUBLE PRECISION   ONE, ZERO
00182       PARAMETER          ( ONE = 1.0D+0, ZERO = 0.0D+0 )
00183 *     ..
00184 *     .. Local Scalars ..
00185       INTEGER            I, J, KD
00186       DOUBLE PRECISION   BIGNUM, RCMAX, RCMIN, SMLNUM, RADIX,
00187      $                   LOGRDX
00188       COMPLEX*16         ZDUM
00189 *     ..
00190 *     .. External Functions ..
00191       DOUBLE PRECISION   DLAMCH
00192       EXTERNAL           DLAMCH
00193 *     ..
00194 *     .. External Subroutines ..
00195       EXTERNAL           XERBLA
00196 *     ..
00197 *     .. Intrinsic Functions ..
00198       INTRINSIC          ABS, MAX, MIN, LOG, REAL, DIMAG
00199 *     ..
00200 *     .. Statement Functions ..
00201       DOUBLE PRECISION   CABS1
00202 *     ..
00203 *     .. Statement Function definitions ..
00204       CABS1( ZDUM ) = ABS( DBLE( ZDUM ) ) + ABS( DIMAG( ZDUM ) )
00205 *     ..
00206 *     .. Executable Statements ..
00207 *
00208 *     Test the input parameters.
00209 *
00210       INFO = 0
00211       IF( M.LT.0 ) THEN
00212          INFO = -1
00213       ELSE IF( N.LT.0 ) THEN
00214          INFO = -2
00215       ELSE IF( KL.LT.0 ) THEN
00216          INFO = -3
00217       ELSE IF( KU.LT.0 ) THEN
00218          INFO = -4
00219       ELSE IF( LDAB.LT.KL+KU+1 ) THEN
00220          INFO = -6
00221       END IF
00222       IF( INFO.NE.0 ) THEN
00223          CALL XERBLA( 'ZGBEQUB', -INFO )
00224          RETURN
00225       END IF
00226 *
00227 *     Quick return if possible.
00228 *
00229       IF( M.EQ.0 .OR. N.EQ.0 ) THEN
00230          ROWCND = ONE
00231          COLCND = ONE
00232          AMAX = ZERO
00233          RETURN
00234       END IF
00235 *
00236 *     Get machine constants.  Assume SMLNUM is a power of the radix.
00237 *
00238       SMLNUM = DLAMCH( 'S' )
00239       BIGNUM = ONE / SMLNUM
00240       RADIX = DLAMCH( 'B' )
00241       LOGRDX = LOG(RADIX)
00242 *
00243 *     Compute row scale factors.
00244 *
00245       DO 10 I = 1, M
00246          R( I ) = ZERO
00247    10 CONTINUE
00248 *
00249 *     Find the maximum element in each row.
00250 *
00251       KD = KU + 1
00252       DO 30 J = 1, N
00253          DO 20 I = MAX( J-KU, 1 ), MIN( J+KL, M )
00254             R( I ) = MAX( R( I ), CABS1( AB( KD+I-J, J ) ) )
00255    20    CONTINUE
00256    30 CONTINUE
00257       DO I = 1, M
00258          IF( R( I ).GT.ZERO ) THEN
00259             R( I ) = RADIX**INT( LOG( R( I ) ) / LOGRDX )
00260          END IF
00261       END DO
00262 *
00263 *     Find the maximum and minimum scale factors.
00264 *
00265       RCMIN = BIGNUM
00266       RCMAX = ZERO
00267       DO 40 I = 1, M
00268          RCMAX = MAX( RCMAX, R( I ) )
00269          RCMIN = MIN( RCMIN, R( I ) )
00270    40 CONTINUE
00271       AMAX = RCMAX
00272 *
00273       IF( RCMIN.EQ.ZERO ) THEN
00274 *
00275 *        Find the first zero scale factor and return an error code.
00276 *
00277          DO 50 I = 1, M
00278             IF( R( I ).EQ.ZERO ) THEN
00279                INFO = I
00280                RETURN
00281             END IF
00282    50    CONTINUE
00283       ELSE
00284 *
00285 *        Invert the scale factors.
00286 *
00287          DO 60 I = 1, M
00288             R( I ) = ONE / MIN( MAX( R( I ), SMLNUM ), BIGNUM )
00289    60    CONTINUE
00290 *
00291 *        Compute ROWCND = min(R(I)) / max(R(I)).
00292 *
00293          ROWCND = MAX( RCMIN, SMLNUM ) / MIN( RCMAX, BIGNUM )
00294       END IF
00295 *
00296 *     Compute column scale factors.
00297 *
00298       DO 70 J = 1, N
00299          C( J ) = ZERO
00300    70 CONTINUE
00301 *
00302 *     Find the maximum element in each column,
00303 *     assuming the row scaling computed above.
00304 *
00305       DO 90 J = 1, N
00306          DO 80 I = MAX( J-KU, 1 ), MIN( J+KL, M )
00307             C( J ) = MAX( C( J ), CABS1( AB( KD+I-J, J ) )*R( I ) )
00308    80    CONTINUE
00309          IF( C( J ).GT.ZERO ) THEN
00310             C( J ) = RADIX**INT( LOG( C( J ) ) / LOGRDX )
00311          END IF
00312    90 CONTINUE
00313 *
00314 *     Find the maximum and minimum scale factors.
00315 *
00316       RCMIN = BIGNUM
00317       RCMAX = ZERO
00318       DO 100 J = 1, N
00319          RCMIN = MIN( RCMIN, C( J ) )
00320          RCMAX = MAX( RCMAX, C( J ) )
00321   100 CONTINUE
00322 *
00323       IF( RCMIN.EQ.ZERO ) THEN
00324 *
00325 *        Find the first zero scale factor and return an error code.
00326 *
00327          DO 110 J = 1, N
00328             IF( C( J ).EQ.ZERO ) THEN
00329                INFO = M + J
00330                RETURN
00331             END IF
00332   110    CONTINUE
00333       ELSE
00334 *
00335 *        Invert the scale factors.
00336 *
00337          DO 120 J = 1, N
00338             C( J ) = ONE / MIN( MAX( C( J ), SMLNUM ), BIGNUM )
00339   120    CONTINUE
00340 *
00341 *        Compute COLCND = min(C(J)) / max(C(J)).
00342 *
00343          COLCND = MAX( RCMIN, SMLNUM ) / MIN( RCMAX, BIGNUM )
00344       END IF
00345 *
00346       RETURN
00347 *
00348 *     End of ZGBEQUB
00349 *
00350       END
 All Files Functions