89 SUBROUTINE sget37( RMAX, LMAX, NINFO, KNT, NIN )
99 INTEGER LMAX( 3 ), NINFO( 3 )
107 parameter( zero = 0.0e0, one = 1.0e0, two = 2.0e0 )
109 parameter( epsin = 5.9605e-8 )
111 parameter( ldt = 20, lwork = 2*ldt*( 10+ldt ) )
114 INTEGER I, ICMP, IFND, INFO, ISCL, J, KMIN, M, N
115 REAL BIGNUM, EPS, SMLNUM, TNRM, TOL, TOLIN, V,
116 $ VIMIN, VMAX, VMUL, VRMIN
119 LOGICAL SELECT( LDT )
120 INTEGER IWORK( 2*LDT ), LCMP( 3 )
121 REAL DUM( 1 ), LE( LDT, LDT ), RE( LDT, LDT ),
122 $ S( LDT ), SEP( LDT ), SEPIN( LDT ),
123 $ SEPTMP( LDT ), SIN( LDT ), STMP( LDT ),
124 $ T( LDT, LDT ), TMP( LDT, LDT ), VAL( 3 ),
125 $ WI( LDT ), WIIN( LDT ), WITMP( LDT ),
126 $ WORK( LWORK ), WR( LDT ), WRIN( LDT ),
131 EXTERNAL slamch, slange
138 INTRINSIC max, real, sqrt
143 smlnum = slamch(
'S' ) / eps
144 bignum = one / smlnum
148 eps = max( eps, epsin )
160 val( 1 ) = sqrt( smlnum )
162 val( 3 ) = sqrt( bignum )
169 READ( nin, fmt = * )n
173 READ( nin, fmt = * )( tmp( i, j ), j = 1, n )
176 READ( nin, fmt = * )wrin( i ), wiin( i ), sin( i ), sepin( i )
178 tnrm = slange(
'M', n, n, tmp, ldt, work )
187 CALL slacpy(
'F', n, n, tmp, ldt, t, ldt )
190 CALL sscal( n, vmul, t( 1, i ), 1 )
197 CALL sgehrd( n, 1, n, t, ldt, work( 1 ), work( n+1 ), lwork-n,
201 ninfo( 1 ) = ninfo( 1 ) + 1
212 CALL shseqr(
'S',
'N', n, 1, n, t, ldt, wr, wi, dum, 1, work,
216 ninfo( 2 ) = ninfo( 2 ) + 1
222 CALL strevc(
'Both',
'All',
SELECT, n, t, ldt, le, ldt, re,
223 $ ldt, n, m, work, info )
227 CALL strsna(
'Both',
'All',
SELECT, n, t, ldt, le, ldt, re,
228 $ ldt, s, sep, n, m, work, n, iwork, info )
231 ninfo( 3 ) = ninfo( 3 ) + 1
238 CALL scopy( n, wr, 1, wrtmp, 1 )
239 CALL scopy( n, wi, 1, witmp, 1 )
240 CALL scopy( n, s, 1, stmp, 1 )
241 CALL scopy( n, sep, 1, septmp, 1 )
242 CALL sscal( n, one / vmul, septmp, 1 )
248 IF( wrtmp( j ).LT.vrmin )
THEN
254 wrtmp( kmin ) = wrtmp( i )
255 witmp( kmin ) = witmp( i )
259 stmp( kmin ) = stmp( i )
261 vrmin = septmp( kmin )
262 septmp( kmin ) = septmp( i )
269 v = max( two*real( n )*eps*tnrm, smlnum )
273 IF( v.GT.septmp( i ) )
THEN
276 tol = v / septmp( i )
278 IF( v.GT.sepin( i ) )
THEN
281 tolin = v / sepin( i )
283 tol = max( tol, smlnum / eps )
284 tolin = max( tolin, smlnum / eps )
285 IF( eps*( sin( i )-tolin ).GT.stmp( i )+tol )
THEN
287 ELSE IF( sin( i )-tolin.GT.stmp( i )+tol )
THEN
288 vmax = ( sin( i )-tolin ) / ( stmp( i )+tol )
289 ELSE IF( sin( i )+tolin.LT.eps*( stmp( i )-tol ) )
THEN
291 ELSE IF( sin( i )+tolin.LT.stmp( i )-tol )
THEN
292 vmax = ( stmp( i )-tol ) / ( sin( i )+tolin )
296 IF( vmax.GT.rmax( 2 ) )
THEN
298 IF( ninfo( 2 ).EQ.0 )
307 IF( v.GT.septmp( i )*stmp( i ) )
THEN
312 IF( v.GT.sepin( i )*sin( i ) )
THEN
317 tol = max( tol, smlnum / eps )
318 tolin = max( tolin, smlnum / eps )
319 IF( eps*( sepin( i )-tolin ).GT.septmp( i )+tol )
THEN
321 ELSE IF( sepin( i )-tolin.GT.septmp( i )+tol )
THEN
322 vmax = ( sepin( i )-tolin ) / ( septmp( i )+tol )
323 ELSE IF( sepin( i )+tolin.LT.eps*( septmp( i )-tol ) )
THEN
325 ELSE IF( sepin( i )+tolin.LT.septmp( i )-tol )
THEN
326 vmax = ( septmp( i )-tol ) / ( sepin( i )+tolin )
330 IF( vmax.GT.rmax( 2 ) )
THEN
332 IF( ninfo( 2 ).EQ.0 )
341 IF( sin( i ).LE.real( 2*n )*eps .AND. stmp( i ).LE.
342 $ real( 2*n )*eps )
THEN
344 ELSE IF( eps*sin( i ).GT.stmp( i ) )
THEN
346 ELSE IF( sin( i ).GT.stmp( i ) )
THEN
347 vmax = sin( i ) / stmp( i )
348 ELSE IF( sin( i ).LT.eps*stmp( i ) )
THEN
350 ELSE IF( sin( i ).LT.stmp( i ) )
THEN
351 vmax = stmp( i ) / sin( i )
355 IF( vmax.GT.rmax( 3 ) )
THEN
357 IF( ninfo( 3 ).EQ.0 )
366 IF( sepin( i ).LE.v .AND. septmp( i ).LE.v )
THEN
368 ELSE IF( eps*sepin( i ).GT.septmp( i ) )
THEN
370 ELSE IF( sepin( i ).GT.septmp( i ) )
THEN
371 vmax = sepin( i ) / septmp( i )
372 ELSE IF( sepin( i ).LT.eps*septmp( i ) )
THEN
374 ELSE IF( sepin( i ).LT.septmp( i ) )
THEN
375 vmax = septmp( i ) / sepin( i )
379 IF( vmax.GT.rmax( 3 ) )
THEN
381 IF( ninfo( 3 ).EQ.0 )
390 CALL scopy( n, dum, 0, stmp, 1 )
391 CALL scopy( n, dum, 0, septmp, 1 )
392 CALL strsna(
'Eigcond',
'All',
SELECT, n, t, ldt, le, ldt, re,
393 $ ldt, stmp, septmp, n, m, work, n, iwork, info )
396 ninfo( 3 ) = ninfo( 3 ) + 1
400 IF( stmp( i ).NE.s( i ) )
402 IF( septmp( i ).NE.dum( 1 ) )
408 CALL scopy( n, dum, 0, stmp, 1 )
409 CALL scopy( n, dum, 0, septmp, 1 )
410 CALL strsna(
'Veccond',
'All',
SELECT, n, t, ldt, le, ldt, re,
411 $ ldt, stmp, septmp, n, m, work, n, iwork, info )
414 ninfo( 3 ) = ninfo( 3 ) + 1
418 IF( stmp( i ).NE.dum( 1 ) )
420 IF( septmp( i ).NE.sep( i ) )
429 CALL scopy( n, dum, 0, stmp, 1 )
430 CALL scopy( n, dum, 0, septmp, 1 )
431 CALL strsna(
'Bothcond',
'Some',
SELECT, n, t, ldt, le, ldt,
432 $ re, ldt, stmp, septmp, n, m, work, n, iwork,
436 ninfo( 3 ) = ninfo( 3 ) + 1
440 IF( septmp( i ).NE.sep( i ) )
442 IF( stmp( i ).NE.s( i ) )
448 CALL scopy( n, dum, 0, stmp, 1 )
449 CALL scopy( n, dum, 0, septmp, 1 )
450 CALL strsna(
'Eigcond',
'Some',
SELECT, n, t, ldt, le, ldt, re,
451 $ ldt, stmp, septmp, n, m, work, n, iwork, info )
454 ninfo( 3 ) = ninfo( 3 ) + 1
458 IF( stmp( i ).NE.s( i ) )
460 IF( septmp( i ).NE.dum( 1 ) )
466 CALL scopy( n, dum, 0, stmp, 1 )
467 CALL scopy( n, dum, 0, septmp, 1 )
468 CALL strsna(
'Veccond',
'Some',
SELECT, n, t, ldt, le, ldt, re,
469 $ ldt, stmp, septmp, n, m, work, n, iwork, info )
472 ninfo( 3 ) = ninfo( 3 ) + 1
476 IF( stmp( i ).NE.dum( 1 ) )
478 IF( septmp( i ).NE.sep( i ) )
481 IF( vmax.GT.rmax( 1 ) )
THEN
483 IF( ninfo( 1 ).EQ.0 )
489 IF( wi( 1 ).EQ.zero )
THEN
493 IF( ifnd.EQ.1 .OR. wi( i ).EQ.zero )
THEN
494 SELECT( i ) = .false.
499 CALL scopy( n, re( 1, i ), 1, re( 1, 2 ), 1 )
500 CALL scopy( n, re( 1, i+1 ), 1, re( 1, 3 ), 1 )
501 CALL scopy( n, le( 1, i ), 1, le( 1, 2 ), 1 )
502 CALL scopy( n, le( 1, i+1 ), 1, le( 1, 3 ), 1 )
515 IF( ifnd.EQ.1 .OR. wi( i ).NE.zero )
THEN
516 SELECT( i ) = .false.
520 CALL scopy( n, re( 1, i ), 1, re( 1, 3 ), 1 )
521 CALL scopy( n, le( 1, i ), 1, le( 1, 3 ), 1 )
533 CALL scopy( icmp, dum, 0, stmp, 1 )
534 CALL scopy( icmp, dum, 0, septmp, 1 )
535 CALL strsna(
'Bothcond',
'Some',
SELECT, n, t, ldt, le, ldt,
536 $ re, ldt, stmp, septmp, n, m, work, n, iwork,
540 ninfo( 3 ) = ninfo( 3 ) + 1
545 IF( septmp( i ).NE.sep( j ) )
547 IF( stmp( i ).NE.s( j ) )
553 CALL scopy( icmp, dum, 0, stmp, 1 )
554 CALL scopy( icmp, dum, 0, septmp, 1 )
555 CALL strsna(
'Eigcond',
'Some',
SELECT, n, t, ldt, le, ldt, re,
556 $ ldt, stmp, septmp, n, m, work, n, iwork, info )
559 ninfo( 3 ) = ninfo( 3 ) + 1
564 IF( stmp( i ).NE.s( j ) )
566 IF( septmp( i ).NE.dum( 1 ) )
572 CALL scopy( icmp, dum, 0, stmp, 1 )
573 CALL scopy( icmp, dum, 0, septmp, 1 )
574 CALL strsna(
'Veccond',
'Some',
SELECT, n, t, ldt, le, ldt, re,
575 $ ldt, stmp, septmp, n, m, work, n, iwork, info )
578 ninfo( 3 ) = ninfo( 3 ) + 1
583 IF( stmp( i ).NE.dum( 1 ) )
585 IF( septmp( i ).NE.sep( j ) )
588 IF( vmax.GT.rmax( 1 ) )
THEN
590 IF( ninfo( 1 ).EQ.0 )
subroutine scopy(n, sx, incx, sy, incy)
SCOPY
subroutine sgehrd(n, ilo, ihi, a, lda, tau, work, lwork, info)
SGEHRD
subroutine shseqr(job, compz, n, ilo, ihi, h, ldh, wr, wi, z, ldz, work, lwork, info)
SHSEQR
subroutine slacpy(uplo, m, n, a, lda, b, ldb)
SLACPY copies all or part of one two-dimensional array to another.
subroutine sscal(n, sa, sx, incx)
SSCAL
subroutine strevc(side, howmny, select, n, t, ldt, vl, ldvl, vr, ldvr, mm, m, work, info)
STREVC
subroutine strsna(job, howmny, select, n, t, ldt, vl, ldvl, vr, ldvr, s, sep, mm, m, work, ldwork, iwork, info)
STRSNA
subroutine sget37(rmax, lmax, ninfo, knt, nin)
SGET37