295 SUBROUTINE zdrgvx( NSIZE, THRESH, NIN, NOUT, A, LDA, B, AI, BI,
296 $ alpha, beta, vl, vr, ilo, ihi, lscale, rscale,
297 $ s, dtru, dif, diftru, work, lwork, rwork,
298 $ iwork, liwork, result, bwork, info )
306 INTEGER IHI, ILO, INFO, LDA, LIWORK, LWORK, NIN, NOUT,
308 DOUBLE PRECISION THRESH
313 DOUBLE PRECISION DIF( * ), DIFTRU( * ), DTRU( * ), LSCALE( * ),
314 $ result( 4 ), rscale( * ), rwork( * ), s( * )
315 COMPLEX*16 A( lda, * ), AI( lda, * ), ALPHA( * ),
316 $ b( lda, * ), beta( * ), bi( lda, * ),
317 $ vl( lda, * ), vr( lda, * ), work( * )
323 DOUBLE PRECISION ZERO, ONE, TEN, TNTH, HALF
324 parameter ( zero = 0.0d+0, one = 1.0d+0, ten = 1.0d+1,
325 $ tnth = 1.0d-1, half = 0.5d+0 )
328 INTEGER I, IPTYPE, IWA, IWB, IWX, IWY, J, LINFO,
329 $ maxwrk, minwrk, n, nerrs, nmax, nptknt, ntestt
330 DOUBLE PRECISION ABNORM, ANORM, BNORM, RATIO1, RATIO2, THRSH2,
334 COMPLEX*16 WEIGHT( 5 )
338 DOUBLE PRECISION DLAMCH, ZLANGE
339 EXTERNAL ilaenv, dlamch, zlange
345 INTRINSIC abs, dcmplx, max, sqrt
355 IF( nsize.LT.0 )
THEN
357 ELSE IF( thresh.LT.zero )
THEN
359 ELSE IF( nin.LE.0 )
THEN
361 ELSE IF( nout.LE.0 )
THEN
363 ELSE IF( lda.LT.1 .OR. lda.LT.nmax )
THEN
365 ELSE IF( liwork.LT.nmax+2 )
THEN
377 IF( info.EQ.0 .AND. lwork.GE.1 )
THEN
378 minwrk = 2*nmax*( nmax+1 )
379 maxwrk = nmax*( 1+ilaenv( 1,
'ZGEQRF',
' ', nmax, 1, nmax,
381 maxwrk = max( maxwrk, 2*nmax*( nmax+1 ) )
385 IF( lwork.LT.minwrk )
389 CALL xerbla(
'ZDRGVX', -info )
406 weight( 1 ) = dcmplx( tnth, zero )
407 weight( 2 ) = dcmplx( half, zero )
409 weight( 4 ) = one / weight( 2 )
410 weight( 5 ) = one / weight( 1 )
420 CALL zlatm6( iptype, 5, a, lda, b, vr, lda, vl,
421 $ lda, weight( iwa ), weight( iwb ),
422 $ weight( iwx ), weight( iwy ), dtru,
429 CALL zlacpy(
'F', n, n, a, lda, ai, lda )
430 CALL zlacpy(
'F', n, n, b, lda, bi, lda )
432 CALL zggevx(
'N',
'V',
'V',
'B', n, ai, lda, bi,
433 $ lda, alpha, beta, vl, lda, vr, lda,
434 $ ilo, ihi, lscale, rscale, anorm,
435 $ bnorm, s, dif, work, lwork, rwork,
436 $ iwork, bwork, linfo )
437 IF( linfo.NE.0 )
THEN
438 WRITE( nout, fmt = 9999 )
'ZGGEVX', linfo, n,
439 $ iptype, iwa, iwb, iwx, iwy
445 CALL zlacpy(
'Full', n, n, ai, lda, work, n )
446 CALL zlacpy(
'Full', n, n, bi, lda, work( n*n+1 ),
448 abnorm = zlange(
'Fro', n, 2*n, work, n, rwork )
453 CALL zget52( .true., n, a, lda, b, lda, vl, lda,
454 $ alpha, beta, work, rwork,
456 IF( result( 2 ).GT.thresh )
THEN
457 WRITE( nout, fmt = 9998 )
'Left',
'ZGGEVX',
458 $ result( 2 ), n, iptype, iwa, iwb, iwx, iwy
462 CALL zget52( .false., n, a, lda, b, lda, vr, lda,
463 $ alpha, beta, work, rwork,
465 IF( result( 3 ).GT.thresh )
THEN
466 WRITE( nout, fmt = 9998 )
'Right',
'ZGGEVX',
467 $ result( 3 ), n, iptype, iwa, iwb, iwx, iwy
474 IF( s( i ).EQ.zero )
THEN
475 IF( dtru( i ).GT.abnorm*ulp )
476 $ result( 3 ) = ulpinv
477 ELSE IF( dtru( i ).EQ.zero )
THEN
478 IF( s( i ).GT.abnorm*ulp )
479 $ result( 3 ) = ulpinv
481 rwork( i ) = max( abs( dtru( i ) / s( i ) ),
482 $ abs( s( i ) / dtru( i ) ) )
483 result( 3 ) = max( result( 3 ), rwork( i ) )
490 IF( dif( 1 ).EQ.zero )
THEN
491 IF( diftru( 1 ).GT.abnorm*ulp )
492 $ result( 4 ) = ulpinv
493 ELSE IF( diftru( 1 ).EQ.zero )
THEN
494 IF( dif( 1 ).GT.abnorm*ulp )
495 $ result( 4 ) = ulpinv
496 ELSE IF( dif( 5 ).EQ.zero )
THEN
497 IF( diftru( 5 ).GT.abnorm*ulp )
498 $ result( 4 ) = ulpinv
499 ELSE IF( diftru( 5 ).EQ.zero )
THEN
500 IF( dif( 5 ).GT.abnorm*ulp )
501 $ result( 4 ) = ulpinv
503 ratio1 = max( abs( diftru( 1 ) / dif( 1 ) ),
504 $ abs( dif( 1 ) / diftru( 1 ) ) )
505 ratio2 = max( abs( diftru( 5 ) / dif( 5 ) ),
506 $ abs( dif( 5 ) / diftru( 5 ) ) )
507 result( 4 ) = max( ratio1, ratio2 )
515 IF( ( result( j ).GE.thrsh2 .AND. j.GE.4 ) .OR.
516 $ ( result( j ).GE.thresh .AND. j.LE.3 ) )
522 IF( nerrs.EQ.0 )
THEN
523 WRITE( nout, fmt = 9997 )
'ZXV'
529 WRITE( nout, fmt = 9995 )
530 WRITE( nout, fmt = 9994 )
531 WRITE( nout, fmt = 9993 )
535 WRITE( nout, fmt = 9992 )
'''',
540 IF( result( j ).LT.10000.0d0 )
THEN
541 WRITE( nout, fmt = 9991 )iptype, iwa,
542 $ iwb, iwx, iwy, j, result( j )
544 WRITE( nout, fmt = 9990 )iptype, iwa,
545 $ iwb, iwx, iwy, j, result( j )
565 READ( nin, fmt = *, end = 150 )n
569 READ( nin, fmt = * )( a( i, j ), j = 1, n )
572 READ( nin, fmt = * )( b( i, j ), j = 1, n )
574 READ( nin, fmt = * )( dtru( i ), i = 1, n )
575 READ( nin, fmt = * )( diftru( i ), i = 1, n )
583 CALL zlacpy(
'F', n, n, a, lda, ai, lda )
584 CALL zlacpy(
'F', n, n, b, lda, bi, lda )
586 CALL zggevx(
'N',
'V',
'V',
'B', n, ai, lda, bi, lda, alpha, beta,
587 $ vl, lda, vr, lda, ilo, ihi, lscale, rscale, anorm,
588 $ bnorm, s, dif, work, lwork, rwork, iwork, bwork,
591 IF( linfo.NE.0 )
THEN
592 WRITE( nout, fmt = 9987 )
'ZGGEVX', linfo, n, nptknt
598 CALL zlacpy(
'Full', n, n, ai, lda, work, n )
599 CALL zlacpy(
'Full', n, n, bi, lda, work( n*n+1 ), n )
600 abnorm = zlange(
'Fro', n, 2*n, work, n, rwork )
605 CALL zget52( .true., n, a, lda, b, lda, vl, lda, alpha, beta,
606 $ work, rwork, result( 1 ) )
607 IF( result( 2 ).GT.thresh )
THEN
608 WRITE( nout, fmt = 9986 )
'Left',
'ZGGEVX', result( 2 ), n,
613 CALL zget52( .false., n, a, lda, b, lda, vr, lda, alpha, beta,
614 $ work, rwork, result( 2 ) )
615 IF( result( 3 ).GT.thresh )
THEN
616 WRITE( nout, fmt = 9986 )
'Right',
'ZGGEVX', result( 3 ), n,
624 IF( s( i ).EQ.zero )
THEN
625 IF( dtru( i ).GT.abnorm*ulp )
626 $ result( 3 ) = ulpinv
627 ELSE IF( dtru( i ).EQ.zero )
THEN
628 IF( s( i ).GT.abnorm*ulp )
629 $ result( 3 ) = ulpinv
631 rwork( i ) = max( abs( dtru( i ) / s( i ) ),
632 $ abs( s( i ) / dtru( i ) ) )
633 result( 3 ) = max( result( 3 ), rwork( i ) )
640 IF( dif( 1 ).EQ.zero )
THEN
641 IF( diftru( 1 ).GT.abnorm*ulp )
642 $ result( 4 ) = ulpinv
643 ELSE IF( diftru( 1 ).EQ.zero )
THEN
644 IF( dif( 1 ).GT.abnorm*ulp )
645 $ result( 4 ) = ulpinv
646 ELSE IF( dif( 5 ).EQ.zero )
THEN
647 IF( diftru( 5 ).GT.abnorm*ulp )
648 $ result( 4 ) = ulpinv
649 ELSE IF( diftru( 5 ).EQ.zero )
THEN
650 IF( dif( 5 ).GT.abnorm*ulp )
651 $ result( 4 ) = ulpinv
653 ratio1 = max( abs( diftru( 1 ) / dif( 1 ) ),
654 $ abs( dif( 1 ) / diftru( 1 ) ) )
655 ratio2 = max( abs( diftru( 5 ) / dif( 5 ) ),
656 $ abs( dif( 5 ) / diftru( 5 ) ) )
657 result( 4 ) = max( ratio1, ratio2 )
665 IF( result( j ).GE.thrsh2 )
THEN
670 IF( nerrs.EQ.0 )
THEN
671 WRITE( nout, fmt = 9997 )
'ZXV'
677 WRITE( nout, fmt = 9996 )
681 WRITE( nout, fmt = 9992 )
'''',
'transpose',
''''
685 IF( result( j ).LT.10000.0d0 )
THEN
686 WRITE( nout, fmt = 9989 )nptknt, n, j, result( j )
688 WRITE( nout, fmt = 9988 )nptknt, n, j, result( j )
700 CALL alasvm(
'ZXV', nout, nerrs, ntestt, 0 )
706 9999
FORMAT(
' ZDRGVX: ', a,
' returned INFO=', i6,
'.', / 9x,
'N=',
707 $ i6,
', JTYPE=', i6,
')' )
709 9998
FORMAT(
' ZDRGVX: ', a,
' Eigenvectors from ', a,
' incorrectly ',
710 $
'normalized.', /
' Bits of error=', 0p, g10.3,
',', 9x,
711 $
'N=', i6,
', JTYPE=', i6,
', IWA=', i5,
', IWB=', i5,
712 $
', IWX=', i5,
', IWY=', i5 )
714 9997
FORMAT( / 1x, a3,
' -- Complex Expert Eigenvalue/vector',
715 $
' problem driver' )
717 9996
FORMAT(
'Input Example' )
719 9995
FORMAT(
' Matrix types: ', / )
721 9994
FORMAT(
' TYPE 1: Da is diagonal, Db is identity, ',
722 $ /
' A = Y^(-H) Da X^(-1), B = Y^(-H) Db X^(-1) ',
723 $ /
' YH and X are left and right eigenvectors. ', / )
725 9993
FORMAT(
' TYPE 2: Da is quasi-diagonal, Db is identity, ',
726 $ /
' A = Y^(-H) Da X^(-1), B = Y^(-H) Db X^(-1) ',
727 $ /
' YH and X are left and right eigenvectors. ', / )
729 9992
FORMAT( /
' Tests performed: ', / 4x,
730 $
' a is alpha, b is beta, l is a left eigenvector, ', / 4x,
731 $
' r is a right eigenvector and ', a,
' means ', a,
'.',
732 $ /
' 1 = max | ( b A - a B )', a,
' l | / const.',
733 $ /
' 2 = max | ( b A - a B ) r | / const.',
734 $ /
' 3 = max ( Sest/Stru, Stru/Sest ) ',
735 $
' over all eigenvalues', /
736 $
' 4 = max( DIFest/DIFtru, DIFtru/DIFest ) ',
737 $
' over the 1st and 5th eigenvectors', / )
739 9991
FORMAT(
' Type=', i2,
',',
' IWA=', i2,
', IWB=', i2,
', IWX=',
740 $ i2,
', IWY=', i2,
', result ', i2,
' is', 0p, f8.2 )
742 9990
FORMAT(
' Type=', i2,
',',
' IWA=', i2,
', IWB=', i2,
', IWX=',
743 $ i2,
', IWY=', i2,
', result ', i2,
' is', 1p, d10.3 )
745 9989
FORMAT(
' Input example #', i2,
', matrix order=', i4,
',',
746 $
' result ', i2,
' is', 0p, f8.2 )
748 9988
FORMAT(
' Input example #', i2,
', matrix order=', i4,
',',
749 $
' result ', i2,
' is', 1p, d10.3 )
751 9987
FORMAT(
' ZDRGVX: ', a,
' returned INFO=', i6,
'.', / 9x,
'N=',
752 $ i6,
', Input example #', i2,
')' )
754 9986
FORMAT(
' ZDRGVX: ', a,
' Eigenvectors from ', a,
' incorrectly ',
755 $
'normalized.', /
' Bits of error=', 0p, g10.3,
',', 9x,
756 $
'N=', i6,
', Input Example #', i2,
')' )
subroutine alasvm(TYPE, NOUT, NFAIL, NRUN, NERRS)
ALASVM
subroutine zdrgvx(NSIZE, THRESH, NIN, NOUT, A, LDA, B, AI, BI, ALPHA, BETA, VL, VR, ILO, IHI, LSCALE, RSCALE, S, DTRU, DIF, DIFTRU, WORK, LWORK, RWORK, IWORK, LIWORK, RESULT, BWORK, INFO)
ZDRGVX
subroutine zlatm6(TYPE, N, A, LDA, B, X, LDX, Y, LDY, ALPHA, BETA, WX, WY, S, DIF)
ZLATM6
subroutine zget52(LEFT, N, A, LDA, B, LDB, E, LDE, ALPHA, BETA, WORK, RWORK, RESULT)
ZGET52
subroutine zlacpy(UPLO, M, N, A, LDA, B, LDB)
ZLACPY copies all or part of one two-dimensional array to another.
subroutine zggevx(BALANC, JOBVL, JOBVR, SENSE, N, A, LDA, B, LDB, ALPHA, BETA, VL, LDVL, VR, LDVR, ILO, IHI, LSCALE, RSCALE, ABNRM, BBNRM, RCONDE, RCONDV, WORK, LWORK, RWORK, IWORK, BWORK, INFO)
ZGGEVX computes the eigenvalues and, optionally, the left and/or right eigenvectors for GE matrices ...
subroutine xerbla(SRNAME, INFO)
XERBLA