TESTING/EIG/cdrgev.f

   1 *> \brief \b CDRGEV
   2 *
   3 *  =========== DOCUMENTATION ===========
   4 *
   5 * Online html documentation available at
   6 *            http://www.netlib.org/lapack/explore-html/
   7 *
   8 *  Definition:
   9 *  ===========
  10 *
  11 *       SUBROUTINE CDRGEV( NSIZES, NN, NTYPES, DOTYPE, ISEED, THRESH,
  12 *                          NOUNIT, A, LDA, B, S, T, Q, LDQ, Z, QE, LDQE,
  13 *                          ALPHA, BETA, ALPHA1, BETA1, WORK, LWORK, RWORK,
  14 *                          RESULT, INFO )
  15 *
  16 *       .. Scalar Arguments ..
  17 *       INTEGER            INFO, LDA, LDQ, LDQE, LWORK, NOUNIT, NSIZES,
  18 *      $                   NTYPES
  19 *       REAL               THRESH
  20 *       ..
  21 *       .. Array Arguments ..
  22 *       LOGICAL            DOTYPE( * )
  23 *       INTEGER            ISEED( 4 ), NN( * )
  24 *       REAL               RESULT( * ), RWORK( * )
  25 *       COMPLEX            A( LDA, * ), ALPHA( * ), ALPHA1( * ),
  26 *      $                   B( LDA, * ), BETA( * ), BETA1( * ),
  27 *      $                   Q( LDQ, * ), QE( LDQE, * ), S( LDA, * ),
  28 *      $                   T( LDA, * ), WORK( * ), Z( LDQ, * )
  29 *       ..
  30 *
  31 *
  32 *> \par Purpose:
  33 *  =============
  34 *>
  35 *> \verbatim
  36 *>
  37 *> CDRGEV checks the nonsymmetric generalized eigenvalue problem driver
  38 *> routine CGGEV.
  39 *>
  40 *> CGGEV computes for a pair of n-by-n nonsymmetric matrices (A,B) the
  41 *> generalized eigenvalues and, optionally, the left and right
  42 *> eigenvectors.
  43 *>
  44 *> A generalized eigenvalue for a pair of matrices (A,B) is a scalar w
  45 *> or a ratio  alpha/beta = w, such that A - w*B is singular.  It is
  46 *> usually represented as the pair (alpha,beta), as there is reasonable
  47 *> interpretation for beta=0, and even for both being zero.
  48 *>
  49 *> A right generalized eigenvector corresponding to a generalized
  50 *> eigenvalue  w  for a pair of matrices (A,B) is a vector r  such that
  51 *> (A - wB) * r = 0.  A left generalized eigenvector is a vector l such
  52 *> that l**H * (A - wB) = 0, where l**H is the conjugate-transpose of l.
  53 *>
  54 *> When CDRGEV is called, a number of matrix "sizes" ("n's") and a
  55 *> number of matrix "types" are specified.  For each size ("n")
  56 *> and each type of matrix, a pair of matrices (A, B) will be generated
  57 *> and used for testing.  For each matrix pair, the following tests
  58 *> will be performed and compared with the threshold THRESH.
  59 *>
  60 *> Results from CGGEV:
  61 *>
  62 *> (1)  max over all left eigenvalue/-vector pairs (alpha/beta,l) of
  63 *>
  64 *>      | VL**H * (beta A - alpha B) |/( ulp max(|beta A|, |alpha B|) )
  65 *>
  66 *>      where VL**H is the conjugate-transpose of VL.
  67 *>
  68 *> (2)  | |VL(i)| - 1 | / ulp and whether largest component real
  69 *>
  70 *>      VL(i) denotes the i-th column of VL.
  71 *>
  72 *> (3)  max over all left eigenvalue/-vector pairs (alpha/beta,r) of
  73 *>
  74 *>      | (beta A - alpha B) * VR | / ( ulp max(|beta A|, |alpha B|) )
  75 *>
  76 *> (4)  | |VR(i)| - 1 | / ulp and whether largest component real
  77 *>
  78 *>      VR(i) denotes the i-th column of VR.
  79 *>
  80 *> (5)  W(full) = W(partial)
  81 *>      W(full) denotes the eigenvalues computed when both l and r
  82 *>      are also computed, and W(partial) denotes the eigenvalues
  83 *>      computed when only W, only W and r, or only W and l are
  84 *>      computed.
  85 *>
  86 *> (6)  VL(full) = VL(partial)
  87 *>      VL(full) denotes the left eigenvectors computed when both l
  88 *>      and r are computed, and VL(partial) denotes the result
  89 *>      when only l is computed.
  90 *>
  91 *> (7)  VR(full) = VR(partial)
  92 *>      VR(full) denotes the right eigenvectors computed when both l
  93 *>      and r are also computed, and VR(partial) denotes the result
  94 *>      when only l is computed.
  95 *>
  96 *>
  97 *> Test Matrices
  98 *> ---- --------
  99 *>
 100 *> The sizes of the test matrices are specified by an array
 101 *> NN(1:NSIZES); the value of each element NN(j) specifies one size.
 102 *> The "types" are specified by a logical array DOTYPE( 1:NTYPES ); if
 103 *> DOTYPE(j) is .TRUE., then matrix type "j" will be generated.
 104 *> Currently, the list of possible types is:
 105 *>
 106 *> (1)  ( 0, 0 )         (a pair of zero matrices)
 107 *>
 108 *> (2)  ( I, 0 )         (an identity and a zero matrix)
 109 *>
 110 *> (3)  ( 0, I )         (an identity and a zero matrix)
 111 *>
 112 *> (4)  ( I, I )         (a pair of identity matrices)
 113 *>
 114 *>         t   t
 115 *> (5)  ( J , J  )       (a pair of transposed Jordan blocks)
 116 *>
 117 *>                                     t                ( I   0  )
 118 *> (6)  ( X, Y )         where  X = ( J   0  )  and Y = (      t )
 119 *>                                  ( 0   I  )          ( 0   J  )
 120 *>                       and I is a k x k identity and J a (k+1)x(k+1)
 121 *>                       Jordan block; k=(N-1)/2
 122 *>
 123 *> (7)  ( D, I )         where D is diag( 0, 1,..., N-1 ) (a diagonal
 124 *>                       matrix with those diagonal entries.)
 125 *> (8)  ( I, D )
 126 *>
 127 *> (9)  ( big*D, small*I ) where "big" is near overflow and small=1/big
 128 *>
 129 *> (10) ( small*D, big*I )
 130 *>
 131 *> (11) ( big*I, small*D )
 132 *>
 133 *> (12) ( small*I, big*D )
 134 *>
 135 *> (13) ( big*D, big*I )
 136 *>
 137 *> (14) ( small*D, small*I )
 138 *>
 139 *> (15) ( D1, D2 )        where D1 is diag( 0, 0, 1, ..., N-3, 0 ) and
 140 *>                        D2 is diag( 0, N-3, N-4,..., 1, 0, 0 )
 141 *>           t   t
 142 *> (16) Q ( J , J ) Z     where Q and Z are random orthogonal matrices.
 143 *>
 144 *> (17) Q ( T1, T2 ) Z    where T1 and T2 are upper triangular matrices
 145 *>                        with random O(1) entries above the diagonal
 146 *>                        and diagonal entries diag(T1) =
 147 *>                        ( 0, 0, 1, ..., N-3, 0 ) and diag(T2) =
 148 *>                        ( 0, N-3, N-4,..., 1, 0, 0 )
 149 *>
 150 *> (18) Q ( T1, T2 ) Z    diag(T1) = ( 0, 0, 1, 1, s, ..., s, 0 )
 151 *>                        diag(T2) = ( 0, 1, 0, 1,..., 1, 0 )
 152 *>                        s = machine precision.
 153 *>
 154 *> (19) Q ( T1, T2 ) Z    diag(T1)=( 0,0,1,1, 1-d, ..., 1-(N-5)*d=s, 0 )
 155 *>                        diag(T2) = ( 0, 1, 0, 1, ..., 1, 0 )
 156 *>
 157 *>                                                        N-5
 158 *> (20) Q ( T1, T2 ) Z    diag(T1)=( 0, 0, 1, 1, a, ..., a   =s, 0 )
 159 *>                        diag(T2) = ( 0, 1, 0, 1, ..., 1, 0, 0 )
 160 *>
 161 *> (21) Q ( T1, T2 ) Z    diag(T1)=( 0, 0, 1, r1, r2, ..., r(N-4), 0 )
 162 *>                        diag(T2) = ( 0, 1, 0, 1, ..., 1, 0, 0 )
 163 *>                        where r1,..., r(N-4) are random.
 164 *>
 165 *> (22) Q ( big*T1, small*T2 ) Z    diag(T1) = ( 0, 0, 1, ..., N-3, 0 )
 166 *>                                  diag(T2) = ( 0, 1, ..., 1, 0, 0 )
 167 *>
 168 *> (23) Q ( small*T1, big*T2 ) Z    diag(T1) = ( 0, 0, 1, ..., N-3, 0 )
 169 *>                                  diag(T2) = ( 0, 1, ..., 1, 0, 0 )
 170 *>
 171 *> (24) Q ( small*T1, small*T2 ) Z  diag(T1) = ( 0, 0, 1, ..., N-3, 0 )
 172 *>                                  diag(T2) = ( 0, 1, ..., 1, 0, 0 )
 173 *>
 174 *> (25) Q ( big*T1, big*T2 ) Z      diag(T1) = ( 0, 0, 1, ..., N-3, 0 )
 175 *>                                  diag(T2) = ( 0, 1, ..., 1, 0, 0 )
 176 *>
 177 *> (26) Q ( T1, T2 ) Z     where T1 and T2 are random upper-triangular
 178 *>                         matrices.
 179 *>
 180 *> \endverbatim
 181 *
 182 *  Arguments:
 183 *  ==========
 184 *
 185 *> \param[in] NSIZES
 186 *> \verbatim
 187 *>          NSIZES is INTEGER
 188 *>          The number of sizes of matrices to use.  If it is zero,
 189 *>          CDRGES does nothing.  NSIZES >= 0.
 190 *> \endverbatim
 191 *>
 192 *> \param[in] NN
 193 *> \verbatim
 194 *>          NN is INTEGER array, dimension (NSIZES)
 195 *>          An array containing the sizes to be used for the matrices.
 196 *>          Zero values will be skipped.  NN >= 0.
 197 *> \endverbatim
 198 *>
 199 *> \param[in] NTYPES
 200 *> \verbatim
 201 *>          NTYPES is INTEGER
 202 *>          The number of elements in DOTYPE.   If it is zero, CDRGEV
 203 *>          does nothing.  It must be at least zero.  If it is MAXTYP+1
 204 *>          and NSIZES is 1, then an additional type, MAXTYP+1 is
 205 *>          defined, which is to use whatever matrix is in A.  This
 206 *>          is only useful if DOTYPE(1:MAXTYP) is .FALSE. and
 207 *>          DOTYPE(MAXTYP+1) is .TRUE. .
 208 *> \endverbatim
 209 *>
 210 *> \param[in] DOTYPE
 211 *> \verbatim
 212 *>          DOTYPE is LOGICAL array, dimension (NTYPES)
 213 *>          If DOTYPE(j) is .TRUE., then for each size in NN a
 214 *>          matrix of that size and of type j will be generated.
 215 *>          If NTYPES is smaller than the maximum number of types
 216 *>          defined (PARAMETER MAXTYP), then types NTYPES+1 through
 217 *>          MAXTYP will not be generated. If NTYPES is larger
 218 *>          than MAXTYP, DOTYPE(MAXTYP+1) through DOTYPE(NTYPES)
 219 *>          will be ignored.
 220 *> \endverbatim
 221 *>
 222 *> \param[in,out] ISEED
 223 *> \verbatim
 224 *>          ISEED is INTEGER array, dimension (4)
 225 *>          On entry ISEED specifies the seed of the random number
 226 *>          generator. The array elements should be between 0 and 4095;
 227 *>          if not they will be reduced mod 4096. Also, ISEED(4) must
 228 *>          be odd.  The random number generator uses a linear
 229 *>          congruential sequence limited to small integers, and so
 230 *>          should produce machine independent random numbers. The
 231 *>          values of ISEED are changed on exit, and can be used in the
 232 *>          next call to CDRGES to continue the same random number
 233 *>          sequence.
 234 *> \endverbatim
 235 *>
 236 *> \param[in] THRESH
 237 *> \verbatim
 238 *>          THRESH is REAL
 239 *>          A test will count as "failed" if the "error", computed as
 240 *>          described above, exceeds THRESH.  Note that the error is
 241 *>          scaled to be O(1), so THRESH should be a reasonably small
 242 *>          multiple of 1, e.g., 10 or 100.  In particular, it should
 243 *>          not depend on the precision (single vs. double) or the size
 244 *>          of the matrix.  It must be at least zero.
 245 *> \endverbatim
 246 *>
 247 *> \param[in] NOUNIT
 248 *> \verbatim
 249 *>          NOUNIT is INTEGER
 250 *>          The FORTRAN unit number for printing out error messages
 251 *>          (e.g., if a routine returns IERR not equal to 0.)
 252 *> \endverbatim
 253 *>
 254 *> \param[in,out] A
 255 *> \verbatim
 256 *>          A is COMPLEX array, dimension(LDA, max(NN))
 257 *>          Used to hold the original A matrix.  Used as input only
 258 *>          if NTYPES=MAXTYP+1, DOTYPE(1:MAXTYP)=.FALSE., and
 259 *>          DOTYPE(MAXTYP+1)=.TRUE.
 260 *> \endverbatim
 261 *>
 262 *> \param[in] LDA
 263 *> \verbatim
 264 *>          LDA is INTEGER
 265 *>          The leading dimension of A, B, S, and T.
 266 *>          It must be at least 1 and at least max( NN ).
 267 *> \endverbatim
 268 *>
 269 *> \param[in,out] B
 270 *> \verbatim
 271 *>          B is COMPLEX array, dimension(LDA, max(NN))
 272 *>          Used to hold the original B matrix.  Used as input only
 273 *>          if NTYPES=MAXTYP+1, DOTYPE(1:MAXTYP)=.FALSE., and
 274 *>          DOTYPE(MAXTYP+1)=.TRUE.
 275 *> \endverbatim
 276 *>
 277 *> \param[out] S
 278 *> \verbatim
 279 *>          S is COMPLEX array, dimension (LDA, max(NN))
 280 *>          The Schur form matrix computed from A by CGGEV.  On exit, S
 281 *>          contains the Schur form matrix corresponding to the matrix
 282 *>          in A.
 283 *> \endverbatim
 284 *>
 285 *> \param[out] T
 286 *> \verbatim
 287 *>          T is COMPLEX array, dimension (LDA, max(NN))
 288 *>          The upper triangular matrix computed from B by CGGEV.
 289 *> \endverbatim
 290 *>
 291 *> \param[out] Q
 292 *> \verbatim
 293 *>          Q is COMPLEX array, dimension (LDQ, max(NN))
 294 *>          The (left) eigenvectors matrix computed by CGGEV.
 295 *> \endverbatim
 296 *>
 297 *> \param[in] LDQ
 298 *> \verbatim
 299 *>          LDQ is INTEGER
 300 *>          The leading dimension of Q and Z. It must
 301 *>          be at least 1 and at least max( NN ).
 302 *> \endverbatim
 303 *>
 304 *> \param[out] Z
 305 *> \verbatim
 306 *>          Z is COMPLEX array, dimension( LDQ, max(NN) )
 307 *>          The (right) orthogonal matrix computed by CGGEV.
 308 *> \endverbatim
 309 *>
 310 *> \param[out] QE
 311 *> \verbatim
 312 *>          QE is COMPLEX array, dimension( LDQ, max(NN) )
 313 *>          QE holds the computed right or left eigenvectors.
 314 *> \endverbatim
 315 *>
 316 *> \param[in] LDQE
 317 *> \verbatim
 318 *>          LDQE is INTEGER
 319 *>          The leading dimension of QE. LDQE >= max(1,max(NN)).
 320 *> \endverbatim
 321 *>
 322 *> \param[out] ALPHA
 323 *> \verbatim
 324 *>          ALPHA is COMPLEX array, dimension (max(NN))
 325 *> \endverbatim
 326 *>
 327 *> \param[out] BETA
 328 *> \verbatim
 329 *>          BETA is COMPLEX array, dimension (max(NN))
 330 *>
 331 *>          The generalized eigenvalues of (A,B) computed by CGGEV.
 332 *>          ( ALPHAR(k)+ALPHAI(k)*i ) / BETA(k) is the k-th
 333 *>          generalized eigenvalue of A and B.
 334 *> \endverbatim
 335 *>
 336 *> \param[out] ALPHA1
 337 *> \verbatim
 338 *>          ALPHA1 is COMPLEX array, dimension (max(NN))
 339 *> \endverbatim
 340 *>
 341 *> \param[out] BETA1
 342 *> \verbatim
 343 *>          BETA1 is COMPLEX array, dimension (max(NN))
 344 *>
 345 *>          Like ALPHAR, ALPHAI, BETA, these arrays contain the
 346 *>          eigenvalues of A and B, but those computed when CGGEV only
 347 *>          computes a partial eigendecomposition, i.e. not the
 348 *>          eigenvalues and left and right eigenvectors.
 349 *> \endverbatim
 350 *>
 351 *> \param[out] WORK
 352 *> \verbatim
 353 *>          WORK is COMPLEX array, dimension (LWORK)
 354 *> \endverbatim
 355 *>
 356 *> \param[in] LWORK
 357 *> \verbatim
 358 *>          LWORK is INTEGER
 359 *>          The number of entries in WORK.  LWORK >= N*(N+1)
 360 *> \endverbatim
 361 *>
 362 *> \param[out] RWORK
 363 *> \verbatim
 364 *>          RWORK is REAL array, dimension (8*N)
 365 *>          Real workspace.
 366 *> \endverbatim
 367 *>
 368 *> \param[out] RESULT
 369 *> \verbatim
 370 *>          RESULT is REAL array, dimension (2)
 371 *>          The values computed by the tests described above.
 372 *>          The values are currently limited to 1/ulp, to avoid overflow.
 373 *> \endverbatim
 374 *>
 375 *> \param[out] INFO
 376 *> \verbatim
 377 *>          INFO is INTEGER
 378 *>          = 0:  successful exit
 379 *>          < 0:  if INFO = -i, the i-th argument had an illegal value.
 380 *>          > 0:  A routine returned an error code.  INFO is the
 381 *>                absolute value of the INFO value returned.
 382 *> \endverbatim
 383 *
 384 *  Authors:
 385 *  ========
 386 *
 387 *> \author Univ. of Tennessee
 388 *> \author Univ. of California Berkeley
 389 *> \author Univ. of Colorado Denver
 390 *> \author NAG Ltd.
 391 *
 392 *> \date June 2016
 393 *
 394 *> \ingroup complex_eig
 395 *
 396 *  =====================================================================
 397       SUBROUTINE CDRGEV( NSIZES, NN, NTYPES, DOTYPE, ISEED, THRESH,
 398      $                   NOUNIT, A, LDA, B, S, T, Q, LDQ, Z, QE, LDQE,
 399      $                   ALPHA, BETA, ALPHA1, BETA1, WORK, LWORK, RWORK,
 400      $                   RESULT, INFO )
 401 *
 402 *  -- LAPACK test routine (version 3.6.1) --
 403 *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
 404 *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
 405 *     June 2016
 406 *
 407 *     .. Scalar Arguments ..
 408       INTEGER            INFO, LDA, LDQ, LDQE, LWORK, NOUNIT, NSIZES,
 409      $                   NTYPES
 410       REAL               THRESH
 411 *     ..
 412 *     .. Array Arguments ..
 413       LOGICAL            DOTYPE( * )
 414       INTEGER            ISEED( 4 ), NN( * )
 415       REAL               RESULT( * ), RWORK( * )
 416       COMPLEX            A( LDA, * ), ALPHA( * ), ALPHA1( * ),
 417      $                   B( LDA, * ), BETA( * ), BETA1( * ),
 418      $                   Q( LDQ, * ), QE( LDQE, * ), S( LDA, * ),
 419      $                   T( LDA, * ), WORK( * ), Z( LDQ, * )
 420 *     ..
 421 *
 422 *  =====================================================================
 423 *
 424 *     .. Parameters ..
 425       REAL               ZERO, ONE
 426       PARAMETER          ( ZERO = 0.0E+0, ONE = 1.0E+0 )
 427       COMPLEX            CZERO, CONE
 428       PARAMETER          ( CZERO = ( 0.0E+0, 0.0E+0 ),
 429      $                   CONE = ( 1.0E+0, 0.0E+0 ) )
 430       INTEGER            MAXTYP
 431       PARAMETER          ( MAXTYP = 26 )
 432 *     ..
 433 *     .. Local Scalars ..
 434       LOGICAL            BADNN
 435       INTEGER            I, IADD, IERR, IN, J, JC, JR, JSIZE, JTYPE,
 436      $                   MAXWRK, MINWRK, MTYPES, N, N1, NB, NERRS,
 437      $                   NMATS, NMAX, NTESTT
 438       REAL               SAFMAX, SAFMIN, ULP, ULPINV
 439       COMPLEX            CTEMP
 440 *     ..
 441 *     .. Local Arrays ..
 442       LOGICAL            LASIGN( MAXTYP ), LBSIGN( MAXTYP )
 443       INTEGER            IOLDSD( 4 ), KADD( 6 ), KAMAGN( MAXTYP ),
 444      $                   KATYPE( MAXTYP ), KAZERO( MAXTYP ),
 445      $                   KBMAGN( MAXTYP ), KBTYPE( MAXTYP ),
 446      $                   KBZERO( MAXTYP ), KCLASS( MAXTYP ),
 447      $                   KTRIAN( MAXTYP ), KZ1( 6 ), KZ2( 6 )
 448       REAL               RMAGN( 0: 3 )
 449 *     ..
 450 *     .. External Functions ..
 451       INTEGER            ILAENV
 452       REAL               SLAMCH
 453       COMPLEX            CLARND
 454       EXTERNAL           ILAENV, SLAMCH, CLARND
 455 *     ..
 456 *     .. External Subroutines ..
 457       EXTERNAL           ALASVM, CGET52, CGGEV, CLACPY, CLARFG, CLASET,
 458      $                   CLATM4, CUNM2R, SLABAD, XERBLA
 459 *     ..
 460 *     .. Intrinsic Functions ..
 461       INTRINSIC          ABS, CONJG, MAX, MIN, REAL, SIGN
 462 *     ..
 463 *     .. Data statements ..
 464       DATA               KCLASS / 15*1, 10*2, 1*3 /
 465       DATA               KZ1 / 0, 1, 2, 1, 3, 3 /
 466       DATA               KZ2 / 0, 0, 1, 2, 1, 1 /
 467       DATA               KADD / 0, 0, 0, 0, 3, 2 /
 468       DATA               KATYPE / 0, 1, 0, 1, 2, 3, 4, 1, 4, 4, 1, 1, 4,
 469      $                   4, 4, 2, 4, 5, 8, 7, 9, 4*4, 0 /
 470       DATA               KBTYPE / 0, 0, 1, 1, 2, -3, 1, 4, 1, 1, 4, 4,
 471      $                   1, 1, -4, 2, -4, 8*8, 0 /
 472       DATA               KAZERO / 6*1, 2, 1, 2*2, 2*1, 2*2, 3, 1, 3,
 473      $                   4*5, 4*3, 1 /
 474       DATA               KBZERO / 6*1, 1, 2, 2*1, 2*2, 2*1, 4, 1, 4,
 475      $                   4*6, 4*4, 1 /
 476       DATA               KAMAGN / 8*1, 2, 3, 2, 3, 2, 3, 7*1, 2, 3, 3,
 477      $                   2, 1 /
 478       DATA               KBMAGN / 8*1, 3, 2, 3, 2, 2, 3, 7*1, 3, 2, 3,
 479      $                   2, 1 /
 480       DATA               KTRIAN / 16*0, 10*1 /
 481       DATA               LASIGN / 6*.FALSE., .TRUE., .FALSE., 2*.TRUE.,
 482      $                   2*.FALSE., 3*.TRUE., .FALSE., .TRUE.,
 483      $                   3*.FALSE., 5*.TRUE., .FALSE. /
 484       DATA               LBSIGN / 7*.FALSE., .TRUE., 2*.FALSE.,
 485      $                   2*.TRUE., 2*.FALSE., .TRUE., .FALSE., .TRUE.,
 486      $                   9*.FALSE. /
 487 *     ..
 488 *     .. Executable Statements ..
 489 *
 490 *     Check for errors
 491 *
 492       INFO = 0
 493 *
 494       BADNN = .FALSE.
 495       NMAX = 1
 496       DO 10 J = 1, NSIZES
 497          NMAX = MAX( NMAX, NN( J ) )
 498          IF( NN( J ).LT.0 )
 499      $      BADNN = .TRUE.
 500    10 CONTINUE
 501 *
 502       IF( NSIZES.LT.0 ) THEN
 503          INFO = -1
 504       ELSE IF( BADNN ) THEN
 505          INFO = -2
 506       ELSE IF( NTYPES.LT.0 ) THEN
 507          INFO = -3
 508       ELSE IF( THRESH.LT.ZERO ) THEN
 509          INFO = -6
 510       ELSE IF( LDA.LE.1 .OR. LDA.LT.NMAX ) THEN
 511          INFO = -9
 512       ELSE IF( LDQ.LE.1 .OR. LDQ.LT.NMAX ) THEN
 513          INFO = -14
 514       ELSE IF( LDQE.LE.1 .OR. LDQE.LT.NMAX ) THEN
 515          INFO = -17
 516       END IF
 517 *
 518 *     Compute workspace
 519 *      (Note: Comments in the code beginning "Workspace:" describe the
 520 *       minimal amount of workspace needed at that point in the code,
 521 *       as well as the preferred amount for good performance.
 522 *       NB refers to the optimal block size for the immediately
 523 *       following subroutine, as returned by ILAENV.
 524 *
 525       MINWRK = 1
 526       IF( INFO.EQ.0 .AND. LWORK.GE.1 ) THEN
 527          MINWRK = NMAX*( NMAX+1 )
 528          NB = MAX( 1, ILAENV( 1, 'CGEQRF', ' ', NMAX, NMAX, -1, -1 ),
 529      $        ILAENV( 1, 'CUNMQR', 'LC', NMAX, NMAX, NMAX, -1 ),
 530      $        ILAENV( 1, 'CUNGQR', ' ', NMAX, NMAX, NMAX, -1 ) )
 531          MAXWRK = MAX( 2*NMAX, NMAX*( NB+1 ), NMAX*( NMAX+1 ) )
 532          WORK( 1 ) = MAXWRK
 533       END IF
 534 *
 535       IF( LWORK.LT.MINWRK )
 536      $   INFO = -23
 537 *
 538       IF( INFO.NE.0 ) THEN
 539          CALL XERBLA( 'CDRGEV', -INFO )
 540          RETURN
 541       END IF
 542 *
 543 *     Quick return if possible
 544 *
 545       IF( NSIZES.EQ.0 .OR. NTYPES.EQ.0 )
 546      $   RETURN
 547 *
 548       ULP = SLAMCH( 'Precision' )
 549       SAFMIN = SLAMCH( 'Safe minimum' )
 550       SAFMIN = SAFMIN / ULP
 551       SAFMAX = ONE / SAFMIN
 552       CALL SLABAD( SAFMIN, SAFMAX )
 553       ULPINV = ONE / ULP
 554 *
 555 *     The values RMAGN(2:3) depend on N, see below.
 556 *
 557       RMAGN( 0 ) = ZERO
 558       RMAGN( 1 ) = ONE
 559 *
 560 *     Loop over sizes, types
 561 *
 562       NTESTT = 0
 563       NERRS = 0
 564       NMATS = 0
 565 *
 566       DO 220 JSIZE = 1, NSIZES
 567          N = NN( JSIZE )
 568          N1 = MAX( 1, N )
 569          RMAGN( 2 ) = SAFMAX*ULP / REAL( N1 )
 570          RMAGN( 3 ) = SAFMIN*ULPINV*N1
 571 *
 572          IF( NSIZES.NE.1 ) THEN
 573             MTYPES = MIN( MAXTYP, NTYPES )
 574          ELSE
 575             MTYPES = MIN( MAXTYP+1, NTYPES )
 576          END IF
 577 *
 578          DO 210 JTYPE = 1, MTYPES
 579             IF( .NOT.DOTYPE( JTYPE ) )
 580      $         GO TO 210
 581             NMATS = NMATS + 1
 582 *
 583 *           Save ISEED in case of an error.
 584 *
 585             DO 20 J = 1, 4
 586                IOLDSD( J ) = ISEED( J )
 587    20       CONTINUE
 588 *
 589 *           Generate test matrices A and B
 590 *
 591 *           Description of control parameters:
 592 *
 593 *           KCLASS: =1 means w/o rotation, =2 means w/ rotation,
 594 *                   =3 means random.
 595 *           KATYPE: the "type" to be passed to CLATM4 for computing A.
 596 *           KAZERO: the pattern of zeros on the diagonal for A:
 597 *                   =1: ( xxx ), =2: (0, xxx ) =3: ( 0, 0, xxx, 0 ),
 598 *                   =4: ( 0, xxx, 0, 0 ), =5: ( 0, 0, 1, xxx, 0 ),
 599 *                   =6: ( 0, 1, 0, xxx, 0 ).  (xxx means a string of
 600 *                   non-zero entries.)
 601 *           KAMAGN: the magnitude of the matrix: =0: zero, =1: O(1),
 602 *                   =2: large, =3: small.
 603 *           LASIGN: .TRUE. if the diagonal elements of A are to be
 604 *                   multiplied by a random magnitude 1 number.
 605 *           KBTYPE, KBZERO, KBMAGN, LBSIGN: the same, but for B.
 606 *           KTRIAN: =0: don't fill in the upper triangle, =1: do.
 607 *           KZ1, KZ2, KADD: used to implement KAZERO and KBZERO.
 608 *           RMAGN: used to implement KAMAGN and KBMAGN.
 609 *
 610             IF( MTYPES.GT.MAXTYP )
 611      $         GO TO 100
 612             IERR = 0
 613             IF( KCLASS( JTYPE ).LT.3 ) THEN
 614 *
 615 *              Generate A (w/o rotation)
 616 *
 617                IF( ABS( KATYPE( JTYPE ) ).EQ.3 ) THEN
 618                   IN = 2*( ( N-1 ) / 2 ) + 1
 619                   IF( IN.NE.N )
 620      $               CALL CLASET( 'Full', N, N, CZERO, CZERO, A, LDA )
 621                ELSE
 622                   IN = N
 623                END IF
 624                CALL CLATM4( KATYPE( JTYPE ), IN, KZ1( KAZERO( JTYPE ) ),
 625      $                      KZ2( KAZERO( JTYPE ) ), LASIGN( JTYPE ),
 626      $                      RMAGN( KAMAGN( JTYPE ) ), ULP,
 627      $                      RMAGN( KTRIAN( JTYPE )*KAMAGN( JTYPE ) ), 2,
 628      $                      ISEED, A, LDA )
 629                IADD = KADD( KAZERO( JTYPE ) )
 630                IF( IADD.GT.0 .AND. IADD.LE.N )
 631      $            A( IADD, IADD ) = RMAGN( KAMAGN( JTYPE ) )
 632 *
 633 *              Generate B (w/o rotation)
 634 *
 635                IF( ABS( KBTYPE( JTYPE ) ).EQ.3 ) THEN
 636                   IN = 2*( ( N-1 ) / 2 ) + 1
 637                   IF( IN.NE.N )
 638      $               CALL CLASET( 'Full', N, N, CZERO, CZERO, B, LDA )
 639                ELSE
 640                   IN = N
 641                END IF
 642                CALL CLATM4( KBTYPE( JTYPE ), IN, KZ1( KBZERO( JTYPE ) ),
 643      $                      KZ2( KBZERO( JTYPE ) ), LBSIGN( JTYPE ),
 644      $                      RMAGN( KBMAGN( JTYPE ) ), ONE,
 645      $                      RMAGN( KTRIAN( JTYPE )*KBMAGN( JTYPE ) ), 2,
 646      $                      ISEED, B, LDA )
 647                IADD = KADD( KBZERO( JTYPE ) )
 648                IF( IADD.NE.0 .AND. IADD.LE.N )
 649      $            B( IADD, IADD ) = RMAGN( KBMAGN( JTYPE ) )
 650 *
 651                IF( KCLASS( JTYPE ).EQ.2 .AND. N.GT.0 ) THEN
 652 *
 653 *                 Include rotations
 654 *
 655 *                 Generate Q, Z as Householder transformations times
 656 *                 a diagonal matrix.
 657 *
 658                   DO 40 JC = 1, N - 1
 659                      DO 30 JR = JC, N
 660                         Q( JR, JC ) = CLARND( 3, ISEED )
 661                         Z( JR, JC ) = CLARND( 3, ISEED )
 662    30                CONTINUE
 663                      CALL CLARFG( N+1-JC, Q( JC, JC ), Q( JC+1, JC ), 1,
 664      $                            WORK( JC ) )
 665                      WORK( 2*N+JC ) = SIGN( ONE, REAL( Q( JC, JC ) ) )
 666                      Q( JC, JC ) = CONE
 667                      CALL CLARFG( N+1-JC, Z( JC, JC ), Z( JC+1, JC ), 1,
 668      $                            WORK( N+JC ) )
 669                      WORK( 3*N+JC ) = SIGN( ONE, REAL( Z( JC, JC ) ) )
 670                      Z( JC, JC ) = CONE
 671    40             CONTINUE
 672                   CTEMP = CLARND( 3, ISEED )
 673                   Q( N, N ) = CONE
 674                   WORK( N ) = CZERO
 675                   WORK( 3*N ) = CTEMP / ABS( CTEMP )
 676                   CTEMP = CLARND( 3, ISEED )
 677                   Z( N, N ) = CONE
 678                   WORK( 2*N ) = CZERO
 679                   WORK( 4*N ) = CTEMP / ABS( CTEMP )
 680 *
 681 *                 Apply the diagonal matrices
 682 *
 683                   DO 60 JC = 1, N
 684                      DO 50 JR = 1, N
 685                         A( JR, JC ) = WORK( 2*N+JR )*
 686      $                                CONJG( WORK( 3*N+JC ) )*
 687      $                                A( JR, JC )
 688                         B( JR, JC ) = WORK( 2*N+JR )*
 689      $                                CONJG( WORK( 3*N+JC ) )*
 690      $                                B( JR, JC )
 691    50                CONTINUE
 692    60             CONTINUE
 693                   CALL CUNM2R( 'L', 'N', N, N, N-1, Q, LDQ, WORK, A,
 694      $                         LDA, WORK( 2*N+1 ), IERR )
 695                   IF( IERR.NE.0 )
 696      $               GO TO 90
 697                   CALL CUNM2R( 'R', 'C', N, N, N-1, Z, LDQ, WORK( N+1 ),
 698      $                         A, LDA, WORK( 2*N+1 ), IERR )
 699                   IF( IERR.NE.0 )
 700      $               GO TO 90
 701                   CALL CUNM2R( 'L', 'N', N, N, N-1, Q, LDQ, WORK, B,
 702      $                         LDA, WORK( 2*N+1 ), IERR )
 703                   IF( IERR.NE.0 )
 704      $               GO TO 90
 705                   CALL CUNM2R( 'R', 'C', N, N, N-1, Z, LDQ, WORK( N+1 ),
 706      $                         B, LDA, WORK( 2*N+1 ), IERR )
 707                   IF( IERR.NE.0 )
 708      $               GO TO 90
 709                END IF
 710             ELSE
 711 *
 712 *              Random matrices
 713 *
 714                DO 80 JC = 1, N
 715                   DO 70 JR = 1, N
 716                      A( JR, JC ) = RMAGN( KAMAGN( JTYPE ) )*
 717      $                             CLARND( 4, ISEED )
 718                      B( JR, JC ) = RMAGN( KBMAGN( JTYPE ) )*
 719      $                             CLARND( 4, ISEED )
 720    70             CONTINUE
 721    80          CONTINUE
 722             END IF
 723 *
 724    90       CONTINUE
 725 *
 726             IF( IERR.NE.0 ) THEN
 727                WRITE( NOUNIT, FMT = 9999 )'Generator', IERR, N, JTYPE,
 728      $            IOLDSD
 729                INFO = ABS( IERR )
 730                RETURN
 731             END IF
 732 *
 733   100       CONTINUE
 734 *
 735             DO 110 I = 1, 7
 736                RESULT( I ) = -ONE
 737   110       CONTINUE
 738 *
 739 *           Call CGGEV to compute eigenvalues and eigenvectors.
 740 *
 741             CALL CLACPY( ' ', N, N, A, LDA, S, LDA )
 742             CALL CLACPY( ' ', N, N, B, LDA, T, LDA )
 743             CALL CGGEV( 'V', 'V', N, S, LDA, T, LDA, ALPHA, BETA, Q,
 744      $                  LDQ, Z, LDQ, WORK, LWORK, RWORK, IERR )
 745             IF( IERR.NE.0 .AND. IERR.NE.N+1 ) THEN
 746                RESULT( 1 ) = ULPINV
 747                WRITE( NOUNIT, FMT = 9999 )'CGGEV1', IERR, N, JTYPE,
 748      $            IOLDSD
 749                INFO = ABS( IERR )
 750                GO TO 190
 751             END IF
 752 *
 753 *           Do the tests (1) and (2)
 754 *
 755             CALL CGET52( .TRUE., N, A, LDA, B, LDA, Q, LDQ, ALPHA, BETA,
 756      $                   WORK, RWORK, RESULT( 1 ) )
 757             IF( RESULT( 2 ).GT.THRESH ) THEN
 758                WRITE( NOUNIT, FMT = 9998 )'Left', 'CGGEV1',
 759      $            RESULT( 2 ), N, JTYPE, IOLDSD
 760             END IF
 761 *
 762 *           Do the tests (3) and (4)
 763 *
 764             CALL CGET52( .FALSE., N, A, LDA, B, LDA, Z, LDQ, ALPHA,
 765      $                   BETA, WORK, RWORK, RESULT( 3 ) )
 766             IF( RESULT( 4 ).GT.THRESH ) THEN
 767                WRITE( NOUNIT, FMT = 9998 )'Right', 'CGGEV1',
 768      $            RESULT( 4 ), N, JTYPE, IOLDSD
 769             END IF
 770 *
 771 *           Do test (5)
 772 *
 773             CALL CLACPY( ' ', N, N, A, LDA, S, LDA )
 774             CALL CLACPY( ' ', N, N, B, LDA, T, LDA )
 775             CALL CGGEV( 'N', 'N', N, S, LDA, T, LDA, ALPHA1, BETA1, Q,
 776      $                  LDQ, Z, LDQ, WORK, LWORK, RWORK, IERR )
 777             IF( IERR.NE.0 .AND. IERR.NE.N+1 ) THEN
 778                RESULT( 1 ) = ULPINV
 779                WRITE( NOUNIT, FMT = 9999 )'CGGEV2', IERR, N, JTYPE,
 780      $            IOLDSD
 781                INFO = ABS( IERR )
 782                GO TO 190
 783             END IF
 784 *
 785             DO 120 J = 1, N
 786                IF( ALPHA( J ).NE.ALPHA1( J ) .OR. BETA( J ).NE.
 787      $             BETA1( J ) )RESULT( 5 ) = ULPINV
 788   120       CONTINUE
 789 *
 790 *           Do test (6): Compute eigenvalues and left eigenvectors,
 791 *           and test them
 792 *
 793             CALL CLACPY( ' ', N, N, A, LDA, S, LDA )
 794             CALL CLACPY( ' ', N, N, B, LDA, T, LDA )
 795             CALL CGGEV( 'V', 'N', N, S, LDA, T, LDA, ALPHA1, BETA1, QE,
 796      $                  LDQE, Z, LDQ, WORK, LWORK, RWORK, IERR )
 797             IF( IERR.NE.0 .AND. IERR.NE.N+1 ) THEN
 798                RESULT( 1 ) = ULPINV
 799                WRITE( NOUNIT, FMT = 9999 )'CGGEV3', IERR, N, JTYPE,
 800      $            IOLDSD
 801                INFO = ABS( IERR )
 802                GO TO 190
 803             END IF
 804 *
 805             DO 130 J = 1, N
 806                IF( ALPHA( J ).NE.ALPHA1( J ) .OR. BETA( J ).NE.
 807      $             BETA1( J ) )RESULT( 6 ) = ULPINV
 808   130       CONTINUE
 809 *
 810             DO 150 J = 1, N
 811                DO 140 JC = 1, N
 812                   IF( Q( J, JC ).NE.QE( J, JC ) )
 813      $               RESULT( 6 ) = ULPINV
 814   140          CONTINUE
 815   150       CONTINUE
 816 *
 817 *           Do test (7): Compute eigenvalues and right eigenvectors,
 818 *           and test them
 819 *
 820             CALL CLACPY( ' ', N, N, A, LDA, S, LDA )
 821             CALL CLACPY( ' ', N, N, B, LDA, T, LDA )
 822             CALL CGGEV( 'N', 'V', N, S, LDA, T, LDA, ALPHA1, BETA1, Q,
 823      $                  LDQ, QE, LDQE, WORK, LWORK, RWORK, IERR )
 824             IF( IERR.NE.0 .AND. IERR.NE.N+1 ) THEN
 825                RESULT( 1 ) = ULPINV
 826                WRITE( NOUNIT, FMT = 9999 )'CGGEV4', IERR, N, JTYPE,
 827      $            IOLDSD
 828                INFO = ABS( IERR )
 829                GO TO 190
 830             END IF
 831 *
 832             DO 160 J = 1, N
 833                IF( ALPHA( J ).NE.ALPHA1( J ) .OR. BETA( J ).NE.
 834      $             BETA1( J ) )RESULT( 7 ) = ULPINV
 835   160       CONTINUE
 836 *
 837             DO 180 J = 1, N
 838                DO 170 JC = 1, N
 839                   IF( Z( J, JC ).NE.QE( J, JC ) )
 840      $               RESULT( 7 ) = ULPINV
 841   170          CONTINUE
 842   180       CONTINUE
 843 *
 844 *           End of Loop -- Check for RESULT(j) > THRESH
 845 *
 846   190       CONTINUE
 847 *
 848             NTESTT = NTESTT + 7
 849 *
 850 *           Print out tests which fail.
 851 *
 852             DO 200 JR = 1, 7
 853                IF( RESULT( JR ).GE.THRESH ) THEN
 854 *
 855 *                 If this is the first test to fail,
 856 *                 print a header to the data file.
 857 *
 858                   IF( NERRS.EQ.0 ) THEN
 859                      WRITE( NOUNIT, FMT = 9997 )'CGV'
 860 *
 861 *                    Matrix types
 862 *
 863                      WRITE( NOUNIT, FMT = 9996 )
 864                      WRITE( NOUNIT, FMT = 9995 )
 865                      WRITE( NOUNIT, FMT = 9994 )'Orthogonal'
 866 *
 867 *                    Tests performed
 868 *
 869                      WRITE( NOUNIT, FMT = 9993 )
 870 *
 871                   END IF
 872                   NERRS = NERRS + 1
 873                   IF( RESULT( JR ).LT.10000.0 ) THEN
 874                      WRITE( NOUNIT, FMT = 9992 )N, JTYPE, IOLDSD, JR,
 875      $                  RESULT( JR )
 876                   ELSE
 877                      WRITE( NOUNIT, FMT = 9991 )N, JTYPE, IOLDSD, JR,
 878      $                  RESULT( JR )
 879                   END IF
 880                END IF
 881   200       CONTINUE
 882 *
 883   210    CONTINUE
 884   220 CONTINUE
 885 *
 886 *     Summary
 887 *
 888       CALL ALASVM( 'CGV', NOUNIT, NERRS, NTESTT, 0 )
 889 *
 890       WORK( 1 ) = MAXWRK
 891 *
 892       RETURN
 893 *
 894  9999 FORMAT( ' CDRGEV: ', A, ' returned INFO=', I6, '.', / 3X, 'N=',
 895      $      I6, ', JTYPE=', I6, ', ISEED=(', 3( I5, ',' ), I5, ')' )
 896 *
 897  9998 FORMAT( ' CDRGEV: ', A, ' Eigenvectors from ', A, ' incorrectly ',
 898      $      'normalized.', / ' Bits of error=', 0P, G10.3, ',', 3X,
 899      $      'N=', I4, ', JTYPE=', I3, ', ISEED=(', 3( I4, ',' ), I5,
 900      $      ')' )
 901 *
 902  9997 FORMAT( / 1X, A3, ' -- Complex Generalized eigenvalue problem ',
 903      $      'driver' )
 904 *
 905  9996 FORMAT( ' Matrix types (see CDRGEV for details): ' )
 906 *
 907  9995 FORMAT( ' Special Matrices:', 23X,
 908      $      '(J''=transposed Jordan block)',
 909      $      / '   1=(0,0)  2=(I,0)  3=(0,I)  4=(I,I)  5=(J'',J'')  ',
 910      $      '6=(diag(J'',I), diag(I,J''))', / ' Diagonal Matrices:  ( ',
 911      $      'D=diag(0,1,2,...) )', / '   7=(D,I)   9=(large*D, small*I',
 912      $      ')  11=(large*I, small*D)  13=(large*D, large*I)', /
 913      $      '   8=(I,D)  10=(small*D, large*I)  12=(small*I, large*D) ',
 914      $      ' 14=(small*D, small*I)', / '  15=(D, reversed D)' )
 915  9994 FORMAT( ' Matrices Rotated by Random ', A, ' Matrices U, V:',
 916      $      / '  16=Transposed Jordan Blocks             19=geometric ',
 917      $      'alpha, beta=0,1', / '  17=arithm. alpha&beta             ',
 918      $      '      20=arithmetic alpha, beta=0,1', / '  18=clustered ',
 919      $      'alpha, beta=0,1            21=random alpha, beta=0,1',
 920      $      / ' Large & Small Matrices:', / '  22=(large, small)   ',
 921      $      '23=(small,large)    24=(small,small)    25=(large,large)',
 922      $      / '  26=random O(1) matrices.' )
 923 *
 924  9993 FORMAT( / ' Tests performed:    ',
 925      $      / ' 1 = max | ( b A - a B )''*l | / const.,',
 926      $      / ' 2 = | |VR(i)| - 1 | / ulp,',
 927      $      / ' 3 = max | ( b A - a B )*r | / const.',
 928      $      / ' 4 = | |VL(i)| - 1 | / ulp,',
 929      $      / ' 5 = 0 if W same no matter if r or l computed,',
 930      $      / ' 6 = 0 if l same no matter if l computed,',
 931      $      / ' 7 = 0 if r same no matter if r computed,', / 1X )
 932  9992 FORMAT( ' Matrix order=', I5, ', type=', I2, ', seed=',
 933      $      4( I4, ',' ), ' result ', I2, ' is', 0P, F8.2 )
 934  9991 FORMAT( ' Matrix order=', I5, ', type=', I2, ', seed=',
 935      $      4( I4, ',' ), ' result ', I2, ' is', 1P, E10.3 )
 936 *
 937 *     End of CDRGEV
 938 *
 939       END