LAPACK  3.5.0
LAPACK: Linear Algebra PACKage
cgebal.f File Reference

Go to the source code of this file.

Functions/Subroutines

subroutine cgebal (JOB, N, A, LDA, ILO, IHI, SCALE, INFO)
 CGEBAL More...
 

Function/Subroutine Documentation

subroutine cgebal ( character  JOB,
integer  N,
complex, dimension( lda, * )  A,
integer  LDA,
integer  ILO,
integer  IHI,
real, dimension( * )  SCALE,
integer  INFO 
)

CGEBAL

Download CGEBAL + dependencies [TGZ] [ZIP] [TXT]

Purpose:
 CGEBAL balances a general complex matrix A.  This involves, first,
 permuting A by a similarity transformation to isolate eigenvalues
 in the first 1 to ILO-1 and last IHI+1 to N elements on the
 diagonal; and second, applying a diagonal similarity transformation
 to rows and columns ILO to IHI to make the rows and columns as
 close in norm as possible.  Both steps are optional.

 Balancing may reduce the 1-norm of the matrix, and improve the
 accuracy of the computed eigenvalues and/or eigenvectors.
Parameters
[in]JOB
          JOB is CHARACTER*1
          Specifies the operations to be performed on A:
          = 'N':  none:  simply set ILO = 1, IHI = N, SCALE(I) = 1.0
                  for i = 1,...,N;
          = 'P':  permute only;
          = 'S':  scale only;
          = 'B':  both permute and scale.
[in]N
          N is INTEGER
          The order of the matrix A.  N >= 0.
[in,out]A
          A is COMPLEX array, dimension (LDA,N)
          On entry, the input matrix A.
          On exit,  A is overwritten by the balanced matrix.
          If JOB = 'N', A is not referenced.
          See Further Details.
[in]LDA
          LDA is INTEGER
          The leading dimension of the array A.  LDA >= max(1,N).
[out]ILO
          ILO is INTEGER
[out]IHI
          IHI is INTEGER
          ILO and IHI are set to integers such that on exit
          A(i,j) = 0 if i > j and j = 1,...,ILO-1 or I = IHI+1,...,N.
          If JOB = 'N' or 'S', ILO = 1 and IHI = N.
[out]SCALE
          SCALE is REAL array, dimension (N)
          Details of the permutations and scaling factors applied to
          A.  If P(j) is the index of the row and column interchanged
          with row and column j and D(j) is the scaling factor
          applied to row and column j, then
          SCALE(j) = P(j)    for j = 1,...,ILO-1
                   = D(j)    for j = ILO,...,IHI
                   = P(j)    for j = IHI+1,...,N.
          The order in which the interchanges are made is N to IHI+1,
          then 1 to ILO-1.
[out]INFO
          INFO is INTEGER
          = 0:  successful exit.
          < 0:  if INFO = -i, the i-th argument had an illegal value.
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Date
November 2013
Further Details:
  The permutations consist of row and column interchanges which put
  the matrix in the form

             ( T1   X   Y  )
     P A P = (  0   B   Z  )
             (  0   0   T2 )

  where T1 and T2 are upper triangular matrices whose eigenvalues lie
  along the diagonal.  The column indices ILO and IHI mark the starting
  and ending columns of the submatrix B. Balancing consists of applying
  a diagonal similarity transformation inv(D) * B * D to make the
  1-norms of each row of B and its corresponding column nearly equal.
  The output matrix is

     ( T1     X*D          Y    )
     (  0  inv(D)*B*D  inv(D)*Z ).
     (  0      0           T2   )

  Information about the permutations P and the diagonal matrix D is
  returned in the vector SCALE.

  This subroutine is based on the EISPACK routine CBAL.

  Modified by Tzu-Yi Chen, Computer Science Division, University of
    California at Berkeley, USA

Definition at line 163 of file cgebal.f.

163 *
164 * -- LAPACK computational routine (version 3.5.0) --
165 * -- LAPACK is a software package provided by Univ. of Tennessee, --
166 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
167 * November 2013
168 *
169 * .. Scalar Arguments ..
170  CHARACTER job
171  INTEGER ihi, ilo, info, lda, n
172 * ..
173 * .. Array Arguments ..
174  REAL scale( * )
175  COMPLEX a( lda, * )
176 * ..
177 *
178 * =====================================================================
179 *
180 * .. Parameters ..
181  REAL zero, one
182  parameter( zero = 0.0e+0, one = 1.0e+0 )
183  REAL sclfac
184  parameter( sclfac = 2.0e+0 )
185  REAL factor
186  parameter( factor = 0.95e+0 )
187 * ..
188 * .. Local Scalars ..
189  LOGICAL noconv
190  INTEGER i, ica, iexc, ira, j, k, l, m
191  REAL c, ca, f, g, r, ra, s, sfmax1, sfmax2, sfmin1,
192  $ sfmin2
193  COMPLEX cdum
194 * ..
195 * .. External Functions ..
196  LOGICAL sisnan, lsame
197  INTEGER icamax
198  REAL slamch, scnrm2
199  EXTERNAL sisnan, lsame, icamax, slamch, scnrm2
200 * ..
201 * .. External Subroutines ..
202  EXTERNAL csscal, cswap, xerbla
203 * ..
204 * .. Intrinsic Functions ..
205  INTRINSIC abs, aimag, max, min, real
206 * ..
207 * .. Statement Functions ..
208  REAL cabs1
209 * ..
210 * .. Statement Function definitions ..
211  cabs1( cdum ) = abs( REAL( CDUM ) ) + abs( aimag( cdum ) )
212 * ..
213 * .. Executable Statements ..
214 *
215 * Test the input parameters
216 *
217  info = 0
218  IF( .NOT.lsame( job, 'N' ) .AND. .NOT.lsame( job, 'P' ) .AND.
219  $ .NOT.lsame( job, 'S' ) .AND. .NOT.lsame( job, 'B' ) ) THEN
220  info = -1
221  ELSE IF( n.LT.0 ) THEN
222  info = -2
223  ELSE IF( lda.LT.max( 1, n ) ) THEN
224  info = -4
225  END IF
226  IF( info.NE.0 ) THEN
227  CALL xerbla( 'CGEBAL', -info )
228  RETURN
229  END IF
230 *
231  k = 1
232  l = n
233 *
234  IF( n.EQ.0 )
235  $ GO TO 210
236 *
237  IF( lsame( job, 'N' ) ) THEN
238  DO 10 i = 1, n
239  scale( i ) = one
240  10 CONTINUE
241  GO TO 210
242  END IF
243 *
244  IF( lsame( job, 'S' ) )
245  $ GO TO 120
246 *
247 * Permutation to isolate eigenvalues if possible
248 *
249  GO TO 50
250 *
251 * Row and column exchange.
252 *
253  20 CONTINUE
254  scale( m ) = j
255  IF( j.EQ.m )
256  $ GO TO 30
257 *
258  CALL cswap( l, a( 1, j ), 1, a( 1, m ), 1 )
259  CALL cswap( n-k+1, a( j, k ), lda, a( m, k ), lda )
260 *
261  30 CONTINUE
262  GO TO ( 40, 80 )iexc
263 *
264 * Search for rows isolating an eigenvalue and push them down.
265 *
266  40 CONTINUE
267  IF( l.EQ.1 )
268  $ GO TO 210
269  l = l - 1
270 *
271  50 CONTINUE
272  DO 70 j = l, 1, -1
273 *
274  DO 60 i = 1, l
275  IF( i.EQ.j )
276  $ GO TO 60
277  IF( REAL( A( J, I ) ).NE.zero .OR. aimag( a( j, i ) ).NE.
278  $ zero )GO TO 70
279  60 CONTINUE
280 *
281  m = l
282  iexc = 1
283  GO TO 20
284  70 CONTINUE
285 *
286  GO TO 90
287 *
288 * Search for columns isolating an eigenvalue and push them left.
289 *
290  80 CONTINUE
291  k = k + 1
292 *
293  90 CONTINUE
294  DO 110 j = k, l
295 *
296  DO 100 i = k, l
297  IF( i.EQ.j )
298  $ GO TO 100
299  IF( REAL( A( I, J ) ).NE.zero .OR. aimag( a( i, j ) ).NE.
300  $ zero )GO TO 110
301  100 CONTINUE
302 *
303  m = k
304  iexc = 2
305  GO TO 20
306  110 CONTINUE
307 *
308  120 CONTINUE
309  DO 130 i = k, l
310  scale( i ) = one
311  130 CONTINUE
312 *
313  IF( lsame( job, 'P' ) )
314  $ GO TO 210
315 *
316 * Balance the submatrix in rows K to L.
317 *
318 * Iterative loop for norm reduction
319 *
320  sfmin1 = slamch( 'S' ) / slamch( 'P' )
321  sfmax1 = one / sfmin1
322  sfmin2 = sfmin1*sclfac
323  sfmax2 = one / sfmin2
324  140 CONTINUE
325  noconv = .false.
326 *
327  DO 200 i = k, l
328 *
329  c = scnrm2( l-k+1, a( k, i ), 1 )
330  r = scnrm2( l-k+1, a( i , k ), lda )
331  ica = icamax( l, a( 1, i ), 1 )
332  ca = abs( a( ica, i ) )
333  ira = icamax( n-k+1, a( i, k ), lda )
334  ra = abs( a( i, ira+k-1 ) )
335 *
336 * Guard against zero C or R due to underflow.
337 *
338  IF( c.EQ.zero .OR. r.EQ.zero )
339  $ GO TO 200
340  g = r / sclfac
341  f = one
342  s = c + r
343  160 CONTINUE
344  IF( c.GE.g .OR. max( f, c, ca ).GE.sfmax2 .OR.
345  $ min( r, g, ra ).LE.sfmin2 )GO TO 170
346  IF( sisnan( c+f+ca+r+g+ra ) ) THEN
347 *
348 * Exit if NaN to avoid infinite loop
349 *
350  info = -3
351  CALL xerbla( 'CGEBAL', -info )
352  RETURN
353  END IF
354  f = f*sclfac
355  c = c*sclfac
356  ca = ca*sclfac
357  r = r / sclfac
358  g = g / sclfac
359  ra = ra / sclfac
360  GO TO 160
361 *
362  170 CONTINUE
363  g = c / sclfac
364  180 CONTINUE
365  IF( g.LT.r .OR. max( r, ra ).GE.sfmax2 .OR.
366  $ min( f, c, g, ca ).LE.sfmin2 )GO TO 190
367  f = f / sclfac
368  c = c / sclfac
369  g = g / sclfac
370  ca = ca / sclfac
371  r = r*sclfac
372  ra = ra*sclfac
373  GO TO 180
374 *
375 * Now balance.
376 *
377  190 CONTINUE
378  IF( ( c+r ).GE.factor*s )
379  $ GO TO 200
380  IF( f.LT.one .AND. scale( i ).LT.one ) THEN
381  IF( f*scale( i ).LE.sfmin1 )
382  $ GO TO 200
383  END IF
384  IF( f.GT.one .AND. scale( i ).GT.one ) THEN
385  IF( scale( i ).GE.sfmax1 / f )
386  $ GO TO 200
387  END IF
388  g = one / f
389  scale( i ) = scale( i )*f
390  noconv = .true.
391 *
392  CALL csscal( n-k+1, g, a( i, k ), lda )
393  CALL csscal( l, f, a( 1, i ), 1 )
394 *
395  200 CONTINUE
396 *
397  IF( noconv )
398  $ GO TO 140
399 *
400  210 CONTINUE
401  ilo = k
402  ihi = l
403 *
404  RETURN
405 *
406 * End of CGEBAL
407 *
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:62
integer function icamax(N, CX, INCX)
ICAMAX
Definition: icamax.f:53
subroutine cswap(N, CX, INCX, CY, INCY)
CSWAP
Definition: cswap.f:52
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:55
real function scnrm2(N, X, INCX)
SCNRM2
Definition: scnrm2.f:56
real function slamch(CMACH)
SLAMCH
Definition: slamch.f:69
subroutine csscal(N, SA, CX, INCX)
CSSCAL
Definition: csscal.f:54
logical function sisnan(SIN)
SISNAN tests input for NaN.
Definition: sisnan.f:61

Here is the call graph for this function:

Here is the caller graph for this function: