LAPACK 3.12.1
LAPACK: Linear Algebra PACKage
Loading...
Searching...
No Matches

◆ sgebal()

subroutine sgebal ( character job,
integer n,
real, dimension( lda, * ) a,
integer lda,
integer ilo,
integer ihi,
real, dimension( * ) scale,
integer info )

SGEBAL

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

Purpose:
!>
!> SGEBAL balances a general real 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 REAL 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.
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 BALANC.
!>
!>  Modified by Tzu-Yi Chen, Computer Science Division, University of
!>    California at Berkeley, USA
!>
!>  Refactored by Evert Provoost, Department of Computer Science,
!>    KU Leuven, Belgium
!> 

Definition at line 160 of file sgebal.f.

161*
162* -- LAPACK computational routine --
163* -- LAPACK is a software package provided by Univ. of Tennessee, --
164* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
165*
166* .. Scalar Arguments ..
167 CHARACTER JOB
168 INTEGER IHI, ILO, INFO, LDA, N
169* ..
170* .. Array Arguments ..
171 REAL A( LDA, * ), SCALE( * )
172* ..
173*
174* =====================================================================
175*
176* .. Parameters ..
177 REAL ZERO, ONE
178 parameter( zero = 0.0e+0, one = 1.0e+0 )
179 REAL SCLFAC
180 parameter( sclfac = 2.0e+0 )
181 REAL FACTOR
182 parameter( factor = 0.95e+0 )
183* ..
184* .. Local Scalars ..
185 LOGICAL NOCONV, CANSWAP
186 INTEGER I, ICA, IRA, J, K, L
187 REAL C, CA, F, G, R, RA, S, SFMAX1, SFMAX2, SFMIN1,
188 $ SFMIN2
189* ..
190* .. External Functions ..
191 LOGICAL SISNAN, LSAME
192 INTEGER ISAMAX
193 REAL SLAMCH, SNRM2
194 EXTERNAL sisnan, lsame, isamax, slamch,
195 $ snrm2
196* ..
197* .. External Subroutines ..
198 EXTERNAL sscal, sswap, xerbla
199* ..
200* .. Intrinsic Functions ..
201 INTRINSIC abs, max, min
202* ..
203* Test the input parameters
204*
205 info = 0
206 IF( .NOT.lsame( job, 'N' ) .AND.
207 $ .NOT.lsame( job, 'P' ) .AND.
208 $ .NOT.lsame( job, 'S' ) .AND.
209 $ .NOT.lsame( job, 'B' ) ) THEN
210 info = -1
211 ELSE IF( n.LT.0 ) THEN
212 info = -2
213 ELSE IF( lda.LT.max( 1, n ) ) THEN
214 info = -4
215 END IF
216 IF( info.NE.0 ) THEN
217 CALL xerbla( 'SGEBAL', -info )
218 RETURN
219 END IF
220*
221* Quick returns.
222*
223 IF( n.EQ.0 ) THEN
224 ilo = 1
225 ihi = 0
226 RETURN
227 END IF
228*
229 IF( lsame( job, 'N' ) ) THEN
230 DO i = 1, n
231 scale( i ) = one
232 END DO
233 ilo = 1
234 ihi = n
235 RETURN
236 END IF
237*
238* Permutation to isolate eigenvalues if possible.
239*
240 k = 1
241 l = n
242*
243 IF( .NOT.lsame( job, 'S' ) ) THEN
244*
245* Row and column exchange.
246*
247 noconv = .true.
248 DO WHILE( noconv )
249*
250* Search for rows isolating an eigenvalue and push them down.
251*
252 noconv = .false.
253 DO i = l, 1, -1
254 canswap = .true.
255 DO j = 1, l
256 IF( i.NE.j .AND. a( i, j ).NE.zero ) THEN
257 canswap = .false.
258 EXIT
259 END IF
260 END DO
261*
262 IF( canswap ) THEN
263 scale( l ) = real( i )
264 IF( i.NE.l ) THEN
265 CALL sswap( l, a( 1, i ), 1, a( 1, l ), 1 )
266 CALL sswap( n-k+1, a( i, k ), lda, a( l, k ),
267 $ lda )
268 END IF
269 noconv = .true.
270*
271 IF( l.EQ.1 ) THEN
272 ilo = 1
273 ihi = 1
274 RETURN
275 END IF
276*
277 l = l - 1
278 END IF
279 END DO
280*
281 END DO
282
283 noconv = .true.
284 DO WHILE( noconv )
285*
286* Search for columns isolating an eigenvalue and push them left.
287*
288 noconv = .false.
289 DO j = k, l
290 canswap = .true.
291 DO i = k, l
292 IF( i.NE.j .AND. a( i, j ).NE.zero ) THEN
293 canswap = .false.
294 EXIT
295 END IF
296 END DO
297*
298 IF( canswap ) THEN
299 scale( k ) = real( j )
300 IF( j.NE.k ) THEN
301 CALL sswap( l, a( 1, j ), 1, a( 1, k ), 1 )
302 CALL sswap( n-k+1, a( j, k ), lda, a( k, k ),
303 $ lda )
304 END IF
305 noconv = .true.
306*
307 k = k + 1
308 END IF
309 END DO
310*
311 END DO
312*
313 END IF
314*
315* Initialize SCALE for non-permuted submatrix.
316*
317 DO i = k, l
318 scale( i ) = one
319 END DO
320*
321* If we only had to permute, we are done.
322*
323 IF( lsame( job, 'P' ) ) THEN
324 ilo = k
325 ihi = l
326 RETURN
327 END IF
328*
329* Balance the submatrix in rows K to L.
330*
331* Iterative loop for norm reduction.
332*
333 sfmin1 = slamch( 'S' ) / slamch( 'P' )
334 sfmax1 = one / sfmin1
335 sfmin2 = sfmin1*sclfac
336 sfmax2 = one / sfmin2
337*
338 noconv = .true.
339 DO WHILE( noconv )
340 noconv = .false.
341*
342 DO i = k, l
343*
344 c = snrm2( l-k+1, a( k, i ), 1 )
345 r = snrm2( l-k+1, a( i, k ), lda )
346 ica = isamax( l, a( 1, i ), 1 )
347 ca = abs( a( ica, i ) )
348 ira = isamax( n-k+1, a( i, k ), lda )
349 ra = abs( a( i, ira+k-1 ) )
350*
351* Guard against zero C or R due to underflow.
352*
353 IF( c.EQ.zero .OR. r.EQ.zero ) cycle
354*
355* Exit if NaN to avoid infinite loop
356*
357 IF( sisnan( c+ca+r+ra ) ) THEN
358 info = -3
359 CALL xerbla( 'SGEBAL', -info )
360 RETURN
361 END IF
362*
363 g = r / sclfac
364 f = one
365 s = c + r
366*
367 DO WHILE( c.LT.g .AND. max( f, c, ca ).LT.sfmax2 .AND.
368 $ min( r, g, ra ).GT.sfmin2 )
369 f = f*sclfac
370 c = c*sclfac
371 ca = ca*sclfac
372 r = r / sclfac
373 g = g / sclfac
374 ra = ra / sclfac
375 END DO
376*
377 g = c / sclfac
378*
379 DO WHILE( g.GE.r .AND. max( r, ra ).LT.sfmax2 .AND.
380 $ min( f, c, g, ca ).GT.sfmin2 )
381 f = f / sclfac
382 c = c / sclfac
383 g = g / sclfac
384 ca = ca / sclfac
385 r = r*sclfac
386 ra = ra*sclfac
387 END DO
388*
389* Now balance.
390*
391 IF( ( c+r ).GE.factor*s ) cycle
392 IF( f.LT.one .AND. scale( i ).LT.one ) THEN
393 IF( f*scale( i ).LE.sfmin1 ) cycle
394 END IF
395 IF( f.GT.one .AND. scale( i ).GT.one ) THEN
396 IF( scale( i ).GE.sfmax1 / f ) cycle
397 END IF
398 g = one / f
399 scale( i ) = scale( i )*f
400 noconv = .true.
401*
402 CALL sscal( n-k+1, g, a( i, k ), lda )
403 CALL sscal( l, f, a( 1, i ), 1 )
404*
405 END DO
406*
407 END DO
408*
409 ilo = k
410 ihi = l
411*
412 RETURN
413*
414* End of SGEBAL
415*
subroutine xerbla(srname, info)
Definition cblat2.f:3285
integer function isamax(n, sx, incx)
ISAMAX
Definition isamax.f:71
logical function sisnan(sin)
SISNAN tests input for NaN.
Definition sisnan.f:57
real function slamch(cmach)
SLAMCH
Definition slamch.f:68
logical function lsame(ca, cb)
LSAME
Definition lsame.f:48
real(wp) function snrm2(n, x, incx)
SNRM2
Definition snrm2.f90:89
subroutine sscal(n, sa, sx, incx)
SSCAL
Definition sscal.f:79
subroutine sswap(n, sx, incx, sy, incy)
SSWAP
Definition sswap.f:82
Here is the call graph for this function:
Here is the caller graph for this function: