LAPACK  3.8.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.
Date
December 2016

Definition at line 337 of file cget24.f.

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