LAPACK  3.10.0 LAPACK: Linear Algebra PACKage

## ◆ sdrgev3()

 subroutine sdrgev3 ( integer NSIZES, integer, dimension( * ) NN, integer NTYPES, logical, dimension( * ) DOTYPE, integer, dimension( 4 ) ISEED, real THRESH, integer NOUNIT, real, dimension( lda, * ) A, integer LDA, real, dimension( lda, * ) B, real, dimension( lda, * ) S, real, dimension( lda, * ) T, real, dimension( ldq, * ) Q, integer LDQ, real, dimension( ldq, * ) Z, real, dimension( ldqe, * ) QE, integer LDQE, real, dimension( * ) ALPHAR, real, dimension( * ) ALPHAI, real, dimension( * ) BETA, real, dimension( * ) ALPHR1, real, dimension( * ) ALPHI1, real, dimension( * ) BETA1, real, dimension( * ) WORK, integer LWORK, real, dimension( * ) RESULT, integer INFO )

SDRGEV3

Purpose:
``` SDRGEV3 checks the nonsymmetric generalized eigenvalue problem driver
routine SGGEV3.

SGGEV3 computes for a pair of n-by-n nonsymmetric matrices (A,B) the
generalized eigenvalues and, optionally, the left and right
eigenvectors.

A generalized eigenvalue for a pair of matrices (A,B) is a scalar w
or a ratio  alpha/beta = w, such that A - w*B is singular.  It is
usually represented as the pair (alpha,beta), as there is reasonable
interpretation for beta=0, and even for both being zero.

A right generalized eigenvector corresponding to a generalized
eigenvalue  w  for a pair of matrices (A,B) is a vector r  such that
(A - wB) * r = 0.  A left generalized eigenvector is a vector l such
that l**H * (A - wB) = 0, where l**H is the conjugate-transpose of l.

When SDRGEV3 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, a pair of matrices (A, B) will be generated
and used for testing.  For each matrix pair, the following tests
will be performed and compared with the threshold THRESH.

Results from SGGEV3:

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

| VL**H * (beta A - alpha B) |/( ulp max(|beta A|, |alpha B|) )

where VL**H is the conjugate-transpose of VL.

(2)  | |VL(i)| - 1 | / ulp and whether largest component real

VL(i) denotes the i-th column of VL.

(3)  max over all left eigenvalue/-vector pairs (alpha/beta,r) of

| (beta A - alpha B) * VR | / ( ulp max(|beta A|, |alpha B|) )

(4)  | |VR(i)| - 1 | / ulp and whether largest component real

VR(i) denotes the i-th column of VR.

(5)  W(full) = W(partial)
W(full) denotes the eigenvalues computed when both l and r
are also computed, and W(partial) denotes the eigenvalues
computed when only W, only W and r, or only W and l are
computed.

(6)  VL(full) = VL(partial)
VL(full) denotes the left eigenvectors computed when both l
and r are computed, and VL(partial) denotes the result
when only l is computed.

(7)  VR(full) = VR(partial)
VR(full) denotes the right eigenvectors computed when both l
and r are also computed, and VR(partial) denotes the result
when only l is computed.

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) Q ( J , J ) Z     where Q and Z are random orthogonal matrices.

(17) Q ( T1, T2 ) Z    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) Q ( T1, T2 ) Z    diag(T1) = ( 0, 0, 1, 1, s, ..., s, 0 )
diag(T2) = ( 0, 1, 0, 1,..., 1, 0 )
s = machine precision.

(19) Q ( T1, T2 ) Z    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) Q ( T1, T2 ) Z    diag(T1)=( 0, 0, 1, 1, a, ..., a   =s, 0 )
diag(T2) = ( 0, 1, 0, 1, ..., 1, 0, 0 )

(21) Q ( T1, T2 ) Z    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) Q ( big*T1, small*T2 ) Z    diag(T1) = ( 0, 0, 1, ..., N-3, 0 )
diag(T2) = ( 0, 1, ..., 1, 0, 0 )

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

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

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

(26) Q ( T1, T2 ) Z     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, SDRGEV3 does nothing. NSIZES >= 0.``` [in] NN ``` NN is INTEGER array, dimension (NSIZES) An array containing the sizes to be used for the matrices. Zero values will be skipped. NN >= 0.``` [in] NTYPES ``` NTYPES is INTEGER The number of elements in DOTYPE. If it is zero, SDRGEV3 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 SDRGEV3 to continue the same random number sequence.``` [in] THRESH ``` THRESH is REAL 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] NOUNIT ``` NOUNIT is INTEGER The FORTRAN unit number for printing out error messages (e.g., if a routine returns IERR not equal to 0.)``` [in,out] A ``` A is REAL 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, S, and T. It must be at least 1 and at least max( NN ).``` [in,out] B ``` B is REAL 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] S ``` S is REAL array, dimension (LDA, max(NN)) The Schur form matrix computed from A by SGGEV3. On exit, S contains the Schur form matrix corresponding to the matrix in A.``` [out] T ``` T is REAL array, dimension (LDA, max(NN)) The upper triangular matrix computed from B by SGGEV3.``` [out] Q ``` Q is REAL array, dimension (LDQ, max(NN)) The (left) eigenvectors matrix computed by SGGEV3.``` [in] LDQ ``` LDQ is INTEGER The leading dimension of Q and Z. It must be at least 1 and at least max( NN ).``` [out] Z ``` Z is REAL array, dimension( LDQ, max(NN) ) The (right) orthogonal matrix computed by SGGEV3.``` [out] QE ``` QE is REAL array, dimension( LDQ, max(NN) ) QE holds the computed right or left eigenvectors.``` [in] LDQE ``` LDQE is INTEGER The leading dimension of QE. LDQE >= max(1,max(NN)).``` [out] ALPHAR ` ALPHAR is REAL array, dimension (max(NN))` [out] ALPHAI ` ALPHAI is REAL array, dimension (max(NN))` [out] BETA ``` BETA is REAL array, dimension (max(NN)) \verbatim The generalized eigenvalues of (A,B) computed by SGGEV3. ( ALPHAR(k)+ALPHAI(k)*i ) / BETA(k) is the k-th generalized eigenvalue of A and B.``` [out] ALPHR1 ` ALPHR1 is REAL array, dimension (max(NN))` [out] ALPHI1 ` ALPHI1 is REAL array, dimension (max(NN))` [out] BETA1 ``` BETA1 is REAL array, dimension (max(NN)) Like ALPHAR, ALPHAI, BETA, these arrays contain the eigenvalues of A and B, but those computed when SGGEV3 only computes a partial eigendecomposition, i.e. not the eigenvalues and left and right eigenvectors.``` [out] WORK ` WORK is REAL array, dimension (LWORK)` [in] LWORK ``` LWORK is INTEGER The number of entries in WORK. LWORK >= MAX( 8*N, N*(N+1) ).``` [out] RESULT ``` RESULT is REAL array, dimension (2) 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.```

Definition at line 404 of file sdrgev3.f.

408 *
409 * -- LAPACK test routine --
410 * -- LAPACK is a software package provided by Univ. of Tennessee, --
411 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
412 *
413 * .. Scalar Arguments ..
414  INTEGER INFO, LDA, LDQ, LDQE, LWORK, NOUNIT, NSIZES,
415  \$ NTYPES
416  REAL THRESH
417 * ..
418 * .. Array Arguments ..
419  LOGICAL DOTYPE( * )
420  INTEGER ISEED( 4 ), NN( * )
421  REAL A( LDA, * ), ALPHAI( * ), ALPHI1( * ),
422  \$ ALPHAR( * ), ALPHR1( * ), B( LDA, * ),
423  \$ BETA( * ), BETA1( * ), Q( LDQ, * ),
424  \$ QE( LDQE, * ), RESULT( * ), S( LDA, * ),
425  \$ T( LDA, * ), WORK( * ), Z( LDQ, * )
426 * ..
427 *
428 * =====================================================================
429 *
430 * .. Parameters ..
431  REAL ZERO, ONE
432  parameter( zero = 0.0e+0, one = 1.0e+0 )
433  INTEGER MAXTYP
434  parameter( maxtyp = 26 )
435 * ..
436 * .. Local Scalars ..
438  INTEGER I, IADD, IERR, IN, J, JC, JR, JSIZE, JTYPE,
439  \$ MAXWRK, MINWRK, MTYPES, N, N1, NERRS, NMATS,
440  \$ NMAX, NTESTT
441  REAL SAFMAX, SAFMIN, ULP, ULPINV
442 * ..
443 * .. Local Arrays ..
444  INTEGER IASIGN( MAXTYP ), IBSIGN( MAXTYP ),
445  \$ IOLDSD( 4 ), KADD( 6 ), KAMAGN( MAXTYP ),
446  \$ KATYPE( MAXTYP ), KAZERO( MAXTYP ),
447  \$ KBMAGN( MAXTYP ), KBTYPE( MAXTYP ),
448  \$ KBZERO( MAXTYP ), KCLASS( MAXTYP ),
449  \$ KTRIAN( MAXTYP ), KZ1( 6 ), KZ2( 6 )
450  REAL RMAGN( 0: 3 )
451 * ..
452 * .. External Functions ..
453  INTEGER ILAENV
454  REAL SLAMCH, SLARND
455  EXTERNAL ilaenv, slamch, slarnd
456 * ..
457 * .. External Subroutines ..
458  EXTERNAL alasvm, sget52, sggev3, slabad, slacpy, slarfg,
460 * ..
461 * .. Intrinsic Functions ..
462  INTRINSIC abs, max, min, real, sign
463 * ..
464 * .. Data statements ..
465  DATA kclass / 15*1, 10*2, 1*3 /
466  DATA kz1 / 0, 1, 2, 1, 3, 3 /
467  DATA kz2 / 0, 0, 1, 2, 1, 1 /
468  DATA kadd / 0, 0, 0, 0, 3, 2 /
469  DATA katype / 0, 1, 0, 1, 2, 3, 4, 1, 4, 4, 1, 1, 4,
470  \$ 4, 4, 2, 4, 5, 8, 7, 9, 4*4, 0 /
471  DATA kbtype / 0, 0, 1, 1, 2, -3, 1, 4, 1, 1, 4, 4,
472  \$ 1, 1, -4, 2, -4, 8*8, 0 /
473  DATA kazero / 6*1, 2, 1, 2*2, 2*1, 2*2, 3, 1, 3,
474  \$ 4*5, 4*3, 1 /
475  DATA kbzero / 6*1, 1, 2, 2*1, 2*2, 2*1, 4, 1, 4,
476  \$ 4*6, 4*4, 1 /
477  DATA kamagn / 8*1, 2, 3, 2, 3, 2, 3, 7*1, 2, 3, 3,
478  \$ 2, 1 /
479  DATA kbmagn / 8*1, 3, 2, 3, 2, 2, 3, 7*1, 3, 2, 3,
480  \$ 2, 1 /
481  DATA ktrian / 16*0, 10*1 /
482  DATA iasign / 6*0, 2, 0, 2*2, 2*0, 3*2, 0, 2, 3*0,
483  \$ 5*2, 0 /
484  DATA ibsign / 7*0, 2, 2*0, 2*2, 2*0, 2, 0, 2, 9*0 /
485 * ..
486 * .. Executable Statements ..
487 *
488 * Check for errors
489 *
490  info = 0
491 *
493  nmax = 1
494  DO 10 j = 1, nsizes
495  nmax = max( nmax, nn( j ) )
496  IF( nn( j ).LT.0 )
498  10 CONTINUE
499 *
500  IF( nsizes.LT.0 ) THEN
501  info = -1
502  ELSE IF( badnn ) THEN
503  info = -2
504  ELSE IF( ntypes.LT.0 ) THEN
505  info = -3
506  ELSE IF( thresh.LT.zero ) THEN
507  info = -6
508  ELSE IF( lda.LE.1 .OR. lda.LT.nmax ) THEN
509  info = -9
510  ELSE IF( ldq.LE.1 .OR. ldq.LT.nmax ) THEN
511  info = -14
512  ELSE IF( ldqe.LE.1 .OR. ldqe.LT.nmax ) THEN
513  info = -17
514  END IF
515 *
516 * Compute workspace
517 * (Note: Comments in the code beginning "Workspace:" describe the
518 * minimal amount of workspace needed at that point in the code,
519 * as well as the preferred amount for good performance.
520 * NB refers to the optimal block size for the immediately
521 * following subroutine, as returned by ILAENV.
522 *
523  minwrk = 1
524  IF( info.EQ.0 .AND. lwork.GE.1 ) THEN
525  minwrk = max( 1, 8*nmax, nmax*( nmax+1 ) )
526  maxwrk = 7*nmax + nmax*ilaenv( 1, 'SGEQRF', ' ', nmax, 1, nmax,
527  \$ 0 )
528  maxwrk = max( maxwrk, nmax*( nmax+1 ) )
529  work( 1 ) = maxwrk
530  END IF
531 *
532  IF( lwork.LT.minwrk )
533  \$ info = -25
534 *
535  IF( info.NE.0 ) THEN
536  CALL xerbla( 'SDRGEV3', -info )
537  RETURN
538  END IF
539 *
540 * Quick return if possible
541 *
542  IF( nsizes.EQ.0 .OR. ntypes.EQ.0 )
543  \$ RETURN
544 *
545  safmin = slamch( 'Safe minimum' )
546  ulp = slamch( 'Epsilon' )*slamch( 'Base' )
547  safmin = safmin / ulp
548  safmax = one / safmin
549  CALL slabad( safmin, safmax )
550  ulpinv = one / ulp
551 *
552 * The values RMAGN(2:3) depend on N, see below.
553 *
554  rmagn( 0 ) = zero
555  rmagn( 1 ) = one
556 *
557 * Loop over sizes, types
558 *
559  ntestt = 0
560  nerrs = 0
561  nmats = 0
562 *
563  DO 220 jsize = 1, nsizes
564  n = nn( jsize )
565  n1 = max( 1, n )
566  rmagn( 2 ) = safmax*ulp / real( n1 )
567  rmagn( 3 ) = safmin*ulpinv*n1
568 *
569  IF( nsizes.NE.1 ) THEN
570  mtypes = min( maxtyp, ntypes )
571  ELSE
572  mtypes = min( maxtyp+1, ntypes )
573  END IF
574 *
575  DO 210 jtype = 1, mtypes
576  IF( .NOT.dotype( jtype ) )
577  \$ GO TO 210
578  nmats = nmats + 1
579 *
580 * Save ISEED in case of an error.
581 *
582  DO 20 j = 1, 4
583  ioldsd( j ) = iseed( j )
584  20 CONTINUE
585 *
586 * Generate test matrices A and B
587 *
588 * Description of control parameters:
589 *
590 * KCLASS: =1 means w/o rotation, =2 means w/ rotation,
591 * =3 means random.
592 * KATYPE: the "type" to be passed to SLATM4 for computing A.
593 * KAZERO: the pattern of zeros on the diagonal for A:
594 * =1: ( xxx ), =2: (0, xxx ) =3: ( 0, 0, xxx, 0 ),
595 * =4: ( 0, xxx, 0, 0 ), =5: ( 0, 0, 1, xxx, 0 ),
596 * =6: ( 0, 1, 0, xxx, 0 ). (xxx means a string of
597 * non-zero entries.)
598 * KAMAGN: the magnitude of the matrix: =0: zero, =1: O(1),
599 * =2: large, =3: small.
600 * IASIGN: 1 if the diagonal elements of A are to be
601 * multiplied by a random magnitude 1 number, =2 if
602 * randomly chosen diagonal blocks are to be rotated
603 * to form 2x2 blocks.
604 * KBTYPE, KBZERO, KBMAGN, IBSIGN: the same, but for B.
605 * KTRIAN: =0: don't fill in the upper triangle, =1: do.
606 * KZ1, KZ2, KADD: used to implement KAZERO and KBZERO.
607 * RMAGN: used to implement KAMAGN and KBMAGN.
608 *
609  IF( mtypes.GT.maxtyp )
610  \$ GO TO 100
611  ierr = 0
612  IF( kclass( jtype ).LT.3 ) THEN
613 *
614 * Generate A (w/o rotation)
615 *
616  IF( abs( katype( jtype ) ).EQ.3 ) THEN
617  in = 2*( ( n-1 ) / 2 ) + 1
618  IF( in.NE.n )
619  \$ CALL slaset( 'Full', n, n, zero, zero, a, lda )
620  ELSE
621  in = n
622  END IF
623  CALL slatm4( katype( jtype ), in, kz1( kazero( jtype ) ),
624  \$ kz2( kazero( jtype ) ), iasign( jtype ),
625  \$ rmagn( kamagn( jtype ) ), ulp,
626  \$ rmagn( ktrian( jtype )*kamagn( jtype ) ), 2,
627  \$ iseed, a, lda )
631 *
632 * Generate B (w/o rotation)
633 *
634  IF( abs( kbtype( jtype ) ).EQ.3 ) THEN
635  in = 2*( ( n-1 ) / 2 ) + 1
636  IF( in.NE.n )
637  \$ CALL slaset( 'Full', n, n, zero, zero, b, lda )
638  ELSE
639  in = n
640  END IF
641  CALL slatm4( kbtype( jtype ), in, kz1( kbzero( jtype ) ),
642  \$ kz2( kbzero( jtype ) ), ibsign( jtype ),
643  \$ rmagn( kbmagn( jtype ) ), one,
644  \$ rmagn( ktrian( jtype )*kbmagn( jtype ) ), 2,
645  \$ iseed, b, lda )
649 *
650  IF( kclass( jtype ).EQ.2 .AND. n.GT.0 ) THEN
651 *
652 * Include rotations
653 *
654 * Generate Q, Z as Householder transformations times
655 * a diagonal matrix.
656 *
657  DO 40 jc = 1, n - 1
658  DO 30 jr = jc, n
659  q( jr, jc ) = slarnd( 3, iseed )
660  z( jr, jc ) = slarnd( 3, iseed )
661  30 CONTINUE
662  CALL slarfg( n+1-jc, q( jc, jc ), q( jc+1, jc ), 1,
663  \$ work( jc ) )
664  work( 2*n+jc ) = sign( one, q( jc, jc ) )
665  q( jc, jc ) = one
666  CALL slarfg( n+1-jc, z( jc, jc ), z( jc+1, jc ), 1,
667  \$ work( n+jc ) )
668  work( 3*n+jc ) = sign( one, z( jc, jc ) )
669  z( jc, jc ) = one
670  40 CONTINUE
671  q( n, n ) = one
672  work( n ) = zero
673  work( 3*n ) = sign( one, slarnd( 2, iseed ) )
674  z( n, n ) = one
675  work( 2*n ) = zero
676  work( 4*n ) = sign( one, slarnd( 2, iseed ) )
677 *
678 * Apply the diagonal matrices
679 *
680  DO 60 jc = 1, n
681  DO 50 jr = 1, n
682  a( jr, jc ) = work( 2*n+jr )*work( 3*n+jc )*
683  \$ a( jr, jc )
684  b( jr, jc ) = work( 2*n+jr )*work( 3*n+jc )*
685  \$ b( jr, jc )
686  50 CONTINUE
687  60 CONTINUE
688  CALL sorm2r( 'L', 'N', n, n, n-1, q, ldq, work, a,
689  \$ lda, work( 2*n+1 ), ierr )
690  IF( ierr.NE.0 )
691  \$ GO TO 90
692  CALL sorm2r( 'R', 'T', n, n, n-1, z, ldq, work( n+1 ),
693  \$ a, lda, work( 2*n+1 ), ierr )
694  IF( ierr.NE.0 )
695  \$ GO TO 90
696  CALL sorm2r( 'L', 'N', n, n, n-1, q, ldq, work, b,
697  \$ lda, work( 2*n+1 ), ierr )
698  IF( ierr.NE.0 )
699  \$ GO TO 90
700  CALL sorm2r( 'R', 'T', n, n, n-1, z, ldq, work( n+1 ),
701  \$ b, lda, work( 2*n+1 ), ierr )
702  IF( ierr.NE.0 )
703  \$ GO TO 90
704  END IF
705  ELSE
706 *
707 * Random matrices
708 *
709  DO 80 jc = 1, n
710  DO 70 jr = 1, n
711  a( jr, jc ) = rmagn( kamagn( jtype ) )*
712  \$ slarnd( 2, iseed )
713  b( jr, jc ) = rmagn( kbmagn( jtype ) )*
714  \$ slarnd( 2, iseed )
715  70 CONTINUE
716  80 CONTINUE
717  END IF
718 *
719  90 CONTINUE
720 *
721  IF( ierr.NE.0 ) THEN
722  WRITE( nounit, fmt = 9999 )'Generator', ierr, n, jtype,
723  \$ ioldsd
724  info = abs( ierr )
725  RETURN
726  END IF
727 *
728  100 CONTINUE
729 *
730  DO 110 i = 1, 7
731  result( i ) = -one
732  110 CONTINUE
733 *
734 * Call XLAENV to set the parameters used in SLAQZ0
735 *
736  CALL xlaenv( 12, 10 )
737  CALL xlaenv( 13, 12 )
738  CALL xlaenv( 14, 13 )
739  CALL xlaenv( 15, 2 )
740  CALL xlaenv( 17, 10 )
741 *
742 * Call SGGEV3 to compute eigenvalues and eigenvectors.
743 *
744  CALL slacpy( ' ', n, n, a, lda, s, lda )
745  CALL slacpy( ' ', n, n, b, lda, t, lda )
746  CALL sggev3( 'V', 'V', n, s, lda, t, lda, alphar, alphai,
747  \$ beta, q, ldq, z, ldq, work, lwork, ierr )
748  IF( ierr.NE.0 .AND. ierr.NE.n+1 ) THEN
749  result( 1 ) = ulpinv
750  WRITE( nounit, fmt = 9999 )'SGGEV31', ierr, n, jtype,
751  \$ ioldsd
752  info = abs( ierr )
753  GO TO 190
754  END IF
755 *
756 * Do the tests (1) and (2)
757 *
758  CALL sget52( .true., n, a, lda, b, lda, q, ldq, alphar,
759  \$ alphai, beta, work, result( 1 ) )
760  IF( result( 2 ).GT.thresh ) THEN
761  WRITE( nounit, fmt = 9998 )'Left', 'SGGEV31',
762  \$ result( 2 ), n, jtype, ioldsd
763  END IF
764 *
765 * Do the tests (3) and (4)
766 *
767  CALL sget52( .false., n, a, lda, b, lda, z, ldq, alphar,
768  \$ alphai, beta, work, result( 3 ) )
769  IF( result( 4 ).GT.thresh ) THEN
770  WRITE( nounit, fmt = 9998 )'Right', 'SGGEV31',
771  \$ result( 4 ), n, jtype, ioldsd
772  END IF
773 *
774 * Do the test (5)
775 *
776  CALL slacpy( ' ', n, n, a, lda, s, lda )
777  CALL slacpy( ' ', n, n, b, lda, t, lda )
778  CALL sggev3( 'N', 'N', n, s, lda, t, lda, alphr1, alphi1,
779  \$ beta1, q, ldq, z, ldq, work, lwork, ierr )
780  IF( ierr.NE.0 .AND. ierr.NE.n+1 ) THEN
781  result( 1 ) = ulpinv
782  WRITE( nounit, fmt = 9999 )'SGGEV32', ierr, n, jtype,
783  \$ ioldsd
784  info = abs( ierr )
785  GO TO 190
786  END IF
787 *
788  DO 120 j = 1, n
789  IF( alphar( j ).NE.alphr1( j ) .OR.
790  \$ beta( j ).NE. beta1( j ) ) THEN
791  result( 5 ) = ulpinv
792  END IF
793  120 CONTINUE
794 *
795 * Do the test (6): Compute eigenvalues and left eigenvectors,
796 * and test them
797 *
798  CALL slacpy( ' ', n, n, a, lda, s, lda )
799  CALL slacpy( ' ', n, n, b, lda, t, lda )
800  CALL sggev3( 'V', 'N', n, s, lda, t, lda, alphr1, alphi1,
801  \$ beta1, qe, ldqe, z, ldq, work, lwork, ierr )
802  IF( ierr.NE.0 .AND. ierr.NE.n+1 ) THEN
803  result( 1 ) = ulpinv
804  WRITE( nounit, fmt = 9999 )'SGGEV33', ierr, n, jtype,
805  \$ ioldsd
806  info = abs( ierr )
807  GO TO 190
808  END IF
809 *
810  DO 130 j = 1, n
811  IF( alphar( j ).NE.alphr1( j ) .OR. alphai( j ).NE.
812  \$ alphi1( j ) .OR. beta( j ).NE.beta1( j ) )
813  \$ result( 6 ) = ulpinv
814  130 CONTINUE
815 *
816  DO 150 j = 1, n
817  DO 140 jc = 1, n
818  IF( q( j, jc ).NE.qe( j, jc ) )
819  \$ result( 6 ) = ulpinv
820  140 CONTINUE
821  150 CONTINUE
822 *
823 * DO the test (7): Compute eigenvalues and right eigenvectors,
824 * and test them
825 *
826  CALL slacpy( ' ', n, n, a, lda, s, lda )
827  CALL slacpy( ' ', n, n, b, lda, t, lda )
828  CALL sggev3( 'N', 'V', n, s, lda, t, lda, alphr1, alphi1,
829  \$ beta1, q, ldq, qe, ldqe, work, lwork, ierr )
830  IF( ierr.NE.0 .AND. ierr.NE.n+1 ) THEN
831  result( 1 ) = ulpinv
832  WRITE( nounit, fmt = 9999 )'SGGEV34', ierr, n, jtype,
833  \$ ioldsd
834  info = abs( ierr )
835  GO TO 190
836  END IF
837 *
838  DO 160 j = 1, n
839  IF( alphar( j ).NE.alphr1( j ) .OR. alphai( j ).NE.
840  \$ alphi1( j ) .OR. beta( j ).NE.beta1( j ) )
841  \$ result( 7 ) = ulpinv
842  160 CONTINUE
843 *
844  DO 180 j = 1, n
845  DO 170 jc = 1, n
846  IF( z( j, jc ).NE.qe( j, jc ) )
847  \$ result( 7 ) = ulpinv
848  170 CONTINUE
849  180 CONTINUE
850 *
851 * End of Loop -- Check for RESULT(j) > THRESH
852 *
853  190 CONTINUE
854 *
855  ntestt = ntestt + 7
856 *
857 * Print out tests which fail.
858 *
859  DO 200 jr = 1, 7
860  IF( result( jr ).GE.thresh ) THEN
861 *
862 * If this is the first test to fail,
863 * print a header to the data file.
864 *
865  IF( nerrs.EQ.0 ) THEN
866  WRITE( nounit, fmt = 9997 )'SGV'
867 *
868 * Matrix types
869 *
870  WRITE( nounit, fmt = 9996 )
871  WRITE( nounit, fmt = 9995 )
872  WRITE( nounit, fmt = 9994 )'Orthogonal'
873 *
874 * Tests performed
875 *
876  WRITE( nounit, fmt = 9993 )
877 *
878  END IF
879  nerrs = nerrs + 1
880  IF( result( jr ).LT.10000.0 ) THEN
881  WRITE( nounit, fmt = 9992 )n, jtype, ioldsd, jr,
882  \$ result( jr )
883  ELSE
884  WRITE( nounit, fmt = 9991 )n, jtype, ioldsd, jr,
885  \$ result( jr )
886  END IF
887  END IF
888  200 CONTINUE
889 *
890  210 CONTINUE
891  220 CONTINUE
892 *
893 * Summary
894 *
895  CALL alasvm( 'SGV', nounit, nerrs, ntestt, 0 )
896 *
897  work( 1 ) = maxwrk
898 *
899  RETURN
900 *
901  9999 FORMAT( ' SDRGEV3: ', a, ' returned INFO=', i6, '.', / 3x, 'N=',
902  \$ i6, ', JTYPE=', i6, ', ISEED=(', 4( i4, ',' ), i5, ')' )
903 *
904  9998 FORMAT( ' SDRGEV3: ', a, ' Eigenvectors from ', a,
905  \$ ' incorrectly normalized.', / ' Bits of error=', 0p, g10.3,
906  \$ ',', 3x, 'N=', i4, ', JTYPE=', i3, ', ISEED=(',
907  \$ 4( i4, ',' ), i5, ')' )
908 *
909  9997 FORMAT( / 1x, a3, ' -- Real Generalized eigenvalue problem driver'
910  \$ )
911 *
912  9996 FORMAT( ' Matrix types (see SDRGEV3 for details): ' )
913 *
914  9995 FORMAT( ' Special Matrices:', 23x,
915  \$ '(J''=transposed Jordan block)',
916  \$ / ' 1=(0,0) 2=(I,0) 3=(0,I) 4=(I,I) 5=(J'',J'') ',
917  \$ '6=(diag(J'',I), diag(I,J''))', / ' Diagonal Matrices: ( ',
918  \$ 'D=diag(0,1,2,...) )', / ' 7=(D,I) 9=(large*D, small*I',
919  \$ ') 11=(large*I, small*D) 13=(large*D, large*I)', /
920  \$ ' 8=(I,D) 10=(small*D, large*I) 12=(small*I, large*D) ',
921  \$ ' 14=(small*D, small*I)', / ' 15=(D, reversed D)' )
922  9994 FORMAT( ' Matrices Rotated by Random ', a, ' Matrices U, V:',
923  \$ / ' 16=Transposed Jordan Blocks 19=geometric ',
924  \$ 'alpha, beta=0,1', / ' 17=arithm. alpha&beta ',
925  \$ ' 20=arithmetic alpha, beta=0,1', / ' 18=clustered ',
926  \$ 'alpha, beta=0,1 21=random alpha, beta=0,1',
927  \$ / ' Large & Small Matrices:', / ' 22=(large, small) ',
928  \$ '23=(small,large) 24=(small,small) 25=(large,large)',
929  \$ / ' 26=random O(1) matrices.' )
930 *
931  9993 FORMAT( / ' Tests performed: ',
932  \$ / ' 1 = max | ( b A - a B )''*l | / const.,',
933  \$ / ' 2 = | |VR(i)| - 1 | / ulp,',
934  \$ / ' 3 = max | ( b A - a B )*r | / const.',
935  \$ / ' 4 = | |VL(i)| - 1 | / ulp,',
936  \$ / ' 5 = 0 if W same no matter if r or l computed,',
937  \$ / ' 6 = 0 if l same no matter if l computed,',
938  \$ / ' 7 = 0 if r same no matter if r computed,', / 1x )
939  9992 FORMAT( ' Matrix order=', i5, ', type=', i2, ', seed=',
940  \$ 4( i4, ',' ), ' result ', i2, ' is', 0p, f8.2 )
941  9991 FORMAT( ' Matrix order=', i5, ', type=', i2, ', seed=',
942  \$ 4( i4, ',' ), ' result ', i2, ' is', 1p, e10.3 )
943 *
944 * End of SDRGEV3
945 *
subroutine slaset(UPLO, M, N, ALPHA, BETA, A, LDA)
SLASET initializes the off-diagonal elements and the diagonal elements of a matrix to given values.
Definition: slaset.f:110
subroutine slacpy(UPLO, M, N, A, LDA, B, LDB)
SLACPY copies all or part of one two-dimensional array to another.
Definition: slacpy.f:103
integer function ilaenv(ISPEC, NAME, OPTS, N1, N2, N3, N4)
ILAENV
Definition: ilaenv.f:162
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:60
subroutine alasvm(TYPE, NOUT, NFAIL, NRUN, NERRS)
ALASVM
Definition: alasvm.f:73
subroutine xlaenv(ISPEC, NVALUE)
XLAENV
Definition: xlaenv.f:81
real function slarnd(IDIST, ISEED)
SLARND
Definition: slarnd.f:73
subroutine sggev3(JOBVL, JOBVR, N, A, LDA, B, LDB, ALPHAR, ALPHAI, BETA, VL, LDVL, VR, LDVR, WORK, LWORK, INFO)
SGGEV3 computes the eigenvalues and, optionally, the left and/or right eigenvectors for GE matrices (...
Definition: sggev3.f:226
subroutine slarfg(N, ALPHA, X, INCX, TAU)
SLARFG generates an elementary reflector (Householder matrix).
Definition: slarfg.f:106
subroutine sorm2r(SIDE, TRANS, M, N, K, A, LDA, TAU, C, LDC, WORK, INFO)
SORM2R multiplies a general matrix by the orthogonal matrix from a QR factorization determined by sge...
Definition: sorm2r.f:159
subroutine sget52(LEFT, N, A, LDA, B, LDB, E, LDE, ALPHAR, ALPHAI, BETA, WORK, RESULT)
SGET52
Definition: sget52.f:199
subroutine slatm4(ITYPE, N, NZ1, NZ2, ISIGN, AMAGN, RCOND, TRIANG, IDIST, ISEED, A, LDA)
SLATM4
Definition: slatm4.f:175
real function slamch(CMACH)
SLAMCH
Definition: slamch.f:68
Here is the call graph for this function:
Here is the caller graph for this function: