LAPACK  3.10.0
LAPACK: Linear Algebra PACKage

◆ cget24()

subroutine cget24 ( logical  COMP,
integer  JTYPE,
real  THRESH,
integer, dimension( 4 )  ISEED,
integer  NOUNIT,
integer  N,
complex, dimension( lda, * )  A,
integer  LDA,
complex, dimension( lda, * )  H,
complex, dimension( lda, * )  HT,
complex, dimension( * )  W,
complex, dimension( * )  WT,
complex, dimension( * )  WTMP,
complex, dimension( ldvs, * )  VS,
integer  LDVS,
complex, dimension( ldvs, * )  VS1,
real  RCDEIN,
real  RCDVIN,
integer  NSLCT,
integer, dimension( * )  ISLCT,
integer  ISRT,
real, dimension( 17 )  RESULT,
complex, dimension( * )  WORK,
integer  LWORK,
real, dimension( * )  RWORK,
logical, dimension( * )  BWORK,
integer  INFO 
)

CGET24

Purpose:
    CGET24 checks the nonsymmetric eigenvalue (Schur form) problem
    expert driver CGEESX.

    If COMP = .FALSE., the first 13 of the following tests will be
    be performed on the input matrix A, and also tests 14 and 15
    if LWORK is sufficiently large.
    If COMP = .TRUE., all 17 test will be performed.

    (1)     0 if T is in Schur form, 1/ulp otherwise
           (no sorting of eigenvalues)

    (2)     | A - VS T VS' | / ( n |A| ulp )

      Here VS is the matrix of Schur eigenvectors, and T is in Schur
      form  (no sorting of eigenvalues).

    (3)     | I - VS VS' | / ( n ulp ) (no sorting of eigenvalues).

    (4)     0     if W are eigenvalues of T
            1/ulp otherwise
            (no sorting of eigenvalues)

    (5)     0     if T(with VS) = T(without VS),
            1/ulp otherwise
            (no sorting of eigenvalues)

    (6)     0     if eigenvalues(with VS) = eigenvalues(without VS),
            1/ulp otherwise
            (no sorting of eigenvalues)

    (7)     0 if T is in Schur form, 1/ulp otherwise
            (with sorting of eigenvalues)

    (8)     | A - VS T VS' | / ( n |A| ulp )

      Here VS is the matrix of Schur eigenvectors, and T is in Schur
      form  (with sorting of eigenvalues).

    (9)     | I - VS VS' | / ( n ulp ) (with sorting of eigenvalues).

    (10)    0     if W are eigenvalues of T
            1/ulp otherwise
            If workspace sufficient, also compare W with and
            without reciprocal condition numbers
            (with sorting of eigenvalues)

    (11)    0     if T(with VS) = T(without VS),
            1/ulp otherwise
            If workspace sufficient, also compare T with and without
            reciprocal condition numbers
            (with sorting of eigenvalues)

    (12)    0     if eigenvalues(with VS) = eigenvalues(without VS),
            1/ulp otherwise
            If workspace sufficient, also compare VS with and without
            reciprocal condition numbers
            (with sorting of eigenvalues)

    (13)    if sorting worked and SDIM is the number of
            eigenvalues which were SELECTed
            If workspace sufficient, also compare SDIM with and
            without reciprocal condition numbers

    (14)    if RCONDE the same no matter if VS and/or RCONDV computed

    (15)    if RCONDV the same no matter if VS and/or RCONDE computed

    (16)  |RCONDE - RCDEIN| / cond(RCONDE)

       RCONDE is the reciprocal average eigenvalue condition number
       computed by CGEESX and RCDEIN (the precomputed true value)
       is supplied as input.  cond(RCONDE) is the condition number
       of RCONDE, and takes errors in computing RCONDE into account,
       so that the resulting quantity should be O(ULP). cond(RCONDE)
       is essentially given by norm(A)/RCONDV.

    (17)  |RCONDV - RCDVIN| / cond(RCONDV)

       RCONDV is the reciprocal right invariant subspace condition
       number computed by CGEESX and RCDVIN (the precomputed true
       value) is supplied as input. cond(RCONDV) is the condition
       number of RCONDV, and takes errors in computing RCONDV into
       account, so that the resulting quantity should be O(ULP).
       cond(RCONDV) is essentially given by norm(A)/RCONDE.
Parameters
[in]COMP
          COMP is LOGICAL
          COMP describes which input tests to perform:
            = .FALSE. if the computed condition numbers are not to
                      be tested against RCDVIN and RCDEIN
            = .TRUE.  if they are to be compared
[in]JTYPE
          JTYPE is INTEGER
          Type of input matrix. Used to label output if error occurs.
[in]ISEED
          ISEED is INTEGER array, dimension (4)
          If COMP = .FALSE., the random number generator seed
          used to produce matrix.
          If COMP = .TRUE., ISEED(1) = the number of the example.
          Used to label output if error occurs.
[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 INFO not equal to 0.)
[in]N
          N is INTEGER
          The dimension of A. N must be at least 0.
[in,out]A
          A is COMPLEX array, dimension (LDA, N)
          Used to hold the matrix whose eigenvalues are to be
          computed.
[in]LDA
          LDA is INTEGER
          The leading dimension of A, and H. LDA must be at
          least 1 and at least N.
[out]H
          H is COMPLEX array, dimension (LDA, N)
          Another copy of the test matrix A, modified by CGEESX.
[out]HT
          HT is COMPLEX array, dimension (LDA, N)
          Yet another copy of the test matrix A, modified by CGEESX.
[out]W
          W is COMPLEX array, dimension (N)
          The computed eigenvalues of A.
[out]WT
          WT is COMPLEX array, dimension (N)
          Like W, this array contains the eigenvalues of A,
          but those computed when CGEESX only computes a partial
          eigendecomposition, i.e. not Schur vectors
[out]WTMP
          WTMP is COMPLEX array, dimension (N)
          Like W, this array contains the eigenvalues of A,
          but sorted by increasing real or imaginary part.
[out]VS
          VS is COMPLEX array, dimension (LDVS, N)
          VS holds the computed Schur vectors.
[in]LDVS
          LDVS is INTEGER
          Leading dimension of VS. Must be at least max(1, N).
[out]VS1
          VS1 is COMPLEX array, dimension (LDVS, N)
          VS1 holds another copy of the computed Schur vectors.
[in]RCDEIN
          RCDEIN is REAL
          When COMP = .TRUE. RCDEIN holds the precomputed reciprocal
          condition number for the average of selected eigenvalues.
[in]RCDVIN
          RCDVIN is REAL
          When COMP = .TRUE. RCDVIN holds the precomputed reciprocal
          condition number for the selected right invariant subspace.
[in]NSLCT
          NSLCT is INTEGER
          When COMP = .TRUE. the number of selected eigenvalues
          corresponding to the precomputed values RCDEIN and RCDVIN.
[in]ISLCT
          ISLCT is INTEGER array, dimension (NSLCT)
          When COMP = .TRUE. ISLCT selects the eigenvalues of the
          input matrix corresponding to the precomputed values RCDEIN
          and RCDVIN. For I=1, ... ,NSLCT, if ISLCT(I) = J, then the
          eigenvalue with the J-th largest real or imaginary part is
          selected. The real part is used if ISRT = 0, and the
          imaginary part if ISRT = 1.
          Not referenced if COMP = .FALSE.
[in]ISRT
          ISRT is INTEGER
          When COMP = .TRUE., ISRT describes how ISLCT is used to
          choose a subset of the spectrum.
          Not referenced if COMP = .FALSE.
[out]RESULT
          RESULT is REAL array, dimension (17)
          The values computed by the 17 tests described above.
          The values are currently limited to 1/ulp, to avoid
          overflow.
[out]WORK
          WORK is COMPLEX array, dimension (2*N*N)
[in]LWORK
          LWORK is INTEGER
          The number of entries in WORK to be passed to CGEESX. This
          must be at least 2*N, and N*(N+1)/2 if tests 14--16 are to
          be performed.
[out]RWORK
          RWORK is REAL array, dimension (N)
[out]BWORK
          BWORK is LOGICAL array, dimension (N)
[out]INFO
          INFO is INTEGER
          If 0,  successful exit.
          If <0, input parameter -INFO had an incorrect value.
          If >0, CGEESX returned an error code, the absolute
                 value of which is returned.
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.

Definition at line 331 of file cget24.f.

335 *
336 * -- LAPACK test routine --
337 * -- LAPACK is a software package provided by Univ. of Tennessee, --
338 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
339 *
340 * .. Scalar Arguments ..
341  LOGICAL COMP
342  INTEGER INFO, ISRT, JTYPE, LDA, LDVS, LWORK, N, NOUNIT,
343  $ NSLCT
344  REAL RCDEIN, RCDVIN, THRESH
345 * ..
346 * .. Array Arguments ..
347  LOGICAL BWORK( * )
348  INTEGER ISEED( 4 ), ISLCT( * )
349  REAL RESULT( 17 ), RWORK( * )
350  COMPLEX A( LDA, * ), H( LDA, * ), HT( LDA, * ),
351  $ VS( LDVS, * ), VS1( LDVS, * ), W( * ),
352  $ WORK( * ), WT( * ), WTMP( * )
353 * ..
354 *
355 * =====================================================================
356 *
357 * .. Parameters ..
358  COMPLEX CZERO, CONE
359  parameter( czero = ( 0.0e+0, 0.0e+0 ),
360  $ cone = ( 1.0e+0, 0.0e+0 ) )
361  REAL ZERO, ONE
362  parameter( zero = 0.0e+0, one = 1.0e+0 )
363  REAL EPSIN
364  parameter( epsin = 5.9605e-8 )
365 * ..
366 * .. Local Scalars ..
367  CHARACTER SORT
368  INTEGER I, IINFO, ISORT, ITMP, J, KMIN, KNTEIG, RSUB,
369  $ SDIM, SDIM1
370  REAL ANORM, EPS, RCNDE1, RCNDV1, RCONDE, RCONDV,
371  $ SMLNUM, TOL, TOLIN, ULP, ULPINV, V, VRICMP,
372  $ VRIMIN, WNORM
373  COMPLEX CTMP
374 * ..
375 * .. Local Arrays ..
376  INTEGER IPNT( 20 )
377 * ..
378 * .. External Functions ..
379  LOGICAL CSLECT
380  REAL CLANGE, SLAMCH
381  EXTERNAL cslect, clange, slamch
382 * ..
383 * .. External Subroutines ..
384  EXTERNAL ccopy, cgeesx, cgemm, clacpy, cunt01, xerbla
385 * ..
386 * .. Intrinsic Functions ..
387  INTRINSIC abs, aimag, max, min, real
388 * ..
389 * .. Arrays in Common ..
390  LOGICAL SELVAL( 20 )
391  REAL SELWI( 20 ), SELWR( 20 )
392 * ..
393 * .. Scalars in Common ..
394  INTEGER SELDIM, SELOPT
395 * ..
396 * .. Common blocks ..
397  COMMON / sslct / selopt, seldim, selval, selwr, selwi
398 * ..
399 * .. Executable Statements ..
400 *
401 * Check for errors
402 *
403  info = 0
404  IF( thresh.LT.zero ) THEN
405  info = -3
406  ELSE IF( nounit.LE.0 ) THEN
407  info = -5
408  ELSE IF( n.LT.0 ) THEN
409  info = -6
410  ELSE IF( lda.LT.1 .OR. lda.LT.n ) THEN
411  info = -8
412  ELSE IF( ldvs.LT.1 .OR. ldvs.LT.n ) THEN
413  info = -15
414  ELSE IF( lwork.LT.2*n ) THEN
415  info = -24
416  END IF
417 *
418  IF( info.NE.0 ) THEN
419  CALL xerbla( 'CGET24', -info )
420  RETURN
421  END IF
422 *
423 * Quick return if nothing to do
424 *
425  DO 10 i = 1, 17
426  result( i ) = -one
427  10 CONTINUE
428 *
429  IF( n.EQ.0 )
430  $ RETURN
431 *
432 * Important constants
433 *
434  smlnum = slamch( 'Safe minimum' )
435  ulp = slamch( 'Precision' )
436  ulpinv = one / ulp
437 *
438 * Perform tests (1)-(13)
439 *
440  selopt = 0
441  DO 90 isort = 0, 1
442  IF( isort.EQ.0 ) THEN
443  sort = 'N'
444  rsub = 0
445  ELSE
446  sort = 'S'
447  rsub = 6
448  END IF
449 *
450 * Compute Schur form and Schur vectors, and test them
451 *
452  CALL clacpy( 'F', n, n, a, lda, h, lda )
453  CALL cgeesx( 'V', sort, cslect, 'N', n, h, lda, sdim, w, vs,
454  $ ldvs, rconde, rcondv, work, lwork, rwork, bwork,
455  $ iinfo )
456  IF( iinfo.NE.0 ) THEN
457  result( 1+rsub ) = ulpinv
458  IF( jtype.NE.22 ) THEN
459  WRITE( nounit, fmt = 9998 )'CGEESX1', iinfo, n, jtype,
460  $ iseed
461  ELSE
462  WRITE( nounit, fmt = 9999 )'CGEESX1', iinfo, n,
463  $ iseed( 1 )
464  END IF
465  info = abs( iinfo )
466  RETURN
467  END IF
468  IF( isort.EQ.0 ) THEN
469  CALL ccopy( n, w, 1, wtmp, 1 )
470  END IF
471 *
472 * Do Test (1) or Test (7)
473 *
474  result( 1+rsub ) = zero
475  DO 30 j = 1, n - 1
476  DO 20 i = j + 1, n
477  IF( h( i, j ).NE.czero )
478  $ result( 1+rsub ) = ulpinv
479  20 CONTINUE
480  30 CONTINUE
481 *
482 * Test (2) or (8): Compute norm(A - Q*H*Q') / (norm(A) * N * ULP)
483 *
484 * Copy A to VS1, used as workspace
485 *
486  CALL clacpy( ' ', n, n, a, lda, vs1, ldvs )
487 *
488 * Compute Q*H and store in HT.
489 *
490  CALL cgemm( 'No transpose', 'No transpose', n, n, n, cone, vs,
491  $ ldvs, h, lda, czero, ht, lda )
492 *
493 * Compute A - Q*H*Q'
494 *
495  CALL cgemm( 'No transpose', 'Conjugate transpose', n, n, n,
496  $ -cone, ht, lda, vs, ldvs, cone, vs1, ldvs )
497 *
498  anorm = max( clange( '1', n, n, a, lda, rwork ), smlnum )
499  wnorm = clange( '1', n, n, vs1, ldvs, rwork )
500 *
501  IF( anorm.GT.wnorm ) THEN
502  result( 2+rsub ) = ( wnorm / anorm ) / ( n*ulp )
503  ELSE
504  IF( anorm.LT.one ) THEN
505  result( 2+rsub ) = ( min( wnorm, n*anorm ) / anorm ) /
506  $ ( n*ulp )
507  ELSE
508  result( 2+rsub ) = min( wnorm / anorm, real( n ) ) /
509  $ ( n*ulp )
510  END IF
511  END IF
512 *
513 * Test (3) or (9): Compute norm( I - Q'*Q ) / ( N * ULP )
514 *
515  CALL cunt01( 'Columns', n, n, vs, ldvs, work, lwork, rwork,
516  $ result( 3+rsub ) )
517 *
518 * Do Test (4) or Test (10)
519 *
520  result( 4+rsub ) = zero
521  DO 40 i = 1, n
522  IF( h( i, i ).NE.w( i ) )
523  $ result( 4+rsub ) = ulpinv
524  40 CONTINUE
525 *
526 * Do Test (5) or Test (11)
527 *
528  CALL clacpy( 'F', n, n, a, lda, ht, lda )
529  CALL cgeesx( 'N', sort, cslect, 'N', n, ht, lda, sdim, wt, vs,
530  $ ldvs, rconde, rcondv, work, lwork, rwork, bwork,
531  $ iinfo )
532  IF( iinfo.NE.0 ) THEN
533  result( 5+rsub ) = ulpinv
534  IF( jtype.NE.22 ) THEN
535  WRITE( nounit, fmt = 9998 )'CGEESX2', iinfo, n, jtype,
536  $ iseed
537  ELSE
538  WRITE( nounit, fmt = 9999 )'CGEESX2', iinfo, n,
539  $ iseed( 1 )
540  END IF
541  info = abs( iinfo )
542  GO TO 220
543  END IF
544 *
545  result( 5+rsub ) = zero
546  DO 60 j = 1, n
547  DO 50 i = 1, n
548  IF( h( i, j ).NE.ht( i, j ) )
549  $ result( 5+rsub ) = ulpinv
550  50 CONTINUE
551  60 CONTINUE
552 *
553 * Do Test (6) or Test (12)
554 *
555  result( 6+rsub ) = zero
556  DO 70 i = 1, n
557  IF( w( i ).NE.wt( i ) )
558  $ result( 6+rsub ) = ulpinv
559  70 CONTINUE
560 *
561 * Do Test (13)
562 *
563  IF( isort.EQ.1 ) THEN
564  result( 13 ) = zero
565  knteig = 0
566  DO 80 i = 1, n
567  IF( cslect( w( i ) ) )
568  $ knteig = knteig + 1
569  IF( i.LT.n ) THEN
570  IF( cslect( w( i+1 ) ) .AND.
571  $ ( .NOT.cslect( w( i ) ) ) )result( 13 ) = ulpinv
572  END IF
573  80 CONTINUE
574  IF( sdim.NE.knteig )
575  $ result( 13 ) = ulpinv
576  END IF
577 *
578  90 CONTINUE
579 *
580 * If there is enough workspace, perform tests (14) and (15)
581 * as well as (10) through (13)
582 *
583  IF( lwork.GE.( n*( n+1 ) ) / 2 ) THEN
584 *
585 * Compute both RCONDE and RCONDV with VS
586 *
587  sort = 'S'
588  result( 14 ) = zero
589  result( 15 ) = zero
590  CALL clacpy( 'F', n, n, a, lda, ht, lda )
591  CALL cgeesx( 'V', sort, cslect, 'B', n, ht, lda, sdim1, wt,
592  $ vs1, ldvs, rconde, rcondv, work, lwork, rwork,
593  $ bwork, iinfo )
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 )'CGEESX3', iinfo, n, jtype,
599  $ iseed
600  ELSE
601  WRITE( nounit, fmt = 9999 )'CGEESX3', iinfo, n,
602  $ iseed( 1 )
603  END IF
604  info = abs( iinfo )
605  GO TO 220
606  END IF
607 *
608 * Perform tests (10), (11), (12), and (13)
609 *
610  DO 110 i = 1, n
611  IF( w( i ).NE.wt( i ) )
612  $ result( 10 ) = ulpinv
613  DO 100 j = 1, n
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
618  100 CONTINUE
619  110 CONTINUE
620  IF( sdim.NE.sdim1 )
621  $ result( 13 ) = ulpinv
622 *
623 * Compute both RCONDE and RCONDV without VS, and compare
624 *
625  CALL clacpy( 'F', n, n, a, lda, ht, lda )
626  CALL cgeesx( 'N', sort, cslect, 'B', n, ht, lda, sdim1, wt,
627  $ vs1, ldvs, rcnde1, rcndv1, work, lwork, rwork,
628  $ bwork, iinfo )
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 )'CGEESX4', iinfo, n, jtype,
634  $ iseed
635  ELSE
636  WRITE( nounit, fmt = 9999 )'CGEESX4', iinfo, n,
637  $ iseed( 1 )
638  END IF
639  info = abs( iinfo )
640  GO TO 220
641  END IF
642 *
643 * Perform tests (14) and (15)
644 *
645  IF( rcnde1.NE.rconde )
646  $ result( 14 ) = ulpinv
647  IF( rcndv1.NE.rcondv )
648  $ result( 15 ) = ulpinv
649 *
650 * Perform tests (10), (11), (12), and (13)
651 *
652  DO 130 i = 1, n
653  IF( w( i ).NE.wt( i ) )
654  $ result( 10 ) = ulpinv
655  DO 120 j = 1, n
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
660  120 CONTINUE
661  130 CONTINUE
662  IF( sdim.NE.sdim1 )
663  $ result( 13 ) = ulpinv
664 *
665 * Compute RCONDE with VS, and compare
666 *
667  CALL clacpy( 'F', n, n, a, lda, ht, lda )
668  CALL cgeesx( 'V', sort, cslect, 'E', n, ht, lda, sdim1, wt,
669  $ vs1, ldvs, rcnde1, rcndv1, work, lwork, rwork,
670  $ bwork, iinfo )
671  IF( iinfo.NE.0 ) THEN
672  result( 14 ) = ulpinv
673  IF( jtype.NE.22 ) THEN
674  WRITE( nounit, fmt = 9998 )'CGEESX5', iinfo, n, jtype,
675  $ iseed
676  ELSE
677  WRITE( nounit, fmt = 9999 )'CGEESX5', iinfo, n,
678  $ iseed( 1 )
679  END IF
680  info = abs( iinfo )
681  GO TO 220
682  END IF
683 *
684 * Perform test (14)
685 *
686  IF( rcnde1.NE.rconde )
687  $ result( 14 ) = ulpinv
688 *
689 * Perform tests (10), (11), (12), and (13)
690 *
691  DO 150 i = 1, n
692  IF( w( i ).NE.wt( i ) )
693  $ result( 10 ) = ulpinv
694  DO 140 j = 1, n
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
699  140 CONTINUE
700  150 CONTINUE
701  IF( sdim.NE.sdim1 )
702  $ result( 13 ) = ulpinv
703 *
704 * Compute RCONDE without VS, and compare
705 *
706  CALL clacpy( 'F', n, n, a, lda, ht, lda )
707  CALL cgeesx( 'N', sort, cslect, 'E', n, ht, lda, sdim1, wt,
708  $ vs1, ldvs, rcnde1, rcndv1, work, lwork, rwork,
709  $ bwork, iinfo )
710  IF( iinfo.NE.0 ) THEN
711  result( 14 ) = ulpinv
712  IF( jtype.NE.22 ) THEN
713  WRITE( nounit, fmt = 9998 )'CGEESX6', iinfo, n, jtype,
714  $ iseed
715  ELSE
716  WRITE( nounit, fmt = 9999 )'CGEESX6', iinfo, n,
717  $ iseed( 1 )
718  END IF
719  info = abs( iinfo )
720  GO TO 220
721  END IF
722 *
723 * Perform test (14)
724 *
725  IF( rcnde1.NE.rconde )
726  $ result( 14 ) = ulpinv
727 *
728 * Perform tests (10), (11), (12), and (13)
729 *
730  DO 170 i = 1, n
731  IF( w( i ).NE.wt( i ) )
732  $ result( 10 ) = ulpinv
733  DO 160 j = 1, n
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
738  160 CONTINUE
739  170 CONTINUE
740  IF( sdim.NE.sdim1 )
741  $ result( 13 ) = ulpinv
742 *
743 * Compute RCONDV with VS, and compare
744 *
745  CALL clacpy( 'F', n, n, a, lda, ht, lda )
746  CALL cgeesx( 'V', sort, cslect, 'V', n, ht, lda, sdim1, wt,
747  $ vs1, ldvs, rcnde1, rcndv1, work, lwork, rwork,
748  $ bwork, iinfo )
749  IF( iinfo.NE.0 ) THEN
750  result( 15 ) = ulpinv
751  IF( jtype.NE.22 ) THEN
752  WRITE( nounit, fmt = 9998 )'CGEESX7', iinfo, n, jtype,
753  $ iseed
754  ELSE
755  WRITE( nounit, fmt = 9999 )'CGEESX7', iinfo, n,
756  $ iseed( 1 )
757  END IF
758  info = abs( iinfo )
759  GO TO 220
760  END IF
761 *
762 * Perform test (15)
763 *
764  IF( rcndv1.NE.rcondv )
765  $ result( 15 ) = ulpinv
766 *
767 * Perform tests (10), (11), (12), and (13)
768 *
769  DO 190 i = 1, n
770  IF( w( i ).NE.wt( i ) )
771  $ result( 10 ) = ulpinv
772  DO 180 j = 1, n
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
777  180 CONTINUE
778  190 CONTINUE
779  IF( sdim.NE.sdim1 )
780  $ result( 13 ) = ulpinv
781 *
782 * Compute RCONDV without VS, and compare
783 *
784  CALL clacpy( 'F', n, n, a, lda, ht, lda )
785  CALL cgeesx( 'N', sort, cslect, 'V', n, ht, lda, sdim1, wt,
786  $ vs1, ldvs, rcnde1, rcndv1, work, lwork, rwork,
787  $ bwork, iinfo )
788  IF( iinfo.NE.0 ) THEN
789  result( 15 ) = ulpinv
790  IF( jtype.NE.22 ) THEN
791  WRITE( nounit, fmt = 9998 )'CGEESX8', iinfo, n, jtype,
792  $ iseed
793  ELSE
794  WRITE( nounit, fmt = 9999 )'CGEESX8', iinfo, n,
795  $ iseed( 1 )
796  END IF
797  info = abs( iinfo )
798  GO TO 220
799  END IF
800 *
801 * Perform test (15)
802 *
803  IF( rcndv1.NE.rcondv )
804  $ result( 15 ) = ulpinv
805 *
806 * Perform tests (10), (11), (12), and (13)
807 *
808  DO 210 i = 1, n
809  IF( w( i ).NE.wt( i ) )
810  $ result( 10 ) = ulpinv
811  DO 200 j = 1, n
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
816  200 CONTINUE
817  210 CONTINUE
818  IF( sdim.NE.sdim1 )
819  $ result( 13 ) = ulpinv
820 *
821  END IF
822 *
823  220 CONTINUE
824 *
825 * If there are precomputed reciprocal condition numbers, compare
826 * computed values with them.
827 *
828  IF( comp ) THEN
829 *
830 * First set up SELOPT, SELDIM, SELVAL, SELWR and SELWI so that
831 * the logical function CSLECT selects the eigenvalues specified
832 * by NSLCT, ISLCT and ISRT.
833 *
834  seldim = n
835  selopt = 1
836  eps = max( ulp, epsin )
837  DO 230 i = 1, n
838  ipnt( i ) = i
839  selval( i ) = .false.
840  selwr( i ) = real( wtmp( i ) )
841  selwi( i ) = aimag( wtmp( i ) )
842  230 CONTINUE
843  DO 250 i = 1, n - 1
844  kmin = i
845  IF( isrt.EQ.0 ) THEN
846  vrimin = real( wtmp( i ) )
847  ELSE
848  vrimin = aimag( wtmp( i ) )
849  END IF
850  DO 240 j = i + 1, n
851  IF( isrt.EQ.0 ) THEN
852  vricmp = real( wtmp( j ) )
853  ELSE
854  vricmp = aimag( wtmp( j ) )
855  END IF
856  IF( vricmp.LT.vrimin ) THEN
857  kmin = j
858  vrimin = vricmp
859  END IF
860  240 CONTINUE
861  ctmp = wtmp( kmin )
862  wtmp( kmin ) = wtmp( i )
863  wtmp( i ) = ctmp
864  itmp = ipnt( i )
865  ipnt( i ) = ipnt( kmin )
866  ipnt( kmin ) = itmp
867  250 CONTINUE
868  DO 260 i = 1, nslct
869  selval( ipnt( islct( i ) ) ) = .true.
870  260 CONTINUE
871 *
872 * Compute condition numbers
873 *
874  CALL clacpy( 'F', n, n, a, lda, ht, lda )
875  CALL cgeesx( 'N', 'S', cslect, 'B', n, ht, lda, sdim1, wt, vs1,
876  $ ldvs, rconde, rcondv, work, lwork, rwork, bwork,
877  $ iinfo )
878  IF( iinfo.NE.0 ) THEN
879  result( 16 ) = ulpinv
880  result( 17 ) = ulpinv
881  WRITE( nounit, fmt = 9999 )'CGEESX9', iinfo, n, iseed( 1 )
882  info = abs( iinfo )
883  GO TO 270
884  END IF
885 *
886 * Compare condition number for average of selected eigenvalues
887 * taking its condition number into account
888 *
889  anorm = clange( '1', n, n, a, lda, rwork )
890  v = max( real( n )*eps*anorm, smlnum )
891  IF( anorm.EQ.zero )
892  $ v = one
893  IF( v.GT.rcondv ) THEN
894  tol = one
895  ELSE
896  tol = v / rcondv
897  END IF
898  IF( v.GT.rcdvin ) THEN
899  tolin = one
900  ELSE
901  tolin = v / rcdvin
902  END IF
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 )
913  ELSE
914  result( 16 ) = one
915  END IF
916 *
917 * Compare condition numbers for right invariant subspace
918 * taking its condition number into account
919 *
920  IF( v.GT.rcondv*rconde ) THEN
921  tol = rcondv
922  ELSE
923  tol = v / rconde
924  END IF
925  IF( v.GT.rcdvin*rcdein ) THEN
926  tolin = rcdvin
927  ELSE
928  tolin = v / rcdein
929  END IF
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 )
940  ELSE
941  result( 17 ) = one
942  END IF
943 *
944  270 CONTINUE
945 *
946  END IF
947 *
948  9999 FORMAT( ' CGET24: ', a, ' returned INFO=', i6, '.', / 9x, 'N=',
949  $ i6, ', INPUT EXAMPLE NUMBER = ', i4 )
950  9998 FORMAT( ' CGET24: ', a, ' returned INFO=', i6, '.', / 9x, 'N=',
951  $ i6, ', JTYPE=', i6, ', ISEED=(', 3( i5, ',' ), i5, ')' )
952 *
953  RETURN
954 *
955 * End of CGET24
956 *
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:60
subroutine ccopy(N, CX, INCX, CY, INCY)
CCOPY
Definition: ccopy.f:81
subroutine cgemm(TRANSA, TRANSB, M, N, K, ALPHA, A, LDA, B, LDB, BETA, C, LDC)
CGEMM
Definition: cgemm.f:187
logical function cslect(Z)
CSLECT
Definition: cslect.f:56
subroutine cunt01(ROWCOL, M, N, U, LDU, WORK, LWORK, RWORK, RESID)
CUNT01
Definition: cunt01.f:126
real function clange(NORM, M, N, A, LDA, WORK)
CLANGE returns the value of the 1-norm, Frobenius norm, infinity-norm, or the largest absolute value ...
Definition: clange.f:115
subroutine cgeesx(JOBVS, SORT, SELECT, SENSE, N, A, LDA, SDIM, W, VS, LDVS, RCONDE, RCONDV, WORK, LWORK, RWORK, BWORK, INFO)
CGEESX computes the eigenvalues, the Schur form, and, optionally, the matrix of Schur vectors for GE ...
Definition: cgeesx.f:239
subroutine clacpy(UPLO, M, N, A, LDA, B, LDB)
CLACPY copies all or part of one two-dimensional array to another.
Definition: clacpy.f:103
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: