LAPACK 3.11.0
LAPACK: Linear Algebra PACKage
Loading...
Searching...
No Matches

◆ dchkgg()

subroutine dchkgg ( integer  NSIZES,
integer, dimension( * )  NN,
integer  NTYPES,
logical, dimension( * )  DOTYPE,
integer, dimension( 4 )  ISEED,
double precision  THRESH,
logical  TSTDIF,
double precision  THRSHN,
integer  NOUNIT,
double precision, dimension( lda, * )  A,
integer  LDA,
double precision, dimension( lda, * )  B,
double precision, dimension( lda, * )  H,
double precision, dimension( lda, * )  T,
double precision, dimension( lda, * )  S1,
double precision, dimension( lda, * )  S2,
double precision, dimension( lda, * )  P1,
double precision, dimension( lda, * )  P2,
double precision, dimension( ldu, * )  U,
integer  LDU,
double precision, dimension( ldu, * )  V,
double precision, dimension( ldu, * )  Q,
double precision, dimension( ldu, * )  Z,
double precision, dimension( * )  ALPHR1,
double precision, dimension( * )  ALPHI1,
double precision, dimension( * )  BETA1,
double precision, dimension( * )  ALPHR3,
double precision, dimension( * )  ALPHI3,
double precision, dimension( * )  BETA3,
double precision, dimension( ldu, * )  EVECTL,
double precision, dimension( ldu, * )  EVECTR,
double precision, dimension( * )  WORK,
integer  LWORK,
logical, dimension( * )  LLWORK,
double precision, dimension( 15 )  RESULT,
integer  INFO 
)

DCHKGG

Purpose:
 DCHKGG  checks the nonsymmetric generalized eigenvalue problem
 routines.
                                T          T        T
 DGGHRD factors A and B as U H V  and U T V , where   means
 transpose, H is hessenberg, T is triangular and U and V are
 orthogonal.
                                 T          T
 DHGEQZ factors H and T as  Q S Z  and Q P Z , where P is upper
 triangular, S is in generalized Schur form (block upper triangular,
 with 1x1 and 2x2 blocks on the diagonal, the 2x2 blocks
 corresponding to complex conjugate pairs of generalized
 eigenvalues), and Q and Z are orthogonal.  It also computes the
 generalized eigenvalues (alpha(1),beta(1)),...,(alpha(n),beta(n)),
 where alpha(j)=S(j,j) and beta(j)=P(j,j) -- thus,
 w(j) = alpha(j)/beta(j) is a root of the generalized eigenvalue
 problem

     det( A - w(j) B ) = 0

 and m(j) = beta(j)/alpha(j) is a root of the essentially equivalent
 problem

     det( m(j) A - B ) = 0

 DTGEVC computes the matrix L of left eigenvectors and the matrix R
 of right eigenvectors for the matrix pair ( S, P ).  In the
 description below,  l and r are left and right eigenvectors
 corresponding to the generalized eigenvalues (alpha,beta).

 When DCHKGG is called, a number of matrix "sizes" ("n's") and a
 number of matrix "types" are specified.  For each size ("n")
 and each type of matrix, one matrix will be generated and used
 to test the nonsymmetric eigenroutines.  For each matrix, 15
 tests will be performed.  The first twelve "test ratios" should be
 small -- O(1).  They will be compared with the threshold THRESH:

                  T
 (1)   | A - U H V  | / ( |A| n ulp )

                  T
 (2)   | B - U T V  | / ( |B| n ulp )

               T
 (3)   | I - UU  | / ( n ulp )

               T
 (4)   | I - VV  | / ( n ulp )

                  T
 (5)   | H - Q S Z  | / ( |H| n ulp )

                  T
 (6)   | T - Q P Z  | / ( |T| n ulp )

               T
 (7)   | I - QQ  | / ( n ulp )

               T
 (8)   | I - ZZ  | / ( n ulp )

 (9)   max over all left eigenvalue/-vector pairs (beta/alpha,l) of

    | l**H * (beta S - alpha P) | / ( ulp max( |beta S|, |alpha P| ) )

 (10)  max over all left eigenvalue/-vector pairs (beta/alpha,l') of
                           T
   | l'**H * (beta H - alpha T) | / ( ulp max( |beta H|, |alpha T| ) )

       where the eigenvectors l' are the result of passing Q to
       DTGEVC and back transforming (HOWMNY='B').

 (11)  max over all right eigenvalue/-vector pairs (beta/alpha,r) of

       | (beta S - alpha T) r | / ( ulp max( |beta S|, |alpha T| ) )

 (12)  max over all right eigenvalue/-vector pairs (beta/alpha,r') of

       | (beta H - alpha T) r' | / ( ulp max( |beta H|, |alpha T| ) )

       where the eigenvectors r' are the result of passing Z to
       DTGEVC and back transforming (HOWMNY='B').

 The last three test ratios will usually be small, but there is no
 mathematical requirement that they be so.  They are therefore
 compared with THRESH only if TSTDIF is .TRUE.

 (13)  | S(Q,Z computed) - S(Q,Z not computed) | / ( |S| ulp )

 (14)  | P(Q,Z computed) - P(Q,Z not computed) | / ( |P| ulp )

 (15)  max( |alpha(Q,Z computed) - alpha(Q,Z not computed)|/|S| ,
            |beta(Q,Z computed) - beta(Q,Z not computed)|/|P| ) / ulp

 In addition, the normalization of L and R are checked, and compared
 with the threshold THRSHN.

 Test Matrices
 ---- --------

 The sizes of the test matrices are specified by an array
 NN(1:NSIZES); the value of each element NN(j) specifies one size.
 The "types" are specified by a logical array DOTYPE( 1:NTYPES ); if
 DOTYPE(j) is .TRUE., then matrix type "j" will be generated.
 Currently, the list of possible types is:

 (1)  ( 0, 0 )         (a pair of zero matrices)

 (2)  ( I, 0 )         (an identity and a zero matrix)

 (3)  ( 0, I )         (an identity and a zero matrix)

 (4)  ( I, I )         (a pair of identity matrices)

         t   t
 (5)  ( J , J  )       (a pair of transposed Jordan blocks)

                                     t                ( I   0  )
 (6)  ( X, Y )         where  X = ( J   0  )  and Y = (      t )
                                  ( 0   I  )          ( 0   J  )
                       and I is a k x k identity and J a (k+1)x(k+1)
                       Jordan block; k=(N-1)/2

 (7)  ( D, I )         where D is diag( 0, 1,..., N-1 ) (a diagonal
                       matrix with those diagonal entries.)
 (8)  ( I, D )

 (9)  ( big*D, small*I ) where "big" is near overflow and small=1/big

 (10) ( small*D, big*I )

 (11) ( big*I, small*D )

 (12) ( small*I, big*D )

 (13) ( big*D, big*I )

 (14) ( small*D, small*I )

 (15) ( D1, D2 )        where D1 is diag( 0, 0, 1, ..., N-3, 0 ) and
                        D2 is diag( 0, N-3, N-4,..., 1, 0, 0 )
           t   t
 (16) U ( J , J ) V     where U and V are random orthogonal matrices.

 (17) U ( T1, T2 ) V    where T1 and T2 are upper triangular matrices
                        with random O(1) entries above the diagonal
                        and diagonal entries diag(T1) =
                        ( 0, 0, 1, ..., N-3, 0 ) and diag(T2) =
                        ( 0, N-3, N-4,..., 1, 0, 0 )

 (18) U ( T1, T2 ) V    diag(T1) = ( 0, 0, 1, 1, s, ..., s, 0 )
                        diag(T2) = ( 0, 1, 0, 1,..., 1, 0 )
                        s = machine precision.

 (19) U ( T1, T2 ) V    diag(T1)=( 0,0,1,1, 1-d, ..., 1-(N-5)*d=s, 0 )
                        diag(T2) = ( 0, 1, 0, 1, ..., 1, 0 )

                                                        N-5
 (20) U ( T1, T2 ) V    diag(T1)=( 0, 0, 1, 1, a, ..., a   =s, 0 )
                        diag(T2) = ( 0, 1, 0, 1, ..., 1, 0, 0 )

 (21) U ( T1, T2 ) V    diag(T1)=( 0, 0, 1, r1, r2, ..., r(N-4), 0 )
                        diag(T2) = ( 0, 1, 0, 1, ..., 1, 0, 0 )
                        where r1,..., r(N-4) are random.

 (22) U ( big*T1, small*T2 ) V    diag(T1) = ( 0, 0, 1, ..., N-3, 0 )
                                  diag(T2) = ( 0, 1, ..., 1, 0, 0 )

 (23) U ( small*T1, big*T2 ) V    diag(T1) = ( 0, 0, 1, ..., N-3, 0 )
                                  diag(T2) = ( 0, 1, ..., 1, 0, 0 )

 (24) U ( small*T1, small*T2 ) V  diag(T1) = ( 0, 0, 1, ..., N-3, 0 )
                                  diag(T2) = ( 0, 1, ..., 1, 0, 0 )

 (25) U ( big*T1, big*T2 ) V      diag(T1) = ( 0, 0, 1, ..., N-3, 0 )
                                  diag(T2) = ( 0, 1, ..., 1, 0, 0 )

 (26) U ( T1, T2 ) V     where T1 and T2 are random upper-triangular
                         matrices.
Parameters
[in]NSIZES
          NSIZES is INTEGER
          The number of sizes of matrices to use.  If it is zero,
          DCHKGG does nothing.  It must be at least zero.
[in]NN
          NN is INTEGER array, dimension (NSIZES)
          An array containing the sizes to be used for the matrices.
          Zero values will be skipped.  The values must be at least
          zero.
[in]NTYPES
          NTYPES is INTEGER
          The number of elements in DOTYPE.   If it is zero, DCHKGG
          does nothing.  It must be at least zero.  If it is MAXTYP+1
          and NSIZES is 1, then an additional type, MAXTYP+1 is
          defined, which is to use whatever matrix is in A.  This
          is only useful if DOTYPE(1:MAXTYP) is .FALSE. and
          DOTYPE(MAXTYP+1) is .TRUE. .
[in]DOTYPE
          DOTYPE is LOGICAL array, dimension (NTYPES)
          If DOTYPE(j) is .TRUE., then for each size in NN a
          matrix of that size and of type j will be generated.
          If NTYPES is smaller than the maximum number of types
          defined (PARAMETER MAXTYP), then types NTYPES+1 through
          MAXTYP will not be generated.  If NTYPES is larger
          than MAXTYP, DOTYPE(MAXTYP+1) through DOTYPE(NTYPES)
          will be ignored.
[in,out]ISEED
          ISEED is INTEGER array, dimension (4)
          On entry ISEED specifies the seed of the random number
          generator. The array elements should be between 0 and 4095;
          if not they will be reduced mod 4096.  Also, ISEED(4) must
          be odd.  The random number generator uses a linear
          congruential sequence limited to small integers, and so
          should produce machine independent random numbers. The
          values of ISEED are changed on exit, and can be used in the
          next call to DCHKGG to continue the same random number
          sequence.
[in]THRESH
          THRESH is DOUBLE PRECISION
          A test will count as "failed" if the "error", computed as
          described above, exceeds THRESH.  Note that the error is
          scaled to be O(1), so THRESH should be a reasonably small
          multiple of 1, e.g., 10 or 100.  In particular, it should
          not depend on the precision (single vs. double) or the size
          of the matrix.  It must be at least zero.
[in]TSTDIF
          TSTDIF is LOGICAL
          Specifies whether test ratios 13-15 will be computed and
          compared with THRESH.
          = .FALSE.: Only test ratios 1-12 will be computed and tested.
                     Ratios 13-15 will be set to zero.
          = .TRUE.:  All the test ratios 1-15 will be computed and
                     tested.
[in]THRSHN
          THRSHN is DOUBLE PRECISION
          Threshold for reporting eigenvector normalization error.
          If the normalization of any eigenvector differs from 1 by
          more than THRSHN*ulp, then a special error message will be
          printed.  (This is handled separately from the other tests,
          since only a compiler or programming error should cause an
          error message, at least if THRSHN is at least 5--10.)
[in]NOUNIT
          NOUNIT is INTEGER
          The FORTRAN unit number for printing out error messages
          (e.g., if a routine returns IINFO not equal to 0.)
[in,out]A
          A is DOUBLE PRECISION array, dimension
                            (LDA, max(NN))
          Used to hold the original A matrix.  Used as input only
          if NTYPES=MAXTYP+1, DOTYPE(1:MAXTYP)=.FALSE., and
          DOTYPE(MAXTYP+1)=.TRUE.
[in]LDA
          LDA is INTEGER
          The leading dimension of A, B, H, T, S1, P1, S2, and P2.
          It must be at least 1 and at least max( NN ).
[in,out]B
          B is DOUBLE PRECISION array, dimension
                            (LDA, max(NN))
          Used to hold the original B matrix.  Used as input only
          if NTYPES=MAXTYP+1, DOTYPE(1:MAXTYP)=.FALSE., and
          DOTYPE(MAXTYP+1)=.TRUE.
[out]H
          H is DOUBLE PRECISION array, dimension (LDA, max(NN))
          The upper Hessenberg matrix computed from A by DGGHRD.
[out]T
          T is DOUBLE PRECISION array, dimension (LDA, max(NN))
          The upper triangular matrix computed from B by DGGHRD.
[out]S1
          S1 is DOUBLE PRECISION array, dimension (LDA, max(NN))
          The Schur (block upper triangular) matrix computed from H by
          DHGEQZ when Q and Z are also computed.
[out]S2
          S2 is DOUBLE PRECISION array, dimension (LDA, max(NN))
          The Schur (block upper triangular) matrix computed from H by
          DHGEQZ when Q and Z are not computed.
[out]P1
          P1 is DOUBLE PRECISION array, dimension (LDA, max(NN))
          The upper triangular matrix computed from T by DHGEQZ
          when Q and Z are also computed.
[out]P2
          P2 is DOUBLE PRECISION array, dimension (LDA, max(NN))
          The upper triangular matrix computed from T by DHGEQZ
          when Q and Z are not computed.
[out]U
          U is DOUBLE PRECISION array, dimension (LDU, max(NN))
          The (left) orthogonal matrix computed by DGGHRD.
[in]LDU
          LDU is INTEGER
          The leading dimension of U, V, Q, Z, EVECTL, and EVEZTR.  It
          must be at least 1 and at least max( NN ).
[out]V
          V is DOUBLE PRECISION array, dimension (LDU, max(NN))
          The (right) orthogonal matrix computed by DGGHRD.
[out]Q
          Q is DOUBLE PRECISION array, dimension (LDU, max(NN))
          The (left) orthogonal matrix computed by DHGEQZ.
[out]Z
          Z is DOUBLE PRECISION array, dimension (LDU, max(NN))
          The (left) orthogonal matrix computed by DHGEQZ.
[out]ALPHR1
          ALPHR1 is DOUBLE PRECISION array, dimension (max(NN))
[out]ALPHI1
          ALPHI1 is DOUBLE PRECISION array, dimension (max(NN))
[out]BETA1
          BETA1 is DOUBLE PRECISION array, dimension (max(NN))

          The generalized eigenvalues of (A,B) computed by DHGEQZ
          when Q, Z, and the full Schur matrices are computed.
          On exit, ( ALPHR1(k)+ALPHI1(k)*i ) / BETA1(k) is the k-th
          generalized eigenvalue of the matrices in A and B.
[out]ALPHR3
          ALPHR3 is DOUBLE PRECISION array, dimension (max(NN))
[out]ALPHI3
          ALPHI3 is DOUBLE PRECISION array, dimension (max(NN))
[out]BETA3
          BETA3 is DOUBLE PRECISION array, dimension (max(NN))
[out]EVECTL
          EVECTL is DOUBLE PRECISION array, dimension (LDU, max(NN))
          The (block lower triangular) left eigenvector matrix for
          the matrices in S1 and P1.  (See DTGEVC for the format.)
[out]EVECTR
          EVECTR is DOUBLE PRECISION array, dimension (LDU, max(NN))
          The (block upper triangular) right eigenvector matrix for
          the matrices in S1 and P1.  (See DTGEVC for the format.)
[out]WORK
          WORK is DOUBLE PRECISION array, dimension (LWORK)
[in]LWORK
          LWORK is INTEGER
          The number of entries in WORK.  This must be at least
          max( 2 * N**2, 6*N, 1 ), for all N=NN(j).
[out]LLWORK
          LLWORK is LOGICAL array, dimension (max(NN))
[out]RESULT
          RESULT is DOUBLE PRECISION array, dimension (15)
          The values computed by the tests described above.
          The values are currently limited to 1/ulp, to avoid
          overflow.
[out]INFO
          INFO is INTEGER
          = 0:  successful exit
          < 0:  if INFO = -i, the i-th argument had an illegal value
          > 0:  A routine returned an error code.  INFO is the
                absolute value of the INFO value returned.
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.

Definition at line 506 of file dchkgg.f.

511*
512* -- LAPACK test routine --
513* -- LAPACK is a software package provided by Univ. of Tennessee, --
514* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
515*
516* .. Scalar Arguments ..
517 LOGICAL TSTDIF
518 INTEGER INFO, LDA, LDU, LWORK, NOUNIT, NSIZES, NTYPES
519 DOUBLE PRECISION THRESH, THRSHN
520* ..
521* .. Array Arguments ..
522 LOGICAL DOTYPE( * ), LLWORK( * )
523 INTEGER ISEED( 4 ), NN( * )
524 DOUBLE PRECISION A( LDA, * ), ALPHI1( * ), ALPHI3( * ),
525 $ ALPHR1( * ), ALPHR3( * ), B( LDA, * ),
526 $ BETA1( * ), BETA3( * ), EVECTL( LDU, * ),
527 $ EVECTR( LDU, * ), H( LDA, * ), P1( LDA, * ),
528 $ P2( LDA, * ), Q( LDU, * ), RESULT( 15 ),
529 $ S1( LDA, * ), S2( LDA, * ), T( LDA, * ),
530 $ U( LDU, * ), V( LDU, * ), WORK( * ),
531 $ Z( LDU, * )
532* ..
533*
534* =====================================================================
535*
536* .. Parameters ..
537 DOUBLE PRECISION ZERO, ONE
538 parameter( zero = 0.0d0, one = 1.0d0 )
539 INTEGER MAXTYP
540 parameter( maxtyp = 26 )
541* ..
542* .. Local Scalars ..
543 LOGICAL BADNN
544 INTEGER I1, IADD, IINFO, IN, J, JC, JR, JSIZE, JTYPE,
545 $ LWKOPT, MTYPES, N, N1, NERRS, NMATS, NMAX,
546 $ NTEST, NTESTT
547 DOUBLE PRECISION ANORM, BNORM, SAFMAX, SAFMIN, TEMP1, TEMP2,
548 $ ULP, ULPINV
549* ..
550* .. Local Arrays ..
551 INTEGER IASIGN( MAXTYP ), IBSIGN( MAXTYP ),
552 $ IOLDSD( 4 ), KADD( 6 ), KAMAGN( MAXTYP ),
553 $ KATYPE( MAXTYP ), KAZERO( MAXTYP ),
554 $ KBMAGN( MAXTYP ), KBTYPE( MAXTYP ),
555 $ KBZERO( MAXTYP ), KCLASS( MAXTYP ),
556 $ KTRIAN( MAXTYP ), KZ1( 6 ), KZ2( 6 )
557 DOUBLE PRECISION DUMMA( 4 ), RMAGN( 0: 3 )
558* ..
559* .. External Functions ..
560 DOUBLE PRECISION DLAMCH, DLANGE, DLARND
561 EXTERNAL dlamch, dlange, dlarnd
562* ..
563* .. External Subroutines ..
564 EXTERNAL dgeqr2, dget51, dget52, dgghrd, dhgeqz, dlabad,
566 $ dtgevc, xerbla
567* ..
568* .. Intrinsic Functions ..
569 INTRINSIC abs, dble, max, min, sign
570* ..
571* .. Data statements ..
572 DATA kclass / 15*1, 10*2, 1*3 /
573 DATA kz1 / 0, 1, 2, 1, 3, 3 /
574 DATA kz2 / 0, 0, 1, 2, 1, 1 /
575 DATA kadd / 0, 0, 0, 0, 3, 2 /
576 DATA katype / 0, 1, 0, 1, 2, 3, 4, 1, 4, 4, 1, 1, 4,
577 $ 4, 4, 2, 4, 5, 8, 7, 9, 4*4, 0 /
578 DATA kbtype / 0, 0, 1, 1, 2, -3, 1, 4, 1, 1, 4, 4,
579 $ 1, 1, -4, 2, -4, 8*8, 0 /
580 DATA kazero / 6*1, 2, 1, 2*2, 2*1, 2*2, 3, 1, 3,
581 $ 4*5, 4*3, 1 /
582 DATA kbzero / 6*1, 1, 2, 2*1, 2*2, 2*1, 4, 1, 4,
583 $ 4*6, 4*4, 1 /
584 DATA kamagn / 8*1, 2, 3, 2, 3, 2, 3, 7*1, 2, 3, 3,
585 $ 2, 1 /
586 DATA kbmagn / 8*1, 3, 2, 3, 2, 2, 3, 7*1, 3, 2, 3,
587 $ 2, 1 /
588 DATA ktrian / 16*0, 10*1 /
589 DATA iasign / 6*0, 2, 0, 2*2, 2*0, 3*2, 0, 2, 3*0,
590 $ 5*2, 0 /
591 DATA ibsign / 7*0, 2, 2*0, 2*2, 2*0, 2, 0, 2, 9*0 /
592* ..
593* .. Executable Statements ..
594*
595* Check for errors
596*
597 info = 0
598*
599 badnn = .false.
600 nmax = 1
601 DO 10 j = 1, nsizes
602 nmax = max( nmax, nn( j ) )
603 IF( nn( j ).LT.0 )
604 $ badnn = .true.
605 10 CONTINUE
606*
607* Maximum blocksize and shift -- we assume that blocksize and number
608* of shifts are monotone increasing functions of N.
609*
610 lwkopt = max( 6*nmax, 2*nmax*nmax, 1 )
611*
612* Check for errors
613*
614 IF( nsizes.LT.0 ) THEN
615 info = -1
616 ELSE IF( badnn ) THEN
617 info = -2
618 ELSE IF( ntypes.LT.0 ) THEN
619 info = -3
620 ELSE IF( thresh.LT.zero ) THEN
621 info = -6
622 ELSE IF( lda.LE.1 .OR. lda.LT.nmax ) THEN
623 info = -10
624 ELSE IF( ldu.LE.1 .OR. ldu.LT.nmax ) THEN
625 info = -19
626 ELSE IF( lwkopt.GT.lwork ) THEN
627 info = -30
628 END IF
629*
630 IF( info.NE.0 ) THEN
631 CALL xerbla( 'DCHKGG', -info )
632 RETURN
633 END IF
634*
635* Quick return if possible
636*
637 IF( nsizes.EQ.0 .OR. ntypes.EQ.0 )
638 $ RETURN
639*
640 safmin = dlamch( 'Safe minimum' )
641 ulp = dlamch( 'Epsilon' )*dlamch( 'Base' )
642 safmin = safmin / ulp
643 safmax = one / safmin
644 CALL dlabad( safmin, safmax )
645 ulpinv = one / ulp
646*
647* The values RMAGN(2:3) depend on N, see below.
648*
649 rmagn( 0 ) = zero
650 rmagn( 1 ) = one
651*
652* Loop over sizes, types
653*
654 ntestt = 0
655 nerrs = 0
656 nmats = 0
657*
658 DO 240 jsize = 1, nsizes
659 n = nn( jsize )
660 n1 = max( 1, n )
661 rmagn( 2 ) = safmax*ulp / dble( n1 )
662 rmagn( 3 ) = safmin*ulpinv*n1
663*
664 IF( nsizes.NE.1 ) THEN
665 mtypes = min( maxtyp, ntypes )
666 ELSE
667 mtypes = min( maxtyp+1, ntypes )
668 END IF
669*
670 DO 230 jtype = 1, mtypes
671 IF( .NOT.dotype( jtype ) )
672 $ GO TO 230
673 nmats = nmats + 1
674 ntest = 0
675*
676* Save ISEED in case of an error.
677*
678 DO 20 j = 1, 4
679 ioldsd( j ) = iseed( j )
680 20 CONTINUE
681*
682* Initialize RESULT
683*
684 DO 30 j = 1, 15
685 result( j ) = zero
686 30 CONTINUE
687*
688* Compute A and B
689*
690* Description of control parameters:
691*
692* KZLASS: =1 means w/o rotation, =2 means w/ rotation,
693* =3 means random.
694* KATYPE: the "type" to be passed to DLATM4 for computing A.
695* KAZERO: the pattern of zeros on the diagonal for A:
696* =1: ( xxx ), =2: (0, xxx ) =3: ( 0, 0, xxx, 0 ),
697* =4: ( 0, xxx, 0, 0 ), =5: ( 0, 0, 1, xxx, 0 ),
698* =6: ( 0, 1, 0, xxx, 0 ). (xxx means a string of
699* non-zero entries.)
700* KAMAGN: the magnitude of the matrix: =0: zero, =1: O(1),
701* =2: large, =3: small.
702* IASIGN: 1 if the diagonal elements of A are to be
703* multiplied by a random magnitude 1 number, =2 if
704* randomly chosen diagonal blocks are to be rotated
705* to form 2x2 blocks.
706* KBTYPE, KBZERO, KBMAGN, IBSIGN: the same, but for B.
707* KTRIAN: =0: don't fill in the upper triangle, =1: do.
708* KZ1, KZ2, KADD: used to implement KAZERO and KBZERO.
709* RMAGN: used to implement KAMAGN and KBMAGN.
710*
711 IF( mtypes.GT.maxtyp )
712 $ GO TO 110
713 iinfo = 0
714 IF( kclass( jtype ).LT.3 ) THEN
715*
716* Generate A (w/o rotation)
717*
718 IF( abs( katype( jtype ) ).EQ.3 ) THEN
719 in = 2*( ( n-1 ) / 2 ) + 1
720 IF( in.NE.n )
721 $ CALL dlaset( 'Full', n, n, zero, zero, a, lda )
722 ELSE
723 in = n
724 END IF
725 CALL dlatm4( katype( jtype ), in, kz1( kazero( jtype ) ),
726 $ kz2( kazero( jtype ) ), iasign( jtype ),
727 $ rmagn( kamagn( jtype ) ), ulp,
728 $ rmagn( ktrian( jtype )*kamagn( jtype ) ), 2,
729 $ iseed, a, lda )
730 iadd = kadd( kazero( jtype ) )
731 IF( iadd.GT.0 .AND. iadd.LE.n )
732 $ a( iadd, iadd ) = rmagn( kamagn( jtype ) )
733*
734* Generate B (w/o rotation)
735*
736 IF( abs( kbtype( jtype ) ).EQ.3 ) THEN
737 in = 2*( ( n-1 ) / 2 ) + 1
738 IF( in.NE.n )
739 $ CALL dlaset( 'Full', n, n, zero, zero, b, lda )
740 ELSE
741 in = n
742 END IF
743 CALL dlatm4( kbtype( jtype ), in, kz1( kbzero( jtype ) ),
744 $ kz2( kbzero( jtype ) ), ibsign( jtype ),
745 $ rmagn( kbmagn( jtype ) ), one,
746 $ rmagn( ktrian( jtype )*kbmagn( jtype ) ), 2,
747 $ iseed, b, lda )
748 iadd = kadd( kbzero( jtype ) )
749 IF( iadd.NE.0 .AND. iadd.LE.n )
750 $ b( iadd, iadd ) = rmagn( kbmagn( jtype ) )
751*
752 IF( kclass( jtype ).EQ.2 .AND. n.GT.0 ) THEN
753*
754* Include rotations
755*
756* Generate U, V as Householder transformations times
757* a diagonal matrix.
758*
759 DO 50 jc = 1, n - 1
760 DO 40 jr = jc, n
761 u( jr, jc ) = dlarnd( 3, iseed )
762 v( jr, jc ) = dlarnd( 3, iseed )
763 40 CONTINUE
764 CALL dlarfg( n+1-jc, u( jc, jc ), u( jc+1, jc ), 1,
765 $ work( jc ) )
766 work( 2*n+jc ) = sign( one, u( jc, jc ) )
767 u( jc, jc ) = one
768 CALL dlarfg( n+1-jc, v( jc, jc ), v( jc+1, jc ), 1,
769 $ work( n+jc ) )
770 work( 3*n+jc ) = sign( one, v( jc, jc ) )
771 v( jc, jc ) = one
772 50 CONTINUE
773 u( n, n ) = one
774 work( n ) = zero
775 work( 3*n ) = sign( one, dlarnd( 2, iseed ) )
776 v( n, n ) = one
777 work( 2*n ) = zero
778 work( 4*n ) = sign( one, dlarnd( 2, iseed ) )
779*
780* Apply the diagonal matrices
781*
782 DO 70 jc = 1, n
783 DO 60 jr = 1, n
784 a( jr, jc ) = work( 2*n+jr )*work( 3*n+jc )*
785 $ a( jr, jc )
786 b( jr, jc ) = work( 2*n+jr )*work( 3*n+jc )*
787 $ b( jr, jc )
788 60 CONTINUE
789 70 CONTINUE
790 CALL dorm2r( 'L', 'N', n, n, n-1, u, ldu, work, a,
791 $ lda, work( 2*n+1 ), iinfo )
792 IF( iinfo.NE.0 )
793 $ GO TO 100
794 CALL dorm2r( 'R', 'T', n, n, n-1, v, ldu, work( n+1 ),
795 $ a, lda, work( 2*n+1 ), iinfo )
796 IF( iinfo.NE.0 )
797 $ GO TO 100
798 CALL dorm2r( 'L', 'N', n, n, n-1, u, ldu, work, b,
799 $ lda, work( 2*n+1 ), iinfo )
800 IF( iinfo.NE.0 )
801 $ GO TO 100
802 CALL dorm2r( 'R', 'T', n, n, n-1, v, ldu, work( n+1 ),
803 $ b, lda, work( 2*n+1 ), iinfo )
804 IF( iinfo.NE.0 )
805 $ GO TO 100
806 END IF
807 ELSE
808*
809* Random matrices
810*
811 DO 90 jc = 1, n
812 DO 80 jr = 1, n
813 a( jr, jc ) = rmagn( kamagn( jtype ) )*
814 $ dlarnd( 2, iseed )
815 b( jr, jc ) = rmagn( kbmagn( jtype ) )*
816 $ dlarnd( 2, iseed )
817 80 CONTINUE
818 90 CONTINUE
819 END IF
820*
821 anorm = dlange( '1', n, n, a, lda, work )
822 bnorm = dlange( '1', n, n, b, lda, work )
823*
824 100 CONTINUE
825*
826 IF( iinfo.NE.0 ) THEN
827 WRITE( nounit, fmt = 9999 )'Generator', iinfo, n, jtype,
828 $ ioldsd
829 info = abs( iinfo )
830 RETURN
831 END IF
832*
833 110 CONTINUE
834*
835* Call DGEQR2, DORM2R, and DGGHRD to compute H, T, U, and V
836*
837 CALL dlacpy( ' ', n, n, a, lda, h, lda )
838 CALL dlacpy( ' ', n, n, b, lda, t, lda )
839 ntest = 1
840 result( 1 ) = ulpinv
841*
842 CALL dgeqr2( n, n, t, lda, work, work( n+1 ), iinfo )
843 IF( iinfo.NE.0 ) THEN
844 WRITE( nounit, fmt = 9999 )'DGEQR2', iinfo, n, jtype,
845 $ ioldsd
846 info = abs( iinfo )
847 GO TO 210
848 END IF
849*
850 CALL dorm2r( 'L', 'T', n, n, n, t, lda, work, h, lda,
851 $ work( n+1 ), iinfo )
852 IF( iinfo.NE.0 ) THEN
853 WRITE( nounit, fmt = 9999 )'DORM2R', iinfo, n, jtype,
854 $ ioldsd
855 info = abs( iinfo )
856 GO TO 210
857 END IF
858*
859 CALL dlaset( 'Full', n, n, zero, one, u, ldu )
860 CALL dorm2r( 'R', 'N', n, n, n, t, lda, work, u, ldu,
861 $ work( n+1 ), iinfo )
862 IF( iinfo.NE.0 ) THEN
863 WRITE( nounit, fmt = 9999 )'DORM2R', iinfo, n, jtype,
864 $ ioldsd
865 info = abs( iinfo )
866 GO TO 210
867 END IF
868*
869 CALL dgghrd( 'V', 'I', n, 1, n, h, lda, t, lda, u, ldu, v,
870 $ ldu, iinfo )
871 IF( iinfo.NE.0 ) THEN
872 WRITE( nounit, fmt = 9999 )'DGGHRD', iinfo, n, jtype,
873 $ ioldsd
874 info = abs( iinfo )
875 GO TO 210
876 END IF
877 ntest = 4
878*
879* Do tests 1--4
880*
881 CALL dget51( 1, n, a, lda, h, lda, u, ldu, v, ldu, work,
882 $ result( 1 ) )
883 CALL dget51( 1, n, b, lda, t, lda, u, ldu, v, ldu, work,
884 $ result( 2 ) )
885 CALL dget51( 3, n, b, lda, t, lda, u, ldu, u, ldu, work,
886 $ result( 3 ) )
887 CALL dget51( 3, n, b, lda, t, lda, v, ldu, v, ldu, work,
888 $ result( 4 ) )
889*
890* Call DHGEQZ to compute S1, P1, S2, P2, Q, and Z, do tests.
891*
892* Compute T1 and UZ
893*
894* Eigenvalues only
895*
896 CALL dlacpy( ' ', n, n, h, lda, s2, lda )
897 CALL dlacpy( ' ', n, n, t, lda, p2, lda )
898 ntest = 5
899 result( 5 ) = ulpinv
900*
901 CALL dhgeqz( 'E', 'N', 'N', n, 1, n, s2, lda, p2, lda,
902 $ alphr3, alphi3, beta3, q, ldu, z, ldu, work,
903 $ lwork, iinfo )
904 IF( iinfo.NE.0 ) THEN
905 WRITE( nounit, fmt = 9999 )'DHGEQZ(E)', iinfo, n, jtype,
906 $ ioldsd
907 info = abs( iinfo )
908 GO TO 210
909 END IF
910*
911* Eigenvalues and Full Schur Form
912*
913 CALL dlacpy( ' ', n, n, h, lda, s2, lda )
914 CALL dlacpy( ' ', n, n, t, lda, p2, lda )
915*
916 CALL dhgeqz( 'S', 'N', 'N', n, 1, n, s2, lda, p2, lda,
917 $ alphr1, alphi1, beta1, q, ldu, z, ldu, work,
918 $ lwork, iinfo )
919 IF( iinfo.NE.0 ) THEN
920 WRITE( nounit, fmt = 9999 )'DHGEQZ(S)', iinfo, n, jtype,
921 $ ioldsd
922 info = abs( iinfo )
923 GO TO 210
924 END IF
925*
926* Eigenvalues, Schur Form, and Schur Vectors
927*
928 CALL dlacpy( ' ', n, n, h, lda, s1, lda )
929 CALL dlacpy( ' ', n, n, t, lda, p1, lda )
930*
931 CALL dhgeqz( 'S', 'I', 'I', n, 1, n, s1, lda, p1, lda,
932 $ alphr1, alphi1, beta1, q, ldu, z, ldu, work,
933 $ lwork, iinfo )
934 IF( iinfo.NE.0 ) THEN
935 WRITE( nounit, fmt = 9999 )'DHGEQZ(V)', iinfo, n, jtype,
936 $ ioldsd
937 info = abs( iinfo )
938 GO TO 210
939 END IF
940*
941 ntest = 8
942*
943* Do Tests 5--8
944*
945 CALL dget51( 1, n, h, lda, s1, lda, q, ldu, z, ldu, work,
946 $ result( 5 ) )
947 CALL dget51( 1, n, t, lda, p1, lda, q, ldu, z, ldu, work,
948 $ result( 6 ) )
949 CALL dget51( 3, n, t, lda, p1, lda, q, ldu, q, ldu, work,
950 $ result( 7 ) )
951 CALL dget51( 3, n, t, lda, p1, lda, z, ldu, z, ldu, work,
952 $ result( 8 ) )
953*
954* Compute the Left and Right Eigenvectors of (S1,P1)
955*
956* 9: Compute the left eigenvector Matrix without
957* back transforming:
958*
959 ntest = 9
960 result( 9 ) = ulpinv
961*
962* To test "SELECT" option, compute half of the eigenvectors
963* in one call, and half in another
964*
965 i1 = n / 2
966 DO 120 j = 1, i1
967 llwork( j ) = .true.
968 120 CONTINUE
969 DO 130 j = i1 + 1, n
970 llwork( j ) = .false.
971 130 CONTINUE
972*
973 CALL dtgevc( 'L', 'S', llwork, n, s1, lda, p1, lda, evectl,
974 $ ldu, dumma, ldu, n, in, work, iinfo )
975 IF( iinfo.NE.0 ) THEN
976 WRITE( nounit, fmt = 9999 )'DTGEVC(L,S1)', iinfo, n,
977 $ jtype, ioldsd
978 info = abs( iinfo )
979 GO TO 210
980 END IF
981*
982 i1 = in
983 DO 140 j = 1, i1
984 llwork( j ) = .false.
985 140 CONTINUE
986 DO 150 j = i1 + 1, n
987 llwork( j ) = .true.
988 150 CONTINUE
989*
990 CALL dtgevc( 'L', 'S', llwork, n, s1, lda, p1, lda,
991 $ evectl( 1, i1+1 ), ldu, dumma, ldu, n, in,
992 $ work, iinfo )
993 IF( iinfo.NE.0 ) THEN
994 WRITE( nounit, fmt = 9999 )'DTGEVC(L,S2)', iinfo, n,
995 $ jtype, ioldsd
996 info = abs( iinfo )
997 GO TO 210
998 END IF
999*
1000 CALL dget52( .true., n, s1, lda, p1, lda, evectl, ldu,
1001 $ alphr1, alphi1, beta1, work, dumma( 1 ) )
1002 result( 9 ) = dumma( 1 )
1003 IF( dumma( 2 ).GT.thrshn ) THEN
1004 WRITE( nounit, fmt = 9998 )'Left', 'DTGEVC(HOWMNY=S)',
1005 $ dumma( 2 ), n, jtype, ioldsd
1006 END IF
1007*
1008* 10: Compute the left eigenvector Matrix with
1009* back transforming:
1010*
1011 ntest = 10
1012 result( 10 ) = ulpinv
1013 CALL dlacpy( 'F', n, n, q, ldu, evectl, ldu )
1014 CALL dtgevc( 'L', 'B', llwork, n, s1, lda, p1, lda, evectl,
1015 $ ldu, dumma, ldu, n, in, work, iinfo )
1016 IF( iinfo.NE.0 ) THEN
1017 WRITE( nounit, fmt = 9999 )'DTGEVC(L,B)', iinfo, n,
1018 $ jtype, ioldsd
1019 info = abs( iinfo )
1020 GO TO 210
1021 END IF
1022*
1023 CALL dget52( .true., n, h, lda, t, lda, evectl, ldu, alphr1,
1024 $ alphi1, beta1, work, dumma( 1 ) )
1025 result( 10 ) = dumma( 1 )
1026 IF( dumma( 2 ).GT.thrshn ) THEN
1027 WRITE( nounit, fmt = 9998 )'Left', 'DTGEVC(HOWMNY=B)',
1028 $ dumma( 2 ), n, jtype, ioldsd
1029 END IF
1030*
1031* 11: Compute the right eigenvector Matrix without
1032* back transforming:
1033*
1034 ntest = 11
1035 result( 11 ) = ulpinv
1036*
1037* To test "SELECT" option, compute half of the eigenvectors
1038* in one call, and half in another
1039*
1040 i1 = n / 2
1041 DO 160 j = 1, i1
1042 llwork( j ) = .true.
1043 160 CONTINUE
1044 DO 170 j = i1 + 1, n
1045 llwork( j ) = .false.
1046 170 CONTINUE
1047*
1048 CALL dtgevc( 'R', 'S', llwork, n, s1, lda, p1, lda, dumma,
1049 $ ldu, evectr, ldu, n, in, work, iinfo )
1050 IF( iinfo.NE.0 ) THEN
1051 WRITE( nounit, fmt = 9999 )'DTGEVC(R,S1)', iinfo, n,
1052 $ jtype, ioldsd
1053 info = abs( iinfo )
1054 GO TO 210
1055 END IF
1056*
1057 i1 = in
1058 DO 180 j = 1, i1
1059 llwork( j ) = .false.
1060 180 CONTINUE
1061 DO 190 j = i1 + 1, n
1062 llwork( j ) = .true.
1063 190 CONTINUE
1064*
1065 CALL dtgevc( 'R', 'S', llwork, n, s1, lda, p1, lda, dumma,
1066 $ ldu, evectr( 1, i1+1 ), ldu, n, in, work,
1067 $ iinfo )
1068 IF( iinfo.NE.0 ) THEN
1069 WRITE( nounit, fmt = 9999 )'DTGEVC(R,S2)', iinfo, n,
1070 $ jtype, ioldsd
1071 info = abs( iinfo )
1072 GO TO 210
1073 END IF
1074*
1075 CALL dget52( .false., n, s1, lda, p1, lda, evectr, ldu,
1076 $ alphr1, alphi1, beta1, work, dumma( 1 ) )
1077 result( 11 ) = dumma( 1 )
1078 IF( dumma( 2 ).GT.thresh ) THEN
1079 WRITE( nounit, fmt = 9998 )'Right', 'DTGEVC(HOWMNY=S)',
1080 $ dumma( 2 ), n, jtype, ioldsd
1081 END IF
1082*
1083* 12: Compute the right eigenvector Matrix with
1084* back transforming:
1085*
1086 ntest = 12
1087 result( 12 ) = ulpinv
1088 CALL dlacpy( 'F', n, n, z, ldu, evectr, ldu )
1089 CALL dtgevc( 'R', 'B', llwork, n, s1, lda, p1, lda, dumma,
1090 $ ldu, evectr, ldu, n, in, work, iinfo )
1091 IF( iinfo.NE.0 ) THEN
1092 WRITE( nounit, fmt = 9999 )'DTGEVC(R,B)', iinfo, n,
1093 $ jtype, ioldsd
1094 info = abs( iinfo )
1095 GO TO 210
1096 END IF
1097*
1098 CALL dget52( .false., n, h, lda, t, lda, evectr, ldu,
1099 $ alphr1, alphi1, beta1, work, dumma( 1 ) )
1100 result( 12 ) = dumma( 1 )
1101 IF( dumma( 2 ).GT.thresh ) THEN
1102 WRITE( nounit, fmt = 9998 )'Right', 'DTGEVC(HOWMNY=B)',
1103 $ dumma( 2 ), n, jtype, ioldsd
1104 END IF
1105*
1106* Tests 13--15 are done only on request
1107*
1108 IF( tstdif ) THEN
1109*
1110* Do Tests 13--14
1111*
1112 CALL dget51( 2, n, s1, lda, s2, lda, q, ldu, z, ldu,
1113 $ work, result( 13 ) )
1114 CALL dget51( 2, n, p1, lda, p2, lda, q, ldu, z, ldu,
1115 $ work, result( 14 ) )
1116*
1117* Do Test 15
1118*
1119 temp1 = zero
1120 temp2 = zero
1121 DO 200 j = 1, n
1122 temp1 = max( temp1, abs( alphr1( j )-alphr3( j ) )+
1123 $ abs( alphi1( j )-alphi3( j ) ) )
1124 temp2 = max( temp2, abs( beta1( j )-beta3( j ) ) )
1125 200 CONTINUE
1126*
1127 temp1 = temp1 / max( safmin, ulp*max( temp1, anorm ) )
1128 temp2 = temp2 / max( safmin, ulp*max( temp2, bnorm ) )
1129 result( 15 ) = max( temp1, temp2 )
1130 ntest = 15
1131 ELSE
1132 result( 13 ) = zero
1133 result( 14 ) = zero
1134 result( 15 ) = zero
1135 ntest = 12
1136 END IF
1137*
1138* End of Loop -- Check for RESULT(j) > THRESH
1139*
1140 210 CONTINUE
1141*
1142 ntestt = ntestt + ntest
1143*
1144* Print out tests which fail.
1145*
1146 DO 220 jr = 1, ntest
1147 IF( result( jr ).GE.thresh ) THEN
1148*
1149* If this is the first test to fail,
1150* print a header to the data file.
1151*
1152 IF( nerrs.EQ.0 ) THEN
1153 WRITE( nounit, fmt = 9997 )'DGG'
1154*
1155* Matrix types
1156*
1157 WRITE( nounit, fmt = 9996 )
1158 WRITE( nounit, fmt = 9995 )
1159 WRITE( nounit, fmt = 9994 )'Orthogonal'
1160*
1161* Tests performed
1162*
1163 WRITE( nounit, fmt = 9993 )'orthogonal', '''',
1164 $ 'transpose', ( '''', j = 1, 10 )
1165*
1166 END IF
1167 nerrs = nerrs + 1
1168 IF( result( jr ).LT.10000.0d0 ) THEN
1169 WRITE( nounit, fmt = 9992 )n, jtype, ioldsd, jr,
1170 $ result( jr )
1171 ELSE
1172 WRITE( nounit, fmt = 9991 )n, jtype, ioldsd, jr,
1173 $ result( jr )
1174 END IF
1175 END IF
1176 220 CONTINUE
1177*
1178 230 CONTINUE
1179 240 CONTINUE
1180*
1181* Summary
1182*
1183 CALL dlasum( 'DGG', nounit, nerrs, ntestt )
1184 RETURN
1185*
1186 9999 FORMAT( ' DCHKGG: ', a, ' returned INFO=', i6, '.', / 9x, 'N=',
1187 $ i6, ', JTYPE=', i6, ', ISEED=(', 3( i5, ',' ), i5, ')' )
1188*
1189 9998 FORMAT( ' DCHKGG: ', a, ' Eigenvectors from ', a, ' incorrectly ',
1190 $ 'normalized.', / ' Bits of error=', 0p, g10.3, ',', 9x,
1191 $ 'N=', i6, ', JTYPE=', i6, ', ISEED=(', 3( i5, ',' ), i5,
1192 $ ')' )
1193*
1194 9997 FORMAT( / 1x, a3, ' -- Real Generalized eigenvalue problem' )
1195*
1196 9996 FORMAT( ' Matrix types (see DCHKGG for details): ' )
1197*
1198 9995 FORMAT( ' Special Matrices:', 23x,
1199 $ '(J''=transposed Jordan block)',
1200 $ / ' 1=(0,0) 2=(I,0) 3=(0,I) 4=(I,I) 5=(J'',J'') ',
1201 $ '6=(diag(J'',I), diag(I,J''))', / ' Diagonal Matrices: ( ',
1202 $ 'D=diag(0,1,2,...) )', / ' 7=(D,I) 9=(large*D, small*I',
1203 $ ') 11=(large*I, small*D) 13=(large*D, large*I)', /
1204 $ ' 8=(I,D) 10=(small*D, large*I) 12=(small*I, large*D) ',
1205 $ ' 14=(small*D, small*I)', / ' 15=(D, reversed D)' )
1206 9994 FORMAT( ' Matrices Rotated by Random ', a, ' Matrices U, V:',
1207 $ / ' 16=Transposed Jordan Blocks 19=geometric ',
1208 $ 'alpha, beta=0,1', / ' 17=arithm. alpha&beta ',
1209 $ ' 20=arithmetic alpha, beta=0,1', / ' 18=clustered ',
1210 $ 'alpha, beta=0,1 21=random alpha, beta=0,1',
1211 $ / ' Large & Small Matrices:', / ' 22=(large, small) ',
1212 $ '23=(small,large) 24=(small,small) 25=(large,large)',
1213 $ / ' 26=random O(1) matrices.' )
1214*
1215 9993 FORMAT( / ' Tests performed: (H is Hessenberg, S is Schur, B, ',
1216 $ 'T, P are triangular,', / 20x, 'U, V, Q, and Z are ', a,
1217 $ ', l and r are the', / 20x,
1218 $ 'appropriate left and right eigenvectors, resp., a is',
1219 $ / 20x, 'alpha, b is beta, and ', a, ' means ', a, '.)',
1220 $ / ' 1 = | A - U H V', a,
1221 $ ' | / ( |A| n ulp ) 2 = | B - U T V', a,
1222 $ ' | / ( |B| n ulp )', / ' 3 = | I - UU', a,
1223 $ ' | / ( n ulp ) 4 = | I - VV', a,
1224 $ ' | / ( n ulp )', / ' 5 = | H - Q S Z', a,
1225 $ ' | / ( |H| n ulp )', 6x, '6 = | T - Q P Z', a,
1226 $ ' | / ( |T| n ulp )', / ' 7 = | I - QQ', a,
1227 $ ' | / ( n ulp ) 8 = | I - ZZ', a,
1228 $ ' | / ( n ulp )', / ' 9 = max | ( b S - a P )', a,
1229 $ ' l | / const. 10 = max | ( b H - a T )', a,
1230 $ ' l | / const.', /
1231 $ ' 11= max | ( b S - a P ) r | / const. 12 = max | ( b H',
1232 $ ' - a T ) r | / const.', / 1x )
1233*
1234 9992 FORMAT( ' Matrix order=', i5, ', type=', i2, ', seed=',
1235 $ 4( i4, ',' ), ' result ', i2, ' is', 0p, f8.2 )
1236 9991 FORMAT( ' Matrix order=', i5, ', type=', i2, ', seed=',
1237 $ 4( i4, ',' ), ' result ', i2, ' is', 1p, d10.3 )
1238*
1239* End of DCHKGG
1240*
double precision function dlamch(CMACH)
DLAMCH
Definition: dlamch.f:69
subroutine dlabad(SMALL, LARGE)
DLABAD
Definition: dlabad.f:74
subroutine dlacpy(UPLO, M, N, A, LDA, B, LDB)
DLACPY copies all or part of one two-dimensional array to another.
Definition: dlacpy.f:103
subroutine dlaset(UPLO, M, N, ALPHA, BETA, A, LDA)
DLASET initializes the off-diagonal elements and the diagonal elements of a matrix to given values.
Definition: dlaset.f:110
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:60
subroutine dlatm4(ITYPE, N, NZ1, NZ2, ISIGN, AMAGN, RCOND, TRIANG, IDIST, ISEED, A, LDA)
DLATM4
Definition: dlatm4.f:175
subroutine dlasum(TYPE, IOUNIT, IE, NRUN)
DLASUM
Definition: dlasum.f:43
subroutine dget52(LEFT, N, A, LDA, B, LDB, E, LDE, ALPHAR, ALPHAI, BETA, WORK, RESULT)
DGET52
Definition: dget52.f:199
subroutine dget51(ITYPE, N, A, LDA, B, LDB, U, LDU, V, LDV, WORK, RESULT)
DGET51
Definition: dget51.f:149
double precision function dlarnd(IDIST, ISEED)
DLARND
Definition: dlarnd.f:73
double precision function dlange(NORM, M, N, A, LDA, WORK)
DLANGE returns the value of the 1-norm, Frobenius norm, infinity-norm, or the largest absolute value ...
Definition: dlange.f:114
subroutine dhgeqz(JOB, COMPQ, COMPZ, N, ILO, IHI, H, LDH, T, LDT, ALPHAR, ALPHAI, BETA, Q, LDQ, Z, LDZ, WORK, LWORK, INFO)
DHGEQZ
Definition: dhgeqz.f:304
subroutine dtgevc(SIDE, HOWMNY, SELECT, N, S, LDS, P, LDP, VL, LDVL, VR, LDVR, MM, M, WORK, INFO)
DTGEVC
Definition: dtgevc.f:295
subroutine dgeqr2(M, N, A, LDA, TAU, WORK, INFO)
DGEQR2 computes the QR factorization of a general rectangular matrix using an unblocked algorithm.
Definition: dgeqr2.f:130
subroutine dlarfg(N, ALPHA, X, INCX, TAU)
DLARFG generates an elementary reflector (Householder matrix).
Definition: dlarfg.f:106
subroutine dgghrd(COMPQ, COMPZ, N, ILO, IHI, A, LDA, B, LDB, Q, LDQ, Z, LDZ, INFO)
DGGHRD
Definition: dgghrd.f:207
subroutine dorm2r(SIDE, TRANS, M, N, K, A, LDA, TAU, C, LDC, WORK, INFO)
DORM2R multiplies a general matrix by the orthogonal matrix from a QR factorization determined by sge...
Definition: dorm2r.f:159
Here is the call graph for this function:
Here is the caller graph for this function: