91 SUBROUTINE cget37( RMAX, LMAX, NINFO, KNT, NIN )
102 INTEGER lmax( 3 ), ninfo( 3 )
110 parameter( zero = 0.0e0, one = 1.0e0, two = 2.0e0 )
112 parameter( epsin = 5.9605e-8 )
114 parameter( ldt = 20, lwork = 2*ldt*( 10+ldt ) )
117 INTEGER i, icmp, info, iscl, isrt,
j, kmin, m, n
118 REAL bignum, eps, smlnum, tnrm, tol, tolin, v,
119 $ vcmin, vmax, vmin, vmul
122 LOGICAL select( ldt )
124 REAL dum( 1 ), rwork( 2*ldt ), s( ldt ), sep( ldt ),
125 $ sepin( ldt ), septmp( ldt ), sin( ldt ),
126 $ stmp( ldt ), val( 3 ), wiin( ldt ),
127 $ wrin( ldt ), wsrt( ldt )
128 COMPLEX cdum( 1 ), le( ldt, ldt ), re( ldt, ldt ),
129 $ t( ldt, ldt ), tmp( ldt, ldt ), w( ldt ),
130 $ work( lwork ), wtmp( ldt )
141 INTRINSIC aimag, max,
REAL, sqrt
146 smlnum =
slamch(
'S' ) / eps
147 bignum = one / smlnum
148 CALL
slabad( smlnum, bignum )
152 eps = max( eps, epsin )
163 val( 1 ) = sqrt( smlnum )
165 val( 3 ) = sqrt( bignum )
172 READ( nin, fmt = * )n, isrt
176 READ( nin, fmt = * )( tmp( i,
j ),
j = 1, n )
179 READ( nin, fmt = * )wrin( i ), wiin( i ), sin( i ), sepin( i )
181 tnrm =
clange(
'M', n, n, tmp, ldt, rwork )
187 CALL
clacpy(
'F', n, n, tmp, ldt, t, ldt )
190 CALL
csscal( n, vmul, t( 1, i ), 1 )
197 CALL
cgehrd( n, 1, n, t, ldt, work( 1 ), work( n+1 ), lwork-n,
201 ninfo( 1 ) = ninfo( 1 ) + 1
212 CALL
chseqr(
'S',
'N', n, 1, n, t, ldt, w, cdum, 1, work,
216 ninfo( 2 ) = ninfo( 2 ) + 1
225 CALL
ctrevc(
'B',
'A',
SELECT, n, t, ldt, le, ldt, re, ldt, n,
226 $ m, work, rwork, info )
230 CALL
ctrsna(
'B',
'A',
SELECT, n, t, ldt, le, ldt, re, ldt, s,
231 $ sep, n, m, work, n, rwork, info )
234 ninfo( 3 ) = ninfo( 3 ) + 1
241 CALL
ccopy( n, w, 1, wtmp, 1 )
247 wsrt( i ) =
REAL( W( I ) )
254 wsrt( i ) = aimag( w( i ) )
257 CALL
scopy( n, s, 1, stmp, 1 )
258 CALL
scopy( n, sep, 1, septmp, 1 )
259 CALL
sscal( n, one / vmul, septmp, 1 )
264 IF( wsrt(
j ).LT.vmin )
THEN
269 wsrt( kmin ) = wsrt( i )
272 wtmp( i ) = w( kmin )
275 stmp( kmin ) = stmp( i )
277 vmin = septmp( kmin )
278 septmp( kmin ) = septmp( i )
285 v = max( two*
REAL( n )*eps*tnrm, smlnum )
289 IF( v.GT.septmp( i ) )
THEN
292 tol = v / septmp( i )
294 IF( v.GT.sepin( i ) )
THEN
297 tolin = v / sepin( i )
299 tol = max( tol, smlnum / eps )
300 tolin = max( tolin, smlnum / eps )
301 IF( eps*( sin( i )-tolin ).GT.stmp( i )+tol )
THEN
303 ELSE IF( sin( i )-tolin.GT.stmp( i )+tol )
THEN
304 vmax = ( sin( i )-tolin ) / ( stmp( i )+tol )
305 ELSE IF( sin( i )+tolin.LT.eps*( stmp( i )-tol ) )
THEN
307 ELSE IF( sin( i )+tolin.LT.stmp( i )-tol )
THEN
308 vmax = ( stmp( i )-tol ) / ( sin( i )+tolin )
312 IF( vmax.GT.rmax( 2 ) )
THEN
314 IF( ninfo( 2 ).EQ.0 )
323 IF( v.GT.septmp( i )*stmp( i ) )
THEN
328 IF( v.GT.sepin( i )*sin( i ) )
THEN
333 tol = max( tol, smlnum / eps )
334 tolin = max( tolin, smlnum / eps )
335 IF( eps*( sepin( i )-tolin ).GT.septmp( i )+tol )
THEN
337 ELSE IF( sepin( i )-tolin.GT.septmp( i )+tol )
THEN
338 vmax = ( sepin( i )-tolin ) / ( septmp( i )+tol )
339 ELSE IF( sepin( i )+tolin.LT.eps*( septmp( i )-tol ) )
THEN
341 ELSE IF( sepin( i )+tolin.LT.septmp( i )-tol )
THEN
342 vmax = ( septmp( i )-tol ) / ( sepin( i )+tolin )
346 IF( vmax.GT.rmax( 2 ) )
THEN
348 IF( ninfo( 2 ).EQ.0 )
357 IF( sin( i ).LE.
REAL( 2*n )*eps .AND. stmp( i ).LE.
358 $
REAL( 2*n )*eps ) then
360 ELSE IF( eps*sin( i ).GT.stmp( i ) )
THEN
362 ELSE IF( sin( i ).GT.stmp( i ) )
THEN
363 vmax = sin( i ) / stmp( i )
364 ELSE IF( sin( i ).LT.eps*stmp( i ) )
THEN
366 ELSE IF( sin( i ).LT.stmp( i ) )
THEN
367 vmax = stmp( i ) / sin( i )
371 IF( vmax.GT.rmax( 3 ) )
THEN
373 IF( ninfo( 3 ).EQ.0 )
382 IF( sepin( i ).LE.v .AND. septmp( i ).LE.v )
THEN
384 ELSE IF( eps*sepin( i ).GT.septmp( i ) )
THEN
386 ELSE IF( sepin( i ).GT.septmp( i ) )
THEN
387 vmax = sepin( i ) / septmp( i )
388 ELSE IF( sepin( i ).LT.eps*septmp( i ) )
THEN
390 ELSE IF( sepin( i ).LT.septmp( i ) )
THEN
391 vmax = septmp( i ) / sepin( i )
395 IF( vmax.GT.rmax( 3 ) )
THEN
397 IF( ninfo( 3 ).EQ.0 )
406 CALL
scopy( n, dum, 0, stmp, 1 )
407 CALL
scopy( n, dum, 0, septmp, 1 )
408 CALL
ctrsna(
'E',
'A',
SELECT, n, t, ldt, le, ldt, re, ldt,
409 $ stmp, septmp, n, m, work, n, rwork, info )
412 ninfo( 3 ) = ninfo( 3 ) + 1
416 IF( stmp( i ).NE.s( i ) )
418 IF( septmp( i ).NE.dum( 1 ) )
424 CALL
scopy( n, dum, 0, stmp, 1 )
425 CALL
scopy( n, dum, 0, septmp, 1 )
426 CALL
ctrsna(
'V',
'A',
SELECT, n, t, ldt, le, ldt, re, ldt,
427 $ stmp, septmp, n, m, work, n, rwork, info )
430 ninfo( 3 ) = ninfo( 3 ) + 1
434 IF( stmp( i ).NE.dum( 1 ) )
436 IF( septmp( i ).NE.sep( i ) )
445 CALL
scopy( n, dum, 0, stmp, 1 )
446 CALL
scopy( n, dum, 0, septmp, 1 )
447 CALL
ctrsna(
'B',
'S',
SELECT, n, t, ldt, le, ldt, re, ldt,
448 $ stmp, septmp, n, m, work, n, rwork, info )
451 ninfo( 3 ) = ninfo( 3 ) + 1
455 IF( septmp( i ).NE.sep( i ) )
457 IF( stmp( i ).NE.s( i ) )
463 CALL
scopy( n, dum, 0, stmp, 1 )
464 CALL
scopy( n, dum, 0, septmp, 1 )
465 CALL
ctrsna(
'E',
'S',
SELECT, n, t, ldt, le, ldt, re, ldt,
466 $ stmp, septmp, n, m, work, n, rwork, info )
469 ninfo( 3 ) = ninfo( 3 ) + 1
473 IF( stmp( i ).NE.s( i ) )
475 IF( septmp( i ).NE.dum( 1 ) )
481 CALL
scopy( n, dum, 0, stmp, 1 )
482 CALL
scopy( n, dum, 0, septmp, 1 )
483 CALL
ctrsna(
'V',
'S',
SELECT, n, t, ldt, le, ldt, re, ldt,
484 $ stmp, septmp, n, m, work, n, rwork, info )
487 ninfo( 3 ) = ninfo( 3 ) + 1
491 IF( stmp( i ).NE.dum( 1 ) )
493 IF( septmp( i ).NE.sep( i ) )
496 IF( vmax.GT.rmax( 1 ) )
THEN
498 IF( ninfo( 1 ).EQ.0 )
505 SELECT( i ) = .false.
512 CALL
ccopy( n, re( 1, 2 ), 1, re( 1, 1 ), 1 )
513 CALL
ccopy( n, le( 1, 2 ), 1, le( 1, 1 ), 1 )
518 SELECT( n-1 ) = .true.
519 CALL
ccopy( n, re( 1, n-1 ), 1, re( 1, 2 ), 1 )
520 CALL
ccopy( n, le( 1, n-1 ), 1, le( 1, 2 ), 1 )
525 CALL
scopy( icmp, dum, 0, stmp, 1 )
526 CALL
scopy( icmp, dum, 0, septmp, 1 )
527 CALL
ctrsna(
'B',
'S',
SELECT, n, t, ldt, le, ldt, re, ldt,
528 $ stmp, septmp, n, m, work, n, rwork, info )
531 ninfo( 3 ) = ninfo( 3 ) + 1
536 IF( septmp( i ).NE.sep(
j ) )
538 IF( stmp( i ).NE.s(
j ) )
544 CALL
scopy( icmp, dum, 0, stmp, 1 )
545 CALL
scopy( icmp, dum, 0, septmp, 1 )
546 CALL
ctrsna(
'E',
'S',
SELECT, n, t, ldt, le, ldt, re, ldt,
547 $ stmp, septmp, n, m, work, n, rwork, info )
550 ninfo( 3 ) = ninfo( 3 ) + 1
555 IF( stmp( i ).NE.s(
j ) )
557 IF( septmp( i ).NE.dum( 1 ) )
563 CALL
scopy( icmp, dum, 0, stmp, 1 )
564 CALL
scopy( icmp, dum, 0, septmp, 1 )
565 CALL
ctrsna(
'V',
'S',
SELECT, n, t, ldt, le, ldt, re, ldt,
566 $ stmp, septmp, n, m, work, n, rwork, info )
569 ninfo( 3 ) = ninfo( 3 ) + 1
574 IF( stmp( i ).NE.dum( 1 ) )
576 IF( septmp( i ).NE.sep(
j ) )
579 IF( vmax.GT.rmax( 1 ) )
THEN
581 IF( ninfo( 1 ).EQ.0 )
real function clange(NORM, M, N, A, LDA, WORK)
CLANGE returns the value of the 1-norm, Frobenius norm, infinity-norm, or the largest absolute value ...
subroutine scopy(N, SX, INCX, SY, INCY)
SCOPY
subroutine ctrsna(JOB, HOWMNY, SELECT, N, T, LDT, VL, LDVL, VR, LDVR, S, SEP, MM, M, WORK, LDWORK, RWORK, INFO)
CTRSNA
subroutine clacpy(UPLO, M, N, A, LDA, B, LDB)
CLACPY copies all or part of one two-dimensional array to another.
subroutine chseqr(JOB, COMPZ, N, ILO, IHI, H, LDH, W, Z, LDZ, WORK, LWORK, INFO)
CHSEQR
real function slamch(CMACH)
SLAMCH
set ue cd $ADTTMP cat<< EOF > tmp f Program LinearEquations Implicit none Real j
subroutine cgehrd(N, ILO, IHI, A, LDA, TAU, WORK, LWORK, INFO)
CGEHRD
subroutine cget37(RMAX, LMAX, NINFO, KNT, NIN)
CGET37
subroutine slabad(SMALL, LARGE)
SLABAD
subroutine ccopy(N, CX, INCX, CY, INCY)
CCOPY
subroutine csscal(N, SA, CX, INCX)
CSSCAL
subroutine ctrevc(SIDE, HOWMNY, SELECT, N, T, LDT, VL, LDVL, VR, LDVR, MM, M, WORK, RWORK, INFO)
CTREVC
subroutine sscal(N, SA, SX, INCX)
SSCAL