LAPACK  3.4.1
LAPACK: Linear Algebra PACKage
zgebal.f
Go to the documentation of this file.
00001 *> \brief \b ZGEBAL
00002 *
00003 *  =========== DOCUMENTATION ===========
00004 *
00005 * Online html documentation available at 
00006 *            http://www.netlib.org/lapack/explore-html/ 
00007 *
00008 *> \htmlonly
00009 *> Download ZGEBAL + dependencies 
00010 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/zgebal.f"> 
00011 *> [TGZ]</a> 
00012 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/zgebal.f"> 
00013 *> [ZIP]</a> 
00014 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/zgebal.f"> 
00015 *> [TXT]</a>
00016 *> \endhtmlonly 
00017 *
00018 *  Definition:
00019 *  ===========
00020 *
00021 *       SUBROUTINE ZGEBAL( JOB, N, A, LDA, ILO, IHI, SCALE, INFO )
00022 * 
00023 *       .. Scalar Arguments ..
00024 *       CHARACTER          JOB
00025 *       INTEGER            IHI, ILO, INFO, LDA, N
00026 *       ..
00027 *       .. Array Arguments ..
00028 *       DOUBLE PRECISION   SCALE( * )
00029 *       COMPLEX*16         A( LDA, * )
00030 *       ..
00031 *  
00032 *
00033 *> \par Purpose:
00034 *  =============
00035 *>
00036 *> \verbatim
00037 *>
00038 *> ZGEBAL balances a general complex matrix A.  This involves, first,
00039 *> permuting A by a similarity transformation to isolate eigenvalues
00040 *> in the first 1 to ILO-1 and last IHI+1 to N elements on the
00041 *> diagonal; and second, applying a diagonal similarity transformation
00042 *> to rows and columns ILO to IHI to make the rows and columns as
00043 *> close in norm as possible.  Both steps are optional.
00044 *>
00045 *> Balancing may reduce the 1-norm of the matrix, and improve the
00046 *> accuracy of the computed eigenvalues and/or eigenvectors.
00047 *> \endverbatim
00048 *
00049 *  Arguments:
00050 *  ==========
00051 *
00052 *> \param[in] JOB
00053 *> \verbatim
00054 *>          JOB is CHARACTER*1
00055 *>          Specifies the operations to be performed on A:
00056 *>          = 'N':  none:  simply set ILO = 1, IHI = N, SCALE(I) = 1.0
00057 *>                  for i = 1,...,N;
00058 *>          = 'P':  permute only;
00059 *>          = 'S':  scale only;
00060 *>          = 'B':  both permute and scale.
00061 *> \endverbatim
00062 *>
00063 *> \param[in] N
00064 *> \verbatim
00065 *>          N is INTEGER
00066 *>          The order of the matrix A.  N >= 0.
00067 *> \endverbatim
00068 *>
00069 *> \param[in,out] A
00070 *> \verbatim
00071 *>          A is COMPLEX*16 array, dimension (LDA,N)
00072 *>          On entry, the input matrix A.
00073 *>          On exit,  A is overwritten by the balanced matrix.
00074 *>          If JOB = 'N', A is not referenced.
00075 *>          See Further Details.
00076 *> \endverbatim
00077 *>
00078 *> \param[in] LDA
00079 *> \verbatim
00080 *>          LDA is INTEGER
00081 *>          The leading dimension of the array A.  LDA >= max(1,N).
00082 *> \endverbatim
00083 *>
00084 *> \param[out] ILO
00085 *> \verbatim
00086 *> \endverbatim
00087 *>
00088 *> \param[out] IHI
00089 *> \verbatim
00090 *>          ILO and IHI are set to INTEGER such that on exit
00091 *>          A(i,j) = 0 if i > j and j = 1,...,ILO-1 or I = IHI+1,...,N.
00092 *>          If JOB = 'N' or 'S', ILO = 1 and IHI = N.
00093 *> \endverbatim
00094 *>
00095 *> \param[out] SCALE
00096 *> \verbatim
00097 *>          SCALE is DOUBLE PRECISION array, dimension (N)
00098 *>          Details of the permutations and scaling factors applied to
00099 *>          A.  If P(j) is the index of the row and column interchanged
00100 *>          with row and column j and D(j) is the scaling factor
00101 *>          applied to row and column j, then
00102 *>          SCALE(j) = P(j)    for j = 1,...,ILO-1
00103 *>                   = D(j)    for j = ILO,...,IHI
00104 *>                   = P(j)    for j = IHI+1,...,N.
00105 *>          The order in which the interchanges are made is N to IHI+1,
00106 *>          then 1 to ILO-1.
00107 *> \endverbatim
00108 *>
00109 *> \param[out] INFO
00110 *> \verbatim
00111 *>          INFO is INTEGER
00112 *>          = 0:  successful exit.
00113 *>          < 0:  if INFO = -i, the i-th argument had an illegal value.
00114 *> \endverbatim
00115 *
00116 *  Authors:
00117 *  ========
00118 *
00119 *> \author Univ. of Tennessee 
00120 *> \author Univ. of California Berkeley 
00121 *> \author Univ. of Colorado Denver 
00122 *> \author NAG Ltd. 
00123 *
00124 *> \date November 2011
00125 *
00126 *> \ingroup complex16GEcomputational
00127 *
00128 *> \par Further Details:
00129 *  =====================
00130 *>
00131 *> \verbatim
00132 *>
00133 *>  The permutations consist of row and column interchanges which put
00134 *>  the matrix in the form
00135 *>
00136 *>             ( T1   X   Y  )
00137 *>     P A P = (  0   B   Z  )
00138 *>             (  0   0   T2 )
00139 *>
00140 *>  where T1 and T2 are upper triangular matrices whose eigenvalues lie
00141 *>  along the diagonal.  The column indices ILO and IHI mark the starting
00142 *>  and ending columns of the submatrix B. Balancing consists of applying
00143 *>  a diagonal similarity transformation inv(D) * B * D to make the
00144 *>  1-norms of each row of B and its corresponding column nearly equal.
00145 *>  The output matrix is
00146 *>
00147 *>     ( T1     X*D          Y    )
00148 *>     (  0  inv(D)*B*D  inv(D)*Z ).
00149 *>     (  0      0           T2   )
00150 *>
00151 *>  Information about the permutations P and the diagonal matrix D is
00152 *>  returned in the vector SCALE.
00153 *>
00154 *>  This subroutine is based on the EISPACK routine CBAL.
00155 *>
00156 *>  Modified by Tzu-Yi Chen, Computer Science Division, University of
00157 *>    California at Berkeley, USA
00158 *> \endverbatim
00159 *>
00160 *  =====================================================================
00161       SUBROUTINE ZGEBAL( JOB, N, A, LDA, ILO, IHI, SCALE, INFO )
00162 *
00163 *  -- LAPACK computational routine (version 3.4.0) --
00164 *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
00165 *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
00166 *     November 2011
00167 *
00168 *     .. Scalar Arguments ..
00169       CHARACTER          JOB
00170       INTEGER            IHI, ILO, INFO, LDA, N
00171 *     ..
00172 *     .. Array Arguments ..
00173       DOUBLE PRECISION   SCALE( * )
00174       COMPLEX*16         A( LDA, * )
00175 *     ..
00176 *
00177 *  =====================================================================
00178 *
00179 *     .. Parameters ..
00180       DOUBLE PRECISION   ZERO, ONE
00181       PARAMETER          ( ZERO = 0.0D+0, ONE = 1.0D+0 )
00182       DOUBLE PRECISION   SCLFAC
00183       PARAMETER          ( SCLFAC = 2.0D+0 )
00184       DOUBLE PRECISION   FACTOR
00185       PARAMETER          ( FACTOR = 0.95D+0 )
00186 *     ..
00187 *     .. Local Scalars ..
00188       LOGICAL            NOCONV
00189       INTEGER            I, ICA, IEXC, IRA, J, K, L, M
00190       DOUBLE PRECISION   C, CA, F, G, R, RA, S, SFMAX1, SFMAX2, SFMIN1,
00191      $                   SFMIN2
00192       COMPLEX*16         CDUM
00193 *     ..
00194 *     .. External Functions ..
00195       LOGICAL            DISNAN, LSAME
00196       INTEGER            IZAMAX
00197       DOUBLE PRECISION   DLAMCH
00198       EXTERNAL           DISNAN, LSAME, IZAMAX, DLAMCH
00199 *     ..
00200 *     .. External Subroutines ..
00201       EXTERNAL           XERBLA, ZDSCAL, ZSWAP
00202 *     ..
00203 *     .. Intrinsic Functions ..
00204       INTRINSIC          ABS, DBLE, DIMAG, MAX, MIN
00205 *     ..
00206 *     .. Statement Functions ..
00207       DOUBLE PRECISION   CABS1
00208 *     ..
00209 *     .. Statement Function definitions ..
00210       CABS1( CDUM ) = ABS( DBLE( CDUM ) ) + ABS( DIMAG( CDUM ) )
00211 *     ..
00212 *     .. Executable Statements ..
00213 *
00214 *     Test the input parameters
00215 *
00216       INFO = 0
00217       IF( .NOT.LSAME( JOB, 'N' ) .AND. .NOT.LSAME( JOB, 'P' ) .AND.
00218      $    .NOT.LSAME( JOB, 'S' ) .AND. .NOT.LSAME( JOB, 'B' ) ) THEN
00219          INFO = -1
00220       ELSE IF( N.LT.0 ) THEN
00221          INFO = -2
00222       ELSE IF( LDA.LT.MAX( 1, N ) ) THEN
00223          INFO = -4
00224       END IF
00225       IF( INFO.NE.0 ) THEN
00226          CALL XERBLA( 'ZGEBAL', -INFO )
00227          RETURN
00228       END IF
00229 *
00230       K = 1
00231       L = N
00232 *
00233       IF( N.EQ.0 )
00234      $   GO TO 210
00235 *
00236       IF( LSAME( JOB, 'N' ) ) THEN
00237          DO 10 I = 1, N
00238             SCALE( I ) = ONE
00239    10    CONTINUE
00240          GO TO 210
00241       END IF
00242 *
00243       IF( LSAME( JOB, 'S' ) )
00244      $   GO TO 120
00245 *
00246 *     Permutation to isolate eigenvalues if possible
00247 *
00248       GO TO 50
00249 *
00250 *     Row and column exchange.
00251 *
00252    20 CONTINUE
00253       SCALE( M ) = J
00254       IF( J.EQ.M )
00255      $   GO TO 30
00256 *
00257       CALL ZSWAP( L, A( 1, J ), 1, A( 1, M ), 1 )
00258       CALL ZSWAP( N-K+1, A( J, K ), LDA, A( M, K ), LDA )
00259 *
00260    30 CONTINUE
00261       GO TO ( 40, 80 )IEXC
00262 *
00263 *     Search for rows isolating an eigenvalue and push them down.
00264 *
00265    40 CONTINUE
00266       IF( L.EQ.1 )
00267      $   GO TO 210
00268       L = L - 1
00269 *
00270    50 CONTINUE
00271       DO 70 J = L, 1, -1
00272 *
00273          DO 60 I = 1, L
00274             IF( I.EQ.J )
00275      $         GO TO 60
00276             IF( DBLE( A( J, I ) ).NE.ZERO .OR. DIMAG( A( J, I ) ).NE.
00277      $          ZERO )GO TO 70
00278    60    CONTINUE
00279 *
00280          M = L
00281          IEXC = 1
00282          GO TO 20
00283    70 CONTINUE
00284 *
00285       GO TO 90
00286 *
00287 *     Search for columns isolating an eigenvalue and push them left.
00288 *
00289    80 CONTINUE
00290       K = K + 1
00291 *
00292    90 CONTINUE
00293       DO 110 J = K, L
00294 *
00295          DO 100 I = K, L
00296             IF( I.EQ.J )
00297      $         GO TO 100
00298             IF( DBLE( A( I, J ) ).NE.ZERO .OR. DIMAG( A( I, J ) ).NE.
00299      $          ZERO )GO TO 110
00300   100    CONTINUE
00301 *
00302          M = K
00303          IEXC = 2
00304          GO TO 20
00305   110 CONTINUE
00306 *
00307   120 CONTINUE
00308       DO 130 I = K, L
00309          SCALE( I ) = ONE
00310   130 CONTINUE
00311 *
00312       IF( LSAME( JOB, 'P' ) )
00313      $   GO TO 210
00314 *
00315 *     Balance the submatrix in rows K to L.
00316 *
00317 *     Iterative loop for norm reduction
00318 *
00319       SFMIN1 = DLAMCH( 'S' ) / DLAMCH( 'P' )
00320       SFMAX1 = ONE / SFMIN1
00321       SFMIN2 = SFMIN1*SCLFAC
00322       SFMAX2 = ONE / SFMIN2
00323   140 CONTINUE
00324       NOCONV = .FALSE.
00325 *
00326       DO 200 I = K, L
00327          C = ZERO
00328          R = ZERO
00329 *
00330          DO 150 J = K, L
00331             IF( J.EQ.I )
00332      $         GO TO 150
00333             C = C + CABS1( A( J, I ) )
00334             R = R + CABS1( A( I, J ) )
00335   150    CONTINUE
00336          ICA = IZAMAX( L, A( 1, I ), 1 )
00337          CA = ABS( A( ICA, I ) )
00338          IRA = IZAMAX( N-K+1, A( I, K ), LDA )
00339          RA = ABS( A( I, IRA+K-1 ) )
00340 *
00341 *        Guard against zero C or R due to underflow.
00342 *
00343          IF( C.EQ.ZERO .OR. R.EQ.ZERO )
00344      $      GO TO 200
00345          G = R / SCLFAC
00346          F = ONE
00347          S = C + R
00348   160    CONTINUE
00349          IF( C.GE.G .OR. MAX( F, C, CA ).GE.SFMAX2 .OR.
00350      $       MIN( R, G, RA ).LE.SFMIN2 )GO TO 170
00351             IF( DISNAN( C+F+CA+R+G+RA ) ) THEN
00352 *
00353 *           Exit if NaN to avoid infinite loop
00354 *
00355             INFO = -3
00356             CALL XERBLA( 'ZGEBAL', -INFO )
00357             RETURN
00358          END IF
00359          F = F*SCLFAC
00360          C = C*SCLFAC
00361          CA = CA*SCLFAC
00362          R = R / SCLFAC
00363          G = G / SCLFAC
00364          RA = RA / SCLFAC
00365          GO TO 160
00366 *
00367   170    CONTINUE
00368          G = C / SCLFAC
00369   180    CONTINUE
00370          IF( G.LT.R .OR. MAX( R, RA ).GE.SFMAX2 .OR.
00371      $       MIN( F, C, G, CA ).LE.SFMIN2 )GO TO 190
00372          F = F / SCLFAC
00373          C = C / SCLFAC
00374          G = G / SCLFAC
00375          CA = CA / SCLFAC
00376          R = R*SCLFAC
00377          RA = RA*SCLFAC
00378          GO TO 180
00379 *
00380 *        Now balance.
00381 *
00382   190    CONTINUE
00383          IF( ( C+R ).GE.FACTOR*S )
00384      $      GO TO 200
00385          IF( F.LT.ONE .AND. SCALE( I ).LT.ONE ) THEN
00386             IF( F*SCALE( I ).LE.SFMIN1 )
00387      $         GO TO 200
00388          END IF
00389          IF( F.GT.ONE .AND. SCALE( I ).GT.ONE ) THEN
00390             IF( SCALE( I ).GE.SFMAX1 / F )
00391      $         GO TO 200
00392          END IF
00393          G = ONE / F
00394          SCALE( I ) = SCALE( I )*F
00395          NOCONV = .TRUE.
00396 *
00397          CALL ZDSCAL( N-K+1, G, A( I, K ), LDA )
00398          CALL ZDSCAL( L, F, A( 1, I ), 1 )
00399 *
00400   200 CONTINUE
00401 *
00402       IF( NOCONV )
00403      $   GO TO 140
00404 *
00405   210 CONTINUE
00406       ILO = K
00407       IHI = L
00408 *
00409       RETURN
00410 *
00411 *     End of ZGEBAL
00412 *
00413       END
 All Files Functions