333 SUBROUTINE zget24( COMP, JTYPE, THRESH, ISEED, NOUNIT, N, A, LDA,
334 $ h, ht, w, wt, wtmp, vs, ldvs, vs1, rcdein,
335 $ rcdvin, nslct, islct, isrt, result, work,
336 $ lwork, rwork, bwork, info )
345 INTEGER INFO, ISRT, JTYPE, LDA, LDVS, LWORK, N, NOUNIT,
347 DOUBLE PRECISION RCDEIN, RCDVIN, THRESH
351 INTEGER ISEED( 4 ), ISLCT( * )
352 DOUBLE PRECISION RESULT( 17 ), RWORK( * )
353 COMPLEX*16 A( lda, * ), H( lda, * ), HT( lda, * ),
354 $ vs( ldvs, * ), vs1( ldvs, * ), w( * ),
355 $ work( * ), wt( * ), wtmp( * )
361 COMPLEX*16 CZERO, CONE
362 parameter ( czero = ( 0.0d+0, 0.0d+0 ),
363 $ cone = ( 1.0d+0, 0.0d+0 ) )
364 DOUBLE PRECISION ZERO, ONE
365 parameter ( zero = 0.0d+0, one = 1.0d+0 )
366 DOUBLE PRECISION EPSIN
367 parameter ( epsin = 5.9605d-8 )
371 INTEGER I, IINFO, ISORT, ITMP, J, KMIN, KNTEIG, RSUB,
373 DOUBLE PRECISION ANORM, EPS, RCNDE1, RCNDV1, RCONDE, RCONDV,
374 $ smlnum, tol, tolin, ulp, ulpinv, v, vricmp,
383 DOUBLE PRECISION DLAMCH, ZLANGE
384 EXTERNAL zslect, dlamch, zlange
390 INTRINSIC abs, dble, dimag, max, min
394 DOUBLE PRECISION SELWI( 20 ), SELWR( 20 )
397 INTEGER SELDIM, SELOPT
400 COMMON / sslct / selopt, seldim, selval, selwr, selwi
407 IF( thresh.LT.zero )
THEN
409 ELSE IF( nounit.LE.0 )
THEN
411 ELSE IF( n.LT.0 )
THEN
413 ELSE IF( lda.LT.1 .OR. lda.LT.n )
THEN
415 ELSE IF( ldvs.LT.1 .OR. ldvs.LT.n )
THEN
417 ELSE IF( lwork.LT.2*n )
THEN
422 CALL xerbla(
'ZGET24', -info )
437 smlnum = dlamch(
'Safe minimum' )
438 ulp = dlamch(
'Precision' )
445 IF( isort.EQ.0 )
THEN
455 CALL zlacpy(
'F', n, n, a, lda, h, lda )
456 CALL zgeesx(
'V', sort, zslect,
'N', n, h, lda, sdim, w, vs,
457 $ ldvs, rconde, rcondv, work, lwork, rwork, bwork,
459 IF( iinfo.NE.0 )
THEN
460 result( 1+rsub ) = ulpinv
461 IF( jtype.NE.22 )
THEN
462 WRITE( nounit, fmt = 9998 )
'ZGEESX1', iinfo, n, jtype,
465 WRITE( nounit, fmt = 9999 )
'ZGEESX1', iinfo, n,
471 IF( isort.EQ.0 )
THEN
472 CALL zcopy( n, w, 1, wtmp, 1 )
477 result( 1+rsub ) = zero
480 IF( h( i, j ).NE.czero )
481 $ result( 1+rsub ) = ulpinv
489 CALL zlacpy(
' ', n, n, a, lda, vs1, ldvs )
493 CALL zgemm(
'No transpose',
'No transpose', n, n, n, cone, vs,
494 $ ldvs, h, lda, czero, ht, lda )
498 CALL zgemm(
'No transpose',
'Conjugate transpose', n, n, n,
499 $ -cone, ht, lda, vs, ldvs, cone, vs1, ldvs )
501 anorm = max( zlange(
'1', n, n, a, lda, rwork ), smlnum )
502 wnorm = zlange(
'1', n, n, vs1, ldvs, rwork )
504 IF( anorm.GT.wnorm )
THEN
505 result( 2+rsub ) = ( wnorm / anorm ) / ( n*ulp )
507 IF( anorm.LT.one )
THEN
508 result( 2+rsub ) = ( min( wnorm, n*anorm ) / anorm ) /
511 result( 2+rsub ) = min( wnorm / anorm, dble( n ) ) /
518 CALL zunt01(
'Columns', n, n, vs, ldvs, work, lwork, rwork,
523 result( 4+rsub ) = zero
525 IF( h( i, i ).NE.w( i ) )
526 $ result( 4+rsub ) = ulpinv
531 CALL zlacpy(
'F', n, n, a, lda, ht, lda )
532 CALL zgeesx(
'N', sort, zslect,
'N', n, ht, lda, sdim, wt, vs,
533 $ ldvs, rconde, rcondv, work, lwork, rwork, bwork,
535 IF( iinfo.NE.0 )
THEN
536 result( 5+rsub ) = ulpinv
537 IF( jtype.NE.22 )
THEN
538 WRITE( nounit, fmt = 9998 )
'ZGEESX2', iinfo, n, jtype,
541 WRITE( nounit, fmt = 9999 )
'ZGEESX2', iinfo, n,
548 result( 5+rsub ) = zero
551 IF( h( i, j ).NE.ht( i, j ) )
552 $ result( 5+rsub ) = ulpinv
558 result( 6+rsub ) = zero
560 IF( w( i ).NE.wt( i ) )
561 $ result( 6+rsub ) = ulpinv
566 IF( isort.EQ.1 )
THEN
570 IF( zslect( w( i ) ) )
571 $ knteig = knteig + 1
573 IF( zslect( w( i+1 ) ) .AND.
574 $ ( .NOT.zslect( w( i ) ) ) )result( 13 ) = ulpinv
578 $ result( 13 ) = ulpinv
586 IF( lwork.GE.( n*( n+1 ) ) / 2 )
THEN
593 CALL zlacpy(
'F', n, n, a, lda, ht, lda )
594 CALL zgeesx(
'V', sort, zslect,
'B', n, ht, lda, sdim1, wt,
595 $ vs1, ldvs, rconde, rcondv, work, lwork, rwork,
597 IF( iinfo.NE.0 )
THEN
598 result( 14 ) = ulpinv
599 result( 15 ) = ulpinv
600 IF( jtype.NE.22 )
THEN
601 WRITE( nounit, fmt = 9998 )
'ZGEESX3', iinfo, n, jtype,
604 WRITE( nounit, fmt = 9999 )
'ZGEESX3', iinfo, n,
614 IF( w( i ).NE.wt( i ) )
615 $ result( 10 ) = ulpinv
617 IF( h( i, j ).NE.ht( i, j ) )
618 $ result( 11 ) = ulpinv
619 IF( vs( i, j ).NE.vs1( i, j ) )
620 $ result( 12 ) = ulpinv
624 $ result( 13 ) = ulpinv
628 CALL zlacpy(
'F', n, n, a, lda, ht, lda )
629 CALL zgeesx(
'N', sort, zslect,
'B', n, ht, lda, sdim1, wt,
630 $ vs1, ldvs, rcnde1, rcndv1, work, lwork, rwork,
632 IF( iinfo.NE.0 )
THEN
633 result( 14 ) = ulpinv
634 result( 15 ) = ulpinv
635 IF( jtype.NE.22 )
THEN
636 WRITE( nounit, fmt = 9998 )
'ZGEESX4', iinfo, n, jtype,
639 WRITE( nounit, fmt = 9999 )
'ZGEESX4', iinfo, n,
648 IF( rcnde1.NE.rconde )
649 $ result( 14 ) = ulpinv
650 IF( rcndv1.NE.rcondv )
651 $ result( 15 ) = ulpinv
656 IF( w( i ).NE.wt( i ) )
657 $ result( 10 ) = ulpinv
659 IF( h( i, j ).NE.ht( i, j ) )
660 $ result( 11 ) = ulpinv
661 IF( vs( i, j ).NE.vs1( i, j ) )
662 $ result( 12 ) = ulpinv
666 $ result( 13 ) = ulpinv
670 CALL zlacpy(
'F', n, n, a, lda, ht, lda )
671 CALL zgeesx(
'V', sort, zslect,
'E', n, ht, lda, sdim1, wt,
672 $ vs1, ldvs, rcnde1, rcndv1, work, lwork, rwork,
674 IF( iinfo.NE.0 )
THEN
675 result( 14 ) = ulpinv
676 IF( jtype.NE.22 )
THEN
677 WRITE( nounit, fmt = 9998 )
'ZGEESX5', iinfo, n, jtype,
680 WRITE( nounit, fmt = 9999 )
'ZGEESX5', iinfo, n,
689 IF( rcnde1.NE.rconde )
690 $ result( 14 ) = ulpinv
695 IF( w( i ).NE.wt( i ) )
696 $ result( 10 ) = ulpinv
698 IF( h( i, j ).NE.ht( i, j ) )
699 $ result( 11 ) = ulpinv
700 IF( vs( i, j ).NE.vs1( i, j ) )
701 $ result( 12 ) = ulpinv
705 $ result( 13 ) = ulpinv
709 CALL zlacpy(
'F', n, n, a, lda, ht, lda )
710 CALL zgeesx(
'N', sort, zslect,
'E', n, ht, lda, sdim1, wt,
711 $ vs1, ldvs, rcnde1, rcndv1, work, lwork, rwork,
713 IF( iinfo.NE.0 )
THEN
714 result( 14 ) = ulpinv
715 IF( jtype.NE.22 )
THEN
716 WRITE( nounit, fmt = 9998 )
'ZGEESX6', iinfo, n, jtype,
719 WRITE( nounit, fmt = 9999 )
'ZGEESX6', iinfo, n,
728 IF( rcnde1.NE.rconde )
729 $ result( 14 ) = ulpinv
734 IF( w( i ).NE.wt( i ) )
735 $ result( 10 ) = ulpinv
737 IF( h( i, j ).NE.ht( i, j ) )
738 $ result( 11 ) = ulpinv
739 IF( vs( i, j ).NE.vs1( i, j ) )
740 $ result( 12 ) = ulpinv
744 $ result( 13 ) = ulpinv
748 CALL zlacpy(
'F', n, n, a, lda, ht, lda )
749 CALL zgeesx(
'V', sort, zslect,
'V', n, ht, lda, sdim1, wt,
750 $ vs1, ldvs, rcnde1, rcndv1, work, lwork, rwork,
752 IF( iinfo.NE.0 )
THEN
753 result( 15 ) = ulpinv
754 IF( jtype.NE.22 )
THEN
755 WRITE( nounit, fmt = 9998 )
'ZGEESX7', iinfo, n, jtype,
758 WRITE( nounit, fmt = 9999 )
'ZGEESX7', iinfo, n,
767 IF( rcndv1.NE.rcondv )
768 $ result( 15 ) = ulpinv
773 IF( w( i ).NE.wt( i ) )
774 $ result( 10 ) = ulpinv
776 IF( h( i, j ).NE.ht( i, j ) )
777 $ result( 11 ) = ulpinv
778 IF( vs( i, j ).NE.vs1( i, j ) )
779 $ result( 12 ) = ulpinv
783 $ result( 13 ) = ulpinv
787 CALL zlacpy(
'F', n, n, a, lda, ht, lda )
788 CALL zgeesx(
'N', sort, zslect,
'V', n, ht, lda, sdim1, wt,
789 $ vs1, ldvs, rcnde1, rcndv1, work, lwork, rwork,
791 IF( iinfo.NE.0 )
THEN
792 result( 15 ) = ulpinv
793 IF( jtype.NE.22 )
THEN
794 WRITE( nounit, fmt = 9998 )
'ZGEESX8', iinfo, n, jtype,
797 WRITE( nounit, fmt = 9999 )
'ZGEESX8', iinfo, n,
806 IF( rcndv1.NE.rcondv )
807 $ result( 15 ) = ulpinv
812 IF( w( i ).NE.wt( i ) )
813 $ result( 10 ) = ulpinv
815 IF( h( i, j ).NE.ht( i, j ) )
816 $ result( 11 ) = ulpinv
817 IF( vs( i, j ).NE.vs1( i, j ) )
818 $ result( 12 ) = ulpinv
822 $ result( 13 ) = ulpinv
839 eps = max( ulp, epsin )
842 selval( i ) = .false.
843 selwr( i ) = dble( wtmp( i ) )
844 selwi( i ) = dimag( wtmp( i ) )
849 vrimin = dble( wtmp( i ) )
851 vrimin = dimag( wtmp( i ) )
855 vricmp = dble( wtmp( j ) )
857 vricmp = dimag( wtmp( j ) )
859 IF( vricmp.LT.vrimin )
THEN
865 wtmp( kmin ) = wtmp( i )
868 ipnt( i ) = ipnt( kmin )
872 selval( ipnt( islct( i ) ) ) = .true.
877 CALL zlacpy(
'F', n, n, a, lda, ht, lda )
878 CALL zgeesx(
'N',
'S', zslect,
'B', n, ht, lda, sdim1, wt, vs1,
879 $ ldvs, rconde, rcondv, work, lwork, rwork, bwork,
881 IF( iinfo.NE.0 )
THEN
882 result( 16 ) = ulpinv
883 result( 17 ) = ulpinv
884 WRITE( nounit, fmt = 9999 )
'ZGEESX9', iinfo, n, iseed( 1 )
892 anorm = zlange(
'1', n, n, a, lda, rwork )
893 v = max( dble( n )*eps*anorm, smlnum )
896 IF( v.GT.rcondv )
THEN
901 IF( v.GT.rcdvin )
THEN
906 tol = max( tol, smlnum / eps )
907 tolin = max( tolin, smlnum / eps )
908 IF( eps*( rcdein-tolin ).GT.rconde+tol )
THEN
909 result( 16 ) = ulpinv
910 ELSE IF( rcdein-tolin.GT.rconde+tol )
THEN
911 result( 16 ) = ( rcdein-tolin ) / ( rconde+tol )
912 ELSE IF( rcdein+tolin.LT.eps*( rconde-tol ) )
THEN
913 result( 16 ) = ulpinv
914 ELSE IF( rcdein+tolin.LT.rconde-tol )
THEN
915 result( 16 ) = ( rconde-tol ) / ( rcdein+tolin )
923 IF( v.GT.rcondv*rconde )
THEN
928 IF( v.GT.rcdvin*rcdein )
THEN
933 tol = max( tol, smlnum / eps )
934 tolin = max( tolin, smlnum / eps )
935 IF( eps*( rcdvin-tolin ).GT.rcondv+tol )
THEN
936 result( 17 ) = ulpinv
937 ELSE IF( rcdvin-tolin.GT.rcondv+tol )
THEN
938 result( 17 ) = ( rcdvin-tolin ) / ( rcondv+tol )
939 ELSE IF( rcdvin+tolin.LT.eps*( rcondv-tol ) )
THEN
940 result( 17 ) = ulpinv
941 ELSE IF( rcdvin+tolin.LT.rcondv-tol )
THEN
942 result( 17 ) = ( rcondv-tol ) / ( rcdvin+tolin )
951 9999
FORMAT(
' ZGET24: ', a,
' returned INFO=', i6,
'.', / 9x,
'N=',
952 $ i6,
', INPUT EXAMPLE NUMBER = ', i4 )
953 9998
FORMAT(
' ZGET24: ', a,
' returned INFO=', i6,
'.', / 9x,
'N=',
954 $ i6,
', JTYPE=', i6,
', ISEED=(', 3( i5,
',' ), i5,
')' )
subroutine zlacpy(UPLO, M, N, A, LDA, B, LDB)
ZLACPY copies all or part of one two-dimensional array to another.
subroutine zcopy(N, ZX, INCX, ZY, INCY)
ZCOPY
subroutine zgemm(TRANSA, TRANSB, M, N, K, ALPHA, A, LDA, B, LDB, BETA, C, LDC)
ZGEMM
subroutine zget24(COMP, JTYPE, THRESH, ISEED, NOUNIT, N, A, LDA, H, HT, W, WT, WTMP, VS, LDVS, VS1, RCDEIN, RCDVIN, NSLCT, ISLCT, ISRT, RESULT, WORK, LWORK, RWORK, BWORK, INFO)
ZGET24
subroutine xerbla(SRNAME, INFO)
XERBLA
subroutine zunt01(ROWCOL, M, N, U, LDU, WORK, LWORK, RWORK, RESID)
ZUNT01
subroutine zgeesx(JOBVS, SORT, SELECT, SENSE, N, A, LDA, SDIM, W, VS, LDVS, RCONDE, RCONDV, WORK, LWORK, RWORK, BWORK, INFO)
ZGEESX computes the eigenvalues, the Schur form, and, optionally, the matrix of Schur vectors for GE...