331 SUBROUTINE zget24( COMP, JTYPE, THRESH, ISEED, NOUNIT, N, A, LDA,
332 $ H, HT, W, WT, WTMP, VS, LDVS, VS1, RCDEIN,
333 $ RCDVIN, NSLCT, ISLCT, ISRT, RESULT, WORK,
334 $ LWORK, RWORK, BWORK, INFO )
342 INTEGER INFO, ISRT, JTYPE, LDA, LDVS, LWORK, N, NOUNIT,
344 DOUBLE PRECISION RCDEIN, RCDVIN, THRESH
348 INTEGER ISEED( 4 ), ISLCT( * )
349 DOUBLE PRECISION RESULT( 17 ), RWORK( * )
350 COMPLEX*16 A( LDA, * ), H( LDA, * ), HT( LDA, * ),
351 $ vs( ldvs, * ), vs1( ldvs, * ), w( * ),
352 $ work( * ), wt( * ), wtmp( * )
358 COMPLEX*16 CZERO, CONE
359 PARAMETER ( CZERO = ( 0.0d+0, 0.0d+0 ),
360 $ cone = ( 1.0d+0, 0.0d+0 ) )
361 DOUBLE PRECISION ZERO, ONE
362 parameter( zero = 0.0d+0, one = 1.0d+0 )
363 DOUBLE PRECISION EPSIN
364 parameter( epsin = 5.9605d-8 )
368 INTEGER I, IINFO, ISORT, ITMP, J, KMIN, KNTEIG, RSUB,
370 DOUBLE PRECISION ANORM, EPS, RCNDE1, RCNDV1, RCONDE, RCONDV,
371 $ smlnum, tol, tolin, ulp, ulpinv, v, vricmp,
380 DOUBLE PRECISION DLAMCH, ZLANGE
381 EXTERNAL ZSLECT, DLAMCH, ZLANGE
387 INTRINSIC abs, dble, dimag, max, min
391 DOUBLE PRECISION SELWI( 20 ), SELWR( 20 )
394 INTEGER SELDIM, SELOPT
397 COMMON / sslct / selopt, seldim, selval, selwr, selwi
404 IF( thresh.LT.zero )
THEN
406 ELSE IF( nounit.LE.0 )
THEN
408 ELSE IF( n.LT.0 )
THEN
410 ELSE IF( lda.LT.1 .OR. lda.LT.n )
THEN
412 ELSE IF( ldvs.LT.1 .OR. ldvs.LT.n )
THEN
414 ELSE IF( lwork.LT.2*n )
THEN
419 CALL xerbla(
'ZGET24', -info )
434 smlnum = dlamch(
'Safe minimum' )
435 ulp = dlamch(
'Precision' )
442 IF( isort.EQ.0 )
THEN
452 CALL zlacpy(
'F', n, n, a, lda, h, lda )
453 CALL zgeesx(
'V', sort, zslect,
'N', n, h, lda, sdim, w, vs,
454 $ ldvs, rconde, rcondv, work, lwork, rwork, bwork,
456 IF( iinfo.NE.0 )
THEN
457 result( 1+rsub ) = ulpinv
458 IF( jtype.NE.22 )
THEN
459 WRITE( nounit, fmt = 9998 )
'ZGEESX1', iinfo, n, jtype,
462 WRITE( nounit, fmt = 9999 )
'ZGEESX1', iinfo, n,
468 IF( isort.EQ.0 )
THEN
469 CALL zcopy( n, w, 1, wtmp, 1 )
474 result( 1+rsub ) = zero
477 IF( h( i, j ).NE.czero )
478 $ result( 1+rsub ) = ulpinv
486 CALL zlacpy(
' ', n, n, a, lda, vs1, ldvs )
490 CALL zgemm(
'No transpose',
'No transpose', n, n, n, cone, vs,
491 $ ldvs, h, lda, czero, ht, lda )
495 CALL zgemm(
'No transpose',
'Conjugate transpose', n, n, n,
496 $ -cone, ht, lda, vs, ldvs, cone, vs1, ldvs )
498 anorm = max( zlange(
'1', n, n, a, lda, rwork ), smlnum )
499 wnorm = zlange(
'1', n, n, vs1, ldvs, rwork )
501 IF( anorm.GT.wnorm )
THEN
502 result( 2+rsub ) = ( wnorm / anorm ) / ( n*ulp )
504 IF( anorm.LT.one )
THEN
505 result( 2+rsub ) = ( min( wnorm, n*anorm ) / anorm ) /
508 result( 2+rsub ) = min( wnorm / anorm, dble( n ) ) /
515 CALL zunt01(
'Columns', n, n, vs, ldvs, work, lwork, rwork,
520 result( 4+rsub ) = zero
522 IF( h( i, i ).NE.w( i ) )
523 $ result( 4+rsub ) = ulpinv
528 CALL zlacpy(
'F', n, n, a, lda, ht, lda )
529 CALL zgeesx(
'N', sort, zslect,
'N', n, ht, lda, sdim, wt, vs,
530 $ ldvs, rconde, rcondv, work, lwork, rwork, bwork,
532 IF( iinfo.NE.0 )
THEN
533 result( 5+rsub ) = ulpinv
534 IF( jtype.NE.22 )
THEN
535 WRITE( nounit, fmt = 9998 )
'ZGEESX2', iinfo, n, jtype,
538 WRITE( nounit, fmt = 9999 )
'ZGEESX2', iinfo, n,
545 result( 5+rsub ) = zero
548 IF( h( i, j ).NE.ht( i, j ) )
549 $ result( 5+rsub ) = ulpinv
555 result( 6+rsub ) = zero
557 IF( w( i ).NE.wt( i ) )
558 $ result( 6+rsub ) = ulpinv
563 IF( isort.EQ.1 )
THEN
567 IF( zslect( w( i ) ) )
568 $ knteig = knteig + 1
570 IF( zslect( w( i+1 ) ) .AND.
571 $ ( .NOT.zslect( w( i ) ) ) )result( 13 ) = ulpinv
575 $ result( 13 ) = ulpinv
583 IF( lwork.GE.( n*( n+1 ) ) / 2 )
THEN
590 CALL zlacpy(
'F', n, n, a, lda, ht, lda )
591 CALL zgeesx(
'V', sort, zslect,
'B', n, ht, lda, sdim1, wt,
592 $ vs1, ldvs, rconde, rcondv, work, lwork, rwork,
594 IF( iinfo.NE.0 )
THEN
595 result( 14 ) = ulpinv
596 result( 15 ) = ulpinv
597 IF( jtype.NE.22 )
THEN
598 WRITE( nounit, fmt = 9998 )
'ZGEESX3', iinfo, n, jtype,
601 WRITE( nounit, fmt = 9999 )
'ZGEESX3', iinfo, n,
611 IF( w( i ).NE.wt( i ) )
612 $ result( 10 ) = ulpinv
614 IF( h( i, j ).NE.ht( i, j ) )
615 $ result( 11 ) = ulpinv
616 IF( vs( i, j ).NE.vs1( i, j ) )
617 $ result( 12 ) = ulpinv
621 $ result( 13 ) = ulpinv
625 CALL zlacpy(
'F', n, n, a, lda, ht, lda )
626 CALL zgeesx(
'N', sort, zslect,
'B', n, ht, lda, sdim1, wt,
627 $ vs1, ldvs, rcnde1, rcndv1, work, lwork, rwork,
629 IF( iinfo.NE.0 )
THEN
630 result( 14 ) = ulpinv
631 result( 15 ) = ulpinv
632 IF( jtype.NE.22 )
THEN
633 WRITE( nounit, fmt = 9998 )
'ZGEESX4', iinfo, n, jtype,
636 WRITE( nounit, fmt = 9999 )
'ZGEESX4', iinfo, n,
645 IF( rcnde1.NE.rconde )
646 $ result( 14 ) = ulpinv
647 IF( rcndv1.NE.rcondv )
648 $ result( 15 ) = ulpinv
653 IF( w( i ).NE.wt( i ) )
654 $ result( 10 ) = ulpinv
656 IF( h( i, j ).NE.ht( i, j ) )
657 $ result( 11 ) = ulpinv
658 IF( vs( i, j ).NE.vs1( i, j ) )
659 $ result( 12 ) = ulpinv
663 $ result( 13 ) = ulpinv
667 CALL zlacpy(
'F', n, n, a, lda, ht, lda )
668 CALL zgeesx(
'V', sort, zslect,
'E', n, ht, lda, sdim1, wt,
669 $ vs1, ldvs, rcnde1, rcndv1, work, lwork, rwork,
671 IF( iinfo.NE.0 )
THEN
672 result( 14 ) = ulpinv
673 IF( jtype.NE.22 )
THEN
674 WRITE( nounit, fmt = 9998 )
'ZGEESX5', iinfo, n, jtype,
677 WRITE( nounit, fmt = 9999 )
'ZGEESX5', iinfo, n,
686 IF( rcnde1.NE.rconde )
687 $ result( 14 ) = ulpinv
692 IF( w( i ).NE.wt( i ) )
693 $ result( 10 ) = ulpinv
695 IF( h( i, j ).NE.ht( i, j ) )
696 $ result( 11 ) = ulpinv
697 IF( vs( i, j ).NE.vs1( i, j ) )
698 $ result( 12 ) = ulpinv
702 $ result( 13 ) = ulpinv
706 CALL zlacpy(
'F', n, n, a, lda, ht, lda )
707 CALL zgeesx(
'N', sort, zslect,
'E', n, ht, lda, sdim1, wt,
708 $ vs1, ldvs, rcnde1, rcndv1, work, lwork, rwork,
710 IF( iinfo.NE.0 )
THEN
711 result( 14 ) = ulpinv
712 IF( jtype.NE.22 )
THEN
713 WRITE( nounit, fmt = 9998 )
'ZGEESX6', iinfo, n, jtype,
716 WRITE( nounit, fmt = 9999 )
'ZGEESX6', iinfo, n,
725 IF( rcnde1.NE.rconde )
726 $ result( 14 ) = ulpinv
731 IF( w( i ).NE.wt( i ) )
732 $ result( 10 ) = ulpinv
734 IF( h( i, j ).NE.ht( i, j ) )
735 $ result( 11 ) = ulpinv
736 IF( vs( i, j ).NE.vs1( i, j ) )
737 $ result( 12 ) = ulpinv
741 $ result( 13 ) = ulpinv
745 CALL zlacpy(
'F', n, n, a, lda, ht, lda )
746 CALL zgeesx(
'V', sort, zslect,
'V', n, ht, lda, sdim1, wt,
747 $ vs1, ldvs, rcnde1, rcndv1, work, lwork, rwork,
749 IF( iinfo.NE.0 )
THEN
750 result( 15 ) = ulpinv
751 IF( jtype.NE.22 )
THEN
752 WRITE( nounit, fmt = 9998 )
'ZGEESX7', iinfo, n, jtype,
755 WRITE( nounit, fmt = 9999 )
'ZGEESX7', iinfo, n,
764 IF( rcndv1.NE.rcondv )
765 $ result( 15 ) = ulpinv
770 IF( w( i ).NE.wt( i ) )
771 $ result( 10 ) = ulpinv
773 IF( h( i, j ).NE.ht( i, j ) )
774 $ result( 11 ) = ulpinv
775 IF( vs( i, j ).NE.vs1( i, j ) )
776 $ result( 12 ) = ulpinv
780 $ result( 13 ) = ulpinv
784 CALL zlacpy(
'F', n, n, a, lda, ht, lda )
785 CALL zgeesx(
'N', sort, zslect,
'V', n, ht, lda, sdim1, wt,
786 $ vs1, ldvs, rcnde1, rcndv1, work, lwork, rwork,
788 IF( iinfo.NE.0 )
THEN
789 result( 15 ) = ulpinv
790 IF( jtype.NE.22 )
THEN
791 WRITE( nounit, fmt = 9998 )
'ZGEESX8', iinfo, n, jtype,
794 WRITE( nounit, fmt = 9999 )
'ZGEESX8', iinfo, n,
803 IF( rcndv1.NE.rcondv )
804 $ result( 15 ) = ulpinv
809 IF( w( i ).NE.wt( i ) )
810 $ result( 10 ) = ulpinv
812 IF( h( i, j ).NE.ht( i, j ) )
813 $ result( 11 ) = ulpinv
814 IF( vs( i, j ).NE.vs1( i, j ) )
815 $ result( 12 ) = ulpinv
819 $ result( 13 ) = ulpinv
836 eps = max( ulp, epsin )
839 selval( i ) = .false.
840 selwr( i ) = dble( wtmp( i ) )
841 selwi( i ) = dimag( wtmp( i ) )
846 vrimin = dble( wtmp( i ) )
848 vrimin = dimag( wtmp( i ) )
852 vricmp = dble( wtmp( j ) )
854 vricmp = dimag( wtmp( j ) )
856 IF( vricmp.LT.vrimin )
THEN
862 wtmp( kmin ) = wtmp( i )
865 ipnt( i ) = ipnt( kmin )
869 selval( ipnt( islct( i ) ) ) = .true.
874 CALL zlacpy(
'F', n, n, a, lda, ht, lda )
875 CALL zgeesx(
'N',
'S', zslect,
'B', n, ht, lda, sdim1, wt, vs1,
876 $ ldvs, rconde, rcondv, work, lwork, rwork, bwork,
878 IF( iinfo.NE.0 )
THEN
879 result( 16 ) = ulpinv
880 result( 17 ) = ulpinv
881 WRITE( nounit, fmt = 9999 )
'ZGEESX9', iinfo, n, iseed( 1 )
889 anorm = zlange(
'1', n, n, a, lda, rwork )
890 v = max( dble( n )*eps*anorm, smlnum )
893 IF( v.GT.rcondv )
THEN
898 IF( v.GT.rcdvin )
THEN
903 tol = max( tol, smlnum / eps )
904 tolin = max( tolin, smlnum / eps )
905 IF( eps*( rcdein-tolin ).GT.rconde+tol )
THEN
906 result( 16 ) = ulpinv
907 ELSE IF( rcdein-tolin.GT.rconde+tol )
THEN
908 result( 16 ) = ( rcdein-tolin ) / ( rconde+tol )
909 ELSE IF( rcdein+tolin.LT.eps*( rconde-tol ) )
THEN
910 result( 16 ) = ulpinv
911 ELSE IF( rcdein+tolin.LT.rconde-tol )
THEN
912 result( 16 ) = ( rconde-tol ) / ( rcdein+tolin )
920 IF( v.GT.rcondv*rconde )
THEN
925 IF( v.GT.rcdvin*rcdein )
THEN
930 tol = max( tol, smlnum / eps )
931 tolin = max( tolin, smlnum / eps )
932 IF( eps*( rcdvin-tolin ).GT.rcondv+tol )
THEN
933 result( 17 ) = ulpinv
934 ELSE IF( rcdvin-tolin.GT.rcondv+tol )
THEN
935 result( 17 ) = ( rcdvin-tolin ) / ( rcondv+tol )
936 ELSE IF( rcdvin+tolin.LT.eps*( rcondv-tol ) )
THEN
937 result( 17 ) = ulpinv
938 ELSE IF( rcdvin+tolin.LT.rcondv-tol )
THEN
939 result( 17 ) = ( rcondv-tol ) / ( rcdvin+tolin )
948 9999
FORMAT(
' ZGET24: ', a,
' returned INFO=', i6,
'.', / 9x,
'N=',
949 $ i6,
', INPUT EXAMPLE NUMBER = ', i4 )
950 9998
FORMAT(
' ZGET24: ', a,
' returned INFO=', i6,
'.', / 9x,
'N=',
951 $ i6,
', JTYPE=', i6,
', ISEED=(', 3( i5,
',' ), i5,
')' )
subroutine xerbla(srname, info)
subroutine zcopy(n, zx, incx, zy, incy)
ZCOPY
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 ...
subroutine zgemm(transa, transb, m, n, k, alpha, a, lda, b, ldb, beta, c, ldc)
ZGEMM
subroutine zlacpy(uplo, m, n, a, lda, b, ldb)
ZLACPY copies all or part of one two-dimensional array to another.
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 zunt01(rowcol, m, n, u, ldu, work, lwork, rwork, resid)
ZUNT01