zdrgev.c
Go to the documentation of this file.
00001 /* zdrgev.f -- translated by f2c (version 20061008).
00002    You must link the resulting object file with libf2c:
00003         on Microsoft Windows system, link with libf2c.lib;
00004         on Linux or Unix systems, link with .../path/to/libf2c.a -lm
00005         or, if you install libf2c.a in a standard place, with -lf2c -lm
00006         -- in that order, at the end of the command line, as in
00007                 cc *.o -lf2c -lm
00008         Source for libf2c is in /netlib/f2c/libf2c.zip, e.g.,
00009 
00010                 http://www.netlib.org/f2c/libf2c.zip
00011 */
00012 
00013 #include "f2c.h"
00014 #include "blaswrap.h"
00015 
00016 /* Table of constant values */
00017 
00018 static doublecomplex c_b1 = {0.,0.};
00019 static integer c__1 = 1;
00020 static integer c_n1 = -1;
00021 static integer c__2 = 2;
00022 static doublereal c_b28 = 1.;
00023 static integer c__3 = 3;
00024 static integer c__4 = 4;
00025 static logical c_true = TRUE_;
00026 static logical c_false = FALSE_;
00027 static integer c__0 = 0;
00028 
00029 /* Subroutine */ int zdrgev_(integer *nsizes, integer *nn, integer *ntypes, 
00030         logical *dotype, integer *iseed, doublereal *thresh, integer *nounit, 
00031         doublecomplex *a, integer *lda, doublecomplex *b, doublecomplex *s, 
00032         doublecomplex *t, doublecomplex *q, integer *ldq, doublecomplex *z__, 
00033         doublecomplex *qe, integer *ldqe, doublecomplex *alpha, doublecomplex 
00034         *beta, doublecomplex *alpha1, doublecomplex *beta1, doublecomplex *
00035         work, integer *lwork, doublereal *rwork, doublereal *result, integer *
00036         info)
00037 {
00038     /* Initialized data */
00039 
00040     static integer kclass[26] = { 1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,2,2,2,2,2,2,2,
00041             2,2,2,3 };
00042     static integer kbmagn[26] = { 1,1,1,1,1,1,1,1,3,2,3,2,2,3,1,1,1,1,1,1,1,3,
00043             2,3,2,1 };
00044     static integer ktrian[26] = { 0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,1,1,1,1,1,
00045             1,1,1,1 };
00046     static logical lasign[26] = { FALSE_,FALSE_,FALSE_,FALSE_,FALSE_,FALSE_,
00047             TRUE_,FALSE_,TRUE_,TRUE_,FALSE_,FALSE_,TRUE_,TRUE_,TRUE_,FALSE_,
00048             TRUE_,FALSE_,FALSE_,FALSE_,TRUE_,TRUE_,TRUE_,TRUE_,TRUE_,FALSE_ };
00049     static logical lbsign[26] = { FALSE_,FALSE_,FALSE_,FALSE_,FALSE_,FALSE_,
00050             FALSE_,TRUE_,FALSE_,FALSE_,TRUE_,TRUE_,FALSE_,FALSE_,TRUE_,FALSE_,
00051             TRUE_,FALSE_,FALSE_,FALSE_,FALSE_,FALSE_,FALSE_,FALSE_,FALSE_,
00052             FALSE_ };
00053     static integer kz1[6] = { 0,1,2,1,3,3 };
00054     static integer kz2[6] = { 0,0,1,2,1,1 };
00055     static integer kadd[6] = { 0,0,0,0,3,2 };
00056     static integer katype[26] = { 0,1,0,1,2,3,4,1,4,4,1,1,4,4,4,2,4,5,8,7,9,4,
00057             4,4,4,0 };
00058     static integer kbtype[26] = { 0,0,1,1,2,-3,1,4,1,1,4,4,1,1,-4,2,-4,8,8,8,
00059             8,8,8,8,8,0 };
00060     static integer kazero[26] = { 1,1,1,1,1,1,2,1,2,2,1,1,2,2,3,1,3,5,5,5,5,3,
00061             3,3,3,1 };
00062     static integer kbzero[26] = { 1,1,1,1,1,1,1,2,1,1,2,2,1,1,4,1,4,6,6,6,6,4,
00063             4,4,4,1 };
00064     static integer kamagn[26] = { 1,1,1,1,1,1,1,1,2,3,2,3,2,3,1,1,1,1,1,1,1,2,
00065             3,3,2,1 };
00066 
00067     /* Format strings */
00068     static char fmt_9999[] = "(\002 ZDRGEV: \002,a,\002 returned INFO=\002,i"
00069             "6,\002.\002,/3x,\002N=\002,i6,\002, JTYPE=\002,i6,\002, ISEED="
00070             "(\002,3(i5,\002,\002),i5,\002)\002)";
00071     static char fmt_9998[] = "(\002 ZDRGEV: \002,a,\002 Eigenvectors from"
00072             " \002,a,\002 incorrectly \002,\002normalized.\002,/\002 Bits of "
00073             "error=\002,0p,g10.3,\002,\002,3x,\002N=\002,i4,\002, JTYPE=\002,"
00074             "i3,\002, ISEED=(\002,3(i4,\002,\002),i5,\002)\002)";
00075     static char fmt_9997[] = "(/1x,a3,\002 -- Complex Generalized eigenvalue"
00076             " problem \002,\002driver\002)";
00077     static char fmt_9996[] = "(\002 Matrix types (see ZDRGEV for details):"
00078             " \002)";
00079     static char fmt_9995[] = "(\002 Special Matrices:\002,23x,\002(J'=transp"
00080             "osed Jordan block)\002,/\002   1=(0,0)  2=(I,0)  3=(0,I)  4=(I,I"
00081             ")  5=(J',J')  \002,\0026=(diag(J',I), diag(I,J'))\002,/\002 Diag"
00082             "onal Matrices:  ( \002,\002D=diag(0,1,2,...) )\002,/\002   7=(D,"
00083             "I)   9=(large*D, small*I\002,\002)  11=(large*I, small*D)  13=(l"
00084             "arge*D, large*I)\002,/\002   8=(I,D)  10=(small*D, large*I)  12="
00085             "(small*I, large*D) \002,\002 14=(small*D, small*I)\002,/\002  15"
00086             "=(D, reversed D)\002)";
00087     static char fmt_9994[] = "(\002 Matrices Rotated by Random \002,a,\002 M"
00088             "atrices U, V:\002,/\002  16=Transposed Jordan Blocks            "
00089             " 19=geometric \002,\002alpha, beta=0,1\002,/\002  17=arithm. alp"
00090             "ha&beta             \002,\002      20=arithmetic alpha, beta=0,"
00091             "1\002,/\002  18=clustered \002,\002alpha, beta=0,1            21"
00092             "=random alpha, beta=0,1\002,/\002 Large & Small Matrices:\002,"
00093             "/\002  22=(large, small)   \002,\00223=(small,large)    24=(smal"
00094             "l,small)    25=(large,large)\002,/\002  26=random O(1) matrices"
00095             ".\002)";
00096     static char fmt_9993[] = "(/\002 Tests performed:    \002,/\002 1 = max "
00097             "| ( b A - a B )'*l | / const.,\002,/\002 2 = | |VR(i)| - 1 | / u"
00098             "lp,\002,/\002 3 = max | ( b A - a B )*r | / const.\002,/\002 4 ="
00099             " | |VL(i)| - 1 | / ulp,\002,/\002 5 = 0 if W same no matter if r"
00100             " or l computed,\002,/\002 6 = 0 if l same no matter if l compute"
00101             "d,\002,/\002 7 = 0 if r same no matter if r computed,\002,/1x)";
00102     static char fmt_9992[] = "(\002 Matrix order=\002,i5,\002, type=\002,i2"
00103             ",\002, seed=\002,4(i4,\002,\002),\002 result \002,i2,\002 is\002"
00104             ",0p,f8.2)";
00105     static char fmt_9991[] = "(\002 Matrix order=\002,i5,\002, type=\002,i2"
00106             ",\002, seed=\002,4(i4,\002,\002),\002 result \002,i2,\002 is\002"
00107             ",1p,d10.3)";
00108 
00109     /* System generated locals */
00110     integer a_dim1, a_offset, b_dim1, b_offset, q_dim1, q_offset, qe_dim1, 
00111             qe_offset, s_dim1, s_offset, t_dim1, t_offset, z_dim1, z_offset, 
00112             i__1, i__2, i__3, i__4, i__5, i__6, i__7;
00113     doublereal d__1, d__2;
00114     doublecomplex z__1, z__2, z__3;
00115 
00116     /* Builtin functions */
00117     double d_sign(doublereal *, doublereal *), z_abs(doublecomplex *);
00118     void d_cnjg(doublecomplex *, doublecomplex *);
00119     integer s_wsfe(cilist *), do_fio(integer *, char *, ftnlen), e_wsfe(void);
00120 
00121     /* Local variables */
00122     integer i__, j, n, n1, jc, nb, in, jr;
00123     doublereal ulp;
00124     integer iadd, ierr, nmax;
00125     logical badnn;
00126     doublereal rmagn[4];
00127     doublecomplex ctemp;
00128     extern /* Subroutine */ int zget52_(logical *, integer *, doublecomplex *, 
00129              integer *, doublecomplex *, integer *, doublecomplex *, integer *
00130 , doublecomplex *, doublecomplex *, doublecomplex *, doublereal *, 
00131              doublereal *);
00132     integer nmats, jsize;
00133     extern /* Subroutine */ int zggev_(char *, char *, integer *, 
00134             doublecomplex *, integer *, doublecomplex *, integer *, 
00135             doublecomplex *, doublecomplex *, doublecomplex *, integer *, 
00136             doublecomplex *, integer *, doublecomplex *, integer *, 
00137             doublereal *, integer *);
00138     integer nerrs, jtype;
00139     extern /* Subroutine */ int dlabad_(doublereal *, doublereal *), zlatm4_(
00140             integer *, integer *, integer *, integer *, logical *, doublereal 
00141             *, doublereal *, doublereal *, integer *, integer *, 
00142             doublecomplex *, integer *);
00143     extern doublereal dlamch_(char *);
00144     extern /* Subroutine */ int zunm2r_(char *, char *, integer *, integer *, 
00145             integer *, doublecomplex *, integer *, doublecomplex *, 
00146             doublecomplex *, integer *, doublecomplex *, integer *);
00147     doublereal safmin, safmax;
00148     integer ioldsd[4];
00149     extern integer ilaenv_(integer *, char *, char *, integer *, integer *, 
00150             integer *, integer *);
00151     extern /* Subroutine */ int alasvm_(char *, integer *, integer *, integer 
00152             *, integer *), xerbla_(char *, integer *), 
00153             zlarfg_(integer *, doublecomplex *, doublecomplex *, integer *, 
00154             doublecomplex *);
00155     extern /* Double Complex */ VOID zlarnd_(doublecomplex *, integer *, 
00156             integer *);
00157     extern /* Subroutine */ int zlacpy_(char *, integer *, integer *, 
00158             doublecomplex *, integer *, doublecomplex *, integer *), 
00159             zlaset_(char *, integer *, integer *, doublecomplex *, 
00160             doublecomplex *, doublecomplex *, integer *);
00161     integer minwrk, maxwrk;
00162     doublereal ulpinv;
00163     integer mtypes, ntestt;
00164 
00165     /* Fortran I/O blocks */
00166     static cilist io___40 = { 0, 0, 0, fmt_9999, 0 };
00167     static cilist io___42 = { 0, 0, 0, fmt_9999, 0 };
00168     static cilist io___43 = { 0, 0, 0, fmt_9998, 0 };
00169     static cilist io___44 = { 0, 0, 0, fmt_9998, 0 };
00170     static cilist io___45 = { 0, 0, 0, fmt_9999, 0 };
00171     static cilist io___46 = { 0, 0, 0, fmt_9999, 0 };
00172     static cilist io___47 = { 0, 0, 0, fmt_9999, 0 };
00173     static cilist io___48 = { 0, 0, 0, fmt_9997, 0 };
00174     static cilist io___49 = { 0, 0, 0, fmt_9996, 0 };
00175     static cilist io___50 = { 0, 0, 0, fmt_9995, 0 };
00176     static cilist io___51 = { 0, 0, 0, fmt_9994, 0 };
00177     static cilist io___52 = { 0, 0, 0, fmt_9993, 0 };
00178     static cilist io___53 = { 0, 0, 0, fmt_9992, 0 };
00179     static cilist io___54 = { 0, 0, 0, fmt_9991, 0 };
00180 
00181 
00182 
00183 /*  -- LAPACK test routine (version 3.1.1) -- */
00184 /*     Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd.. */
00185 /*     February 2007 */
00186 
00187 /*     .. Scalar Arguments .. */
00188 /*     .. */
00189 /*     .. Array Arguments .. */
00190 /*     .. */
00191 
00192 /*  Purpose */
00193 /*  ======= */
00194 
00195 /*  ZDRGEV checks the nonsymmetric generalized eigenvalue problem driver */
00196 /*  routine ZGGEV. */
00197 
00198 /*  ZGGEV computes for a pair of n-by-n nonsymmetric matrices (A,B) the */
00199 /*  generalized eigenvalues and, optionally, the left and right */
00200 /*  eigenvectors. */
00201 
00202 /*  A generalized eigenvalue for a pair of matrices (A,B) is a scalar w */
00203 /*  or a ratio  alpha/beta = w, such that A - w*B is singular.  It is */
00204 /*  usually represented as the pair (alpha,beta), as there is reasonalbe */
00205 /*  interpretation for beta=0, and even for both being zero. */
00206 
00207 /*  A right generalized eigenvector corresponding to a generalized */
00208 /*  eigenvalue  w  for a pair of matrices (A,B) is a vector r  such that */
00209 /*  (A - wB) * r = 0.  A left generalized eigenvector is a vector l such */
00210 /*  that l**H * (A - wB) = 0, where l**H is the conjugate-transpose of l. */
00211 
00212 /*  When ZDRGEV is called, a number of matrix "sizes" ("n's") and a */
00213 /*  number of matrix "types" are specified.  For each size ("n") */
00214 /*  and each type of matrix, a pair of matrices (A, B) will be generated */
00215 /*  and used for testing.  For each matrix pair, the following tests */
00216 /*  will be performed and compared with the threshhold THRESH. */
00217 
00218 /*  Results from ZGGEV: */
00219 
00220 /*  (1)  max over all left eigenvalue/-vector pairs (alpha/beta,l) of */
00221 
00222 /*       | VL**H * (beta A - alpha B) |/( ulp max(|beta A|, |alpha B|) ) */
00223 
00224 /*       where VL**H is the conjugate-transpose of VL. */
00225 
00226 /*  (2)  | |VL(i)| - 1 | / ulp and whether largest component real */
00227 
00228 /*       VL(i) denotes the i-th column of VL. */
00229 
00230 /*  (3)  max over all left eigenvalue/-vector pairs (alpha/beta,r) of */
00231 
00232 /*       | (beta A - alpha B) * VR | / ( ulp max(|beta A|, |alpha B|) ) */
00233 
00234 /*  (4)  | |VR(i)| - 1 | / ulp and whether largest component real */
00235 
00236 /*       VR(i) denotes the i-th column of VR. */
00237 
00238 /*  (5)  W(full) = W(partial) */
00239 /*       W(full) denotes the eigenvalues computed when both l and r */
00240 /*       are also computed, and W(partial) denotes the eigenvalues */
00241 /*       computed when only W, only W and r, or only W and l are */
00242 /*       computed. */
00243 
00244 /*  (6)  VL(full) = VL(partial) */
00245 /*       VL(full) denotes the left eigenvectors computed when both l */
00246 /*       and r are computed, and VL(partial) denotes the result */
00247 /*       when only l is computed. */
00248 
00249 /*  (7)  VR(full) = VR(partial) */
00250 /*       VR(full) denotes the right eigenvectors computed when both l */
00251 /*       and r are also computed, and VR(partial) denotes the result */
00252 /*       when only l is computed. */
00253 
00254 
00255 /*  Test Matrices */
00256 /*  ---- -------- */
00257 
00258 /*  The sizes of the test matrices are specified by an array */
00259 /*  NN(1:NSIZES); the value of each element NN(j) specifies one size. */
00260 /*  The "types" are specified by a logical array DOTYPE( 1:NTYPES ); if */
00261 /*  DOTYPE(j) is .TRUE., then matrix type "j" will be generated. */
00262 /*  Currently, the list of possible types is: */
00263 
00264 /*  (1)  ( 0, 0 )         (a pair of zero matrices) */
00265 
00266 /*  (2)  ( I, 0 )         (an identity and a zero matrix) */
00267 
00268 /*  (3)  ( 0, I )         (an identity and a zero matrix) */
00269 
00270 /*  (4)  ( I, I )         (a pair of identity matrices) */
00271 
00272 /*          t   t */
00273 /*  (5)  ( J , J  )       (a pair of transposed Jordan blocks) */
00274 
00275 /*                                      t                ( I   0  ) */
00276 /*  (6)  ( X, Y )         where  X = ( J   0  )  and Y = (      t ) */
00277 /*                                   ( 0   I  )          ( 0   J  ) */
00278 /*                        and I is a k x k identity and J a (k+1)x(k+1) */
00279 /*                        Jordan block; k=(N-1)/2 */
00280 
00281 /*  (7)  ( D, I )         where D is diag( 0, 1,..., N-1 ) (a diagonal */
00282 /*                        matrix with those diagonal entries.) */
00283 /*  (8)  ( I, D ) */
00284 
00285 /*  (9)  ( big*D, small*I ) where "big" is near overflow and small=1/big */
00286 
00287 /*  (10) ( small*D, big*I ) */
00288 
00289 /*  (11) ( big*I, small*D ) */
00290 
00291 /*  (12) ( small*I, big*D ) */
00292 
00293 /*  (13) ( big*D, big*I ) */
00294 
00295 /*  (14) ( small*D, small*I ) */
00296 
00297 /*  (15) ( D1, D2 )        where D1 is diag( 0, 0, 1, ..., N-3, 0 ) and */
00298 /*                         D2 is diag( 0, N-3, N-4,..., 1, 0, 0 ) */
00299 /*            t   t */
00300 /*  (16) Q ( J , J ) Z     where Q and Z are random orthogonal matrices. */
00301 
00302 /*  (17) Q ( T1, T2 ) Z    where T1 and T2 are upper triangular matrices */
00303 /*                         with random O(1) entries above the diagonal */
00304 /*                         and diagonal entries diag(T1) = */
00305 /*                         ( 0, 0, 1, ..., N-3, 0 ) and diag(T2) = */
00306 /*                         ( 0, N-3, N-4,..., 1, 0, 0 ) */
00307 
00308 /*  (18) Q ( T1, T2 ) Z    diag(T1) = ( 0, 0, 1, 1, s, ..., s, 0 ) */
00309 /*                         diag(T2) = ( 0, 1, 0, 1,..., 1, 0 ) */
00310 /*                         s = machine precision. */
00311 
00312 /*  (19) Q ( T1, T2 ) Z    diag(T1)=( 0,0,1,1, 1-d, ..., 1-(N-5)*d=s, 0 ) */
00313 /*                         diag(T2) = ( 0, 1, 0, 1, ..., 1, 0 ) */
00314 
00315 /*                                                         N-5 */
00316 /*  (20) Q ( T1, T2 ) Z    diag(T1)=( 0, 0, 1, 1, a, ..., a   =s, 0 ) */
00317 /*                         diag(T2) = ( 0, 1, 0, 1, ..., 1, 0, 0 ) */
00318 
00319 /*  (21) Q ( T1, T2 ) Z    diag(T1)=( 0, 0, 1, r1, r2, ..., r(N-4), 0 ) */
00320 /*                         diag(T2) = ( 0, 1, 0, 1, ..., 1, 0, 0 ) */
00321 /*                         where r1,..., r(N-4) are random. */
00322 
00323 /*  (22) Q ( big*T1, small*T2 ) Z    diag(T1) = ( 0, 0, 1, ..., N-3, 0 ) */
00324 /*                                   diag(T2) = ( 0, 1, ..., 1, 0, 0 ) */
00325 
00326 /*  (23) Q ( small*T1, big*T2 ) Z    diag(T1) = ( 0, 0, 1, ..., N-3, 0 ) */
00327 /*                                   diag(T2) = ( 0, 1, ..., 1, 0, 0 ) */
00328 
00329 /*  (24) Q ( small*T1, small*T2 ) Z  diag(T1) = ( 0, 0, 1, ..., N-3, 0 ) */
00330 /*                                   diag(T2) = ( 0, 1, ..., 1, 0, 0 ) */
00331 
00332 /*  (25) Q ( big*T1, big*T2 ) Z      diag(T1) = ( 0, 0, 1, ..., N-3, 0 ) */
00333 /*                                   diag(T2) = ( 0, 1, ..., 1, 0, 0 ) */
00334 
00335 /*  (26) Q ( T1, T2 ) Z     where T1 and T2 are random upper-triangular */
00336 /*                          matrices. */
00337 
00338 
00339 /*  Arguments */
00340 /*  ========= */
00341 
00342 /*  NSIZES  (input) INTEGER */
00343 /*          The number of sizes of matrices to use.  If it is zero, */
00344 /*          ZDRGES does nothing.  NSIZES >= 0. */
00345 
00346 /*  NN      (input) INTEGER array, dimension (NSIZES) */
00347 /*          An array containing the sizes to be used for the matrices. */
00348 /*          Zero values will be skipped.  NN >= 0. */
00349 
00350 /*  NTYPES  (input) INTEGER */
00351 /*          The number of elements in DOTYPE.   If it is zero, ZDRGEV */
00352 /*          does nothing.  It must be at least zero.  If it is MAXTYP+1 */
00353 /*          and NSIZES is 1, then an additional type, MAXTYP+1 is */
00354 /*          defined, which is to use whatever matrix is in A.  This */
00355 /*          is only useful if DOTYPE(1:MAXTYP) is .FALSE. and */
00356 /*          DOTYPE(MAXTYP+1) is .TRUE. . */
00357 
00358 /*  DOTYPE  (input) LOGICAL array, dimension (NTYPES) */
00359 /*          If DOTYPE(j) is .TRUE., then for each size in NN a */
00360 /*          matrix of that size and of type j will be generated. */
00361 /*          If NTYPES is smaller than the maximum number of types */
00362 /*          defined (PARAMETER MAXTYP), then types NTYPES+1 through */
00363 /*          MAXTYP will not be generated. If NTYPES is larger */
00364 /*          than MAXTYP, DOTYPE(MAXTYP+1) through DOTYPE(NTYPES) */
00365 /*          will be ignored. */
00366 
00367 /*  ISEED   (input/output) INTEGER array, dimension (4) */
00368 /*          On entry ISEED specifies the seed of the random number */
00369 /*          generator. The array elements should be between 0 and 4095; */
00370 /*          if not they will be reduced mod 4096. Also, ISEED(4) must */
00371 /*          be odd.  The random number generator uses a linear */
00372 /*          congruential sequence limited to small integers, and so */
00373 /*          should produce machine independent random numbers. The */
00374 /*          values of ISEED are changed on exit, and can be used in the */
00375 /*          next call to ZDRGES to continue the same random number */
00376 /*          sequence. */
00377 
00378 /*  THRESH  (input) DOUBLE PRECISION */
00379 /*          A test will count as "failed" if the "error", computed as */
00380 /*          described above, exceeds THRESH.  Note that the error is */
00381 /*          scaled to be O(1), so THRESH should be a reasonably small */
00382 /*          multiple of 1, e.g., 10 or 100.  In particular, it should */
00383 /*          not depend on the precision (single vs. double) or the size */
00384 /*          of the matrix.  It must be at least zero. */
00385 
00386 /*  NOUNIT  (input) INTEGER */
00387 /*          The FORTRAN unit number for printing out error messages */
00388 /*          (e.g., if a routine returns IERR not equal to 0.) */
00389 
00390 /*  A       (input/workspace) COMPLEX*16 array, dimension(LDA, max(NN)) */
00391 /*          Used to hold the original A matrix.  Used as input only */
00392 /*          if NTYPES=MAXTYP+1, DOTYPE(1:MAXTYP)=.FALSE., and */
00393 /*          DOTYPE(MAXTYP+1)=.TRUE. */
00394 
00395 /*  LDA     (input) INTEGER */
00396 /*          The leading dimension of A, B, S, and T. */
00397 /*          It must be at least 1 and at least max( NN ). */
00398 
00399 /*  B       (input/workspace) COMPLEX*16 array, dimension(LDA, max(NN)) */
00400 /*          Used to hold the original B matrix.  Used as input only */
00401 /*          if NTYPES=MAXTYP+1, DOTYPE(1:MAXTYP)=.FALSE., and */
00402 /*          DOTYPE(MAXTYP+1)=.TRUE. */
00403 
00404 /*  S       (workspace) COMPLEX*16 array, dimension (LDA, max(NN)) */
00405 /*          The Schur form matrix computed from A by ZGGEV.  On exit, S */
00406 /*          contains the Schur form matrix corresponding to the matrix */
00407 /*          in A. */
00408 
00409 /*  T       (workspace) COMPLEX*16 array, dimension (LDA, max(NN)) */
00410 /*          The upper triangular matrix computed from B by ZGGEV. */
00411 
00412 /*  Q      (workspace) COMPLEX*16 array, dimension (LDQ, max(NN)) */
00413 /*          The (left) eigenvectors matrix computed by ZGGEV. */
00414 
00415 /*  LDQ     (input) INTEGER */
00416 /*          The leading dimension of Q and Z. It must */
00417 /*          be at least 1 and at least max( NN ). */
00418 
00419 /*  Z       (workspace) COMPLEX*16 array, dimension( LDQ, max(NN) ) */
00420 /*          The (right) orthogonal matrix computed by ZGGEV. */
00421 
00422 /*  QE      (workspace) COMPLEX*16 array, dimension( LDQ, max(NN) ) */
00423 /*          QE holds the computed right or left eigenvectors. */
00424 
00425 /*  LDQE    (input) INTEGER */
00426 /*          The leading dimension of QE. LDQE >= max(1,max(NN)). */
00427 
00428 /*  ALPHA   (workspace) COMPLEX*16 array, dimension (max(NN)) */
00429 /*  BETA    (workspace) COMPLEX*16 array, dimension (max(NN)) */
00430 /*          The generalized eigenvalues of (A,B) computed by ZGGEV. */
00431 /*          ( ALPHAR(k)+ALPHAI(k)*i ) / BETA(k) is the k-th */
00432 /*          generalized eigenvalue of A and B. */
00433 
00434 /*  ALPHA1  (workspace) COMPLEX*16 array, dimension (max(NN)) */
00435 /*  BETA1   (workspace) COMPLEX*16 array, dimension (max(NN)) */
00436 /*          Like ALPHAR, ALPHAI, BETA, these arrays contain the */
00437 /*          eigenvalues of A and B, but those computed when ZGGEV only */
00438 /*          computes a partial eigendecomposition, i.e. not the */
00439 /*          eigenvalues and left and right eigenvectors. */
00440 
00441 /*  WORK    (workspace) COMPLEX*16 array, dimension (LWORK) */
00442 
00443 /*  LWORK   (input) INTEGER */
00444 /*          The number of entries in WORK.  LWORK >= N*(N+1) */
00445 
00446 /*  RWORK   (workspace) DOUBLE PRECISION array, dimension (8*N) */
00447 /*          Real workspace. */
00448 
00449 /*  RESULT  (output) DOUBLE PRECISION array, dimension (2) */
00450 /*          The values computed by the tests described above. */
00451 /*          The values are currently limited to 1/ulp, to avoid overflow. */
00452 
00453 /*  INFO    (output) INTEGER */
00454 /*          = 0:  successful exit */
00455 /*          < 0:  if INFO = -i, the i-th argument had an illegal value. */
00456 /*          > 0:  A routine returned an error code.  INFO is the */
00457 /*                absolute value of the INFO value returned. */
00458 
00459 /*  ===================================================================== */
00460 
00461 /*     .. Parameters .. */
00462 /*     .. */
00463 /*     .. Local Scalars .. */
00464 /*     .. */
00465 /*     .. Local Arrays .. */
00466 /*     .. */
00467 /*     .. External Functions .. */
00468 /*     .. */
00469 /*     .. External Subroutines .. */
00470 /*     .. */
00471 /*     .. Intrinsic Functions .. */
00472 /*     .. */
00473 /*     .. Data statements .. */
00474     /* Parameter adjustments */
00475     --nn;
00476     --dotype;
00477     --iseed;
00478     t_dim1 = *lda;
00479     t_offset = 1 + t_dim1;
00480     t -= t_offset;
00481     s_dim1 = *lda;
00482     s_offset = 1 + s_dim1;
00483     s -= s_offset;
00484     b_dim1 = *lda;
00485     b_offset = 1 + b_dim1;
00486     b -= b_offset;
00487     a_dim1 = *lda;
00488     a_offset = 1 + a_dim1;
00489     a -= a_offset;
00490     z_dim1 = *ldq;
00491     z_offset = 1 + z_dim1;
00492     z__ -= z_offset;
00493     q_dim1 = *ldq;
00494     q_offset = 1 + q_dim1;
00495     q -= q_offset;
00496     qe_dim1 = *ldqe;
00497     qe_offset = 1 + qe_dim1;
00498     qe -= qe_offset;
00499     --alpha;
00500     --beta;
00501     --alpha1;
00502     --beta1;
00503     --work;
00504     --rwork;
00505     --result;
00506 
00507     /* Function Body */
00508 /*     .. */
00509 /*     .. Executable Statements .. */
00510 
00511 /*     Check for errors */
00512 
00513     *info = 0;
00514 
00515     badnn = FALSE_;
00516     nmax = 1;
00517     i__1 = *nsizes;
00518     for (j = 1; j <= i__1; ++j) {
00519 /* Computing MAX */
00520         i__2 = nmax, i__3 = nn[j];
00521         nmax = max(i__2,i__3);
00522         if (nn[j] < 0) {
00523             badnn = TRUE_;
00524         }
00525 /* L10: */
00526     }
00527 
00528     if (*nsizes < 0) {
00529         *info = -1;
00530     } else if (badnn) {
00531         *info = -2;
00532     } else if (*ntypes < 0) {
00533         *info = -3;
00534     } else if (*thresh < 0.) {
00535         *info = -6;
00536     } else if (*lda <= 1 || *lda < nmax) {
00537         *info = -9;
00538     } else if (*ldq <= 1 || *ldq < nmax) {
00539         *info = -14;
00540     } else if (*ldqe <= 1 || *ldqe < nmax) {
00541         *info = -17;
00542     }
00543 
00544 /*     Compute workspace */
00545 /*      (Note: Comments in the code beginning "Workspace:" describe the */
00546 /*       minimal amount of workspace needed at that point in the code, */
00547 /*       as well as the preferred amount for good performance. */
00548 /*       NB refers to the optimal block size for the immediately */
00549 /*       following subroutine, as returned by ILAENV. */
00550 
00551     minwrk = 1;
00552     if (*info == 0 && *lwork >= 1) {
00553         minwrk = nmax * (nmax + 1);
00554 /* Computing MAX */
00555         i__1 = 1, i__2 = ilaenv_(&c__1, "ZGEQRF", " ", &nmax, &nmax, &c_n1, &
00556                 c_n1), i__1 = max(i__1,i__2), i__2 = 
00557                 ilaenv_(&c__1, "ZUNMQR", "LC", &nmax, &nmax, &nmax, &c_n1), i__1 = max(i__1,i__2), i__2 = ilaenv_(&
00558                 c__1, "ZUNGQR", " ", &nmax, &nmax, &nmax, &c_n1);
00559         nb = max(i__1,i__2);
00560 /* Computing MAX */
00561         i__1 = nmax << 1, i__2 = nmax * (nb + 1), i__1 = max(i__1,i__2), i__2 
00562                 = nmax * (nmax + 1);
00563         maxwrk = max(i__1,i__2);
00564         work[1].r = (doublereal) maxwrk, work[1].i = 0.;
00565     }
00566 
00567     if (*lwork < minwrk) {
00568         *info = -23;
00569     }
00570 
00571     if (*info != 0) {
00572         i__1 = -(*info);
00573         xerbla_("ZDRGEV", &i__1);
00574         return 0;
00575     }
00576 
00577 /*     Quick return if possible */
00578 
00579     if (*nsizes == 0 || *ntypes == 0) {
00580         return 0;
00581     }
00582 
00583     ulp = dlamch_("Precision");
00584     safmin = dlamch_("Safe minimum");
00585     safmin /= ulp;
00586     safmax = 1. / safmin;
00587     dlabad_(&safmin, &safmax);
00588     ulpinv = 1. / ulp;
00589 
00590 /*     The values RMAGN(2:3) depend on N, see below. */
00591 
00592     rmagn[0] = 0.;
00593     rmagn[1] = 1.;
00594 
00595 /*     Loop over sizes, types */
00596 
00597     ntestt = 0;
00598     nerrs = 0;
00599     nmats = 0;
00600 
00601     i__1 = *nsizes;
00602     for (jsize = 1; jsize <= i__1; ++jsize) {
00603         n = nn[jsize];
00604         n1 = max(1,n);
00605         rmagn[2] = safmax * ulp / (doublereal) n1;
00606         rmagn[3] = safmin * ulpinv * n1;
00607 
00608         if (*nsizes != 1) {
00609             mtypes = min(26,*ntypes);
00610         } else {
00611             mtypes = min(27,*ntypes);
00612         }
00613 
00614         i__2 = mtypes;
00615         for (jtype = 1; jtype <= i__2; ++jtype) {
00616             if (! dotype[jtype]) {
00617                 goto L210;
00618             }
00619             ++nmats;
00620 
00621 /*           Save ISEED in case of an error. */
00622 
00623             for (j = 1; j <= 4; ++j) {
00624                 ioldsd[j - 1] = iseed[j];
00625 /* L20: */
00626             }
00627 
00628 /*           Generate test matrices A and B */
00629 
00630 /*           Description of control parameters: */
00631 
00632 /*           KZLASS: =1 means w/o rotation, =2 means w/ rotation, */
00633 /*                   =3 means random. */
00634 /*           KATYPE: the "type" to be passed to ZLATM4 for computing A. */
00635 /*           KAZERO: the pattern of zeros on the diagonal for A: */
00636 /*                   =1: ( xxx ), =2: (0, xxx ) =3: ( 0, 0, xxx, 0 ), */
00637 /*                   =4: ( 0, xxx, 0, 0 ), =5: ( 0, 0, 1, xxx, 0 ), */
00638 /*                   =6: ( 0, 1, 0, xxx, 0 ).  (xxx means a string of */
00639 /*                   non-zero entries.) */
00640 /*           KAMAGN: the magnitude of the matrix: =0: zero, =1: O(1), */
00641 /*                   =2: large, =3: small. */
00642 /*           LASIGN: .TRUE. if the diagonal elements of A are to be */
00643 /*                   multiplied by a random magnitude 1 number. */
00644 /*           KBTYPE, KBZERO, KBMAGN, LBSIGN: the same, but for B. */
00645 /*           KTRIAN: =0: don't fill in the upper triangle, =1: do. */
00646 /*           KZ1, KZ2, KADD: used to implement KAZERO and KBZERO. */
00647 /*           RMAGN: used to implement KAMAGN and KBMAGN. */
00648 
00649             if (mtypes > 26) {
00650                 goto L100;
00651             }
00652             ierr = 0;
00653             if (kclass[jtype - 1] < 3) {
00654 
00655 /*              Generate A (w/o rotation) */
00656 
00657                 if ((i__3 = katype[jtype - 1], abs(i__3)) == 3) {
00658                     in = ((n - 1) / 2 << 1) + 1;
00659                     if (in != n) {
00660                         zlaset_("Full", &n, &n, &c_b1, &c_b1, &a[a_offset], 
00661                                 lda);
00662                     }
00663                 } else {
00664                     in = n;
00665                 }
00666                 zlatm4_(&katype[jtype - 1], &in, &kz1[kazero[jtype - 1] - 1], 
00667                         &kz2[kazero[jtype - 1] - 1], &lasign[jtype - 1], &
00668                         rmagn[kamagn[jtype - 1]], &ulp, &rmagn[ktrian[jtype - 
00669                         1] * kamagn[jtype - 1]], &c__2, &iseed[1], &a[
00670                         a_offset], lda);
00671                 iadd = kadd[kazero[jtype - 1] - 1];
00672                 if (iadd > 0 && iadd <= n) {
00673                     i__3 = iadd + iadd * a_dim1;
00674                     i__4 = kamagn[jtype - 1];
00675                     a[i__3].r = rmagn[i__4], a[i__3].i = 0.;
00676                 }
00677 
00678 /*              Generate B (w/o rotation) */
00679 
00680                 if ((i__3 = kbtype[jtype - 1], abs(i__3)) == 3) {
00681                     in = ((n - 1) / 2 << 1) + 1;
00682                     if (in != n) {
00683                         zlaset_("Full", &n, &n, &c_b1, &c_b1, &b[b_offset], 
00684                                 lda);
00685                     }
00686                 } else {
00687                     in = n;
00688                 }
00689                 zlatm4_(&kbtype[jtype - 1], &in, &kz1[kbzero[jtype - 1] - 1], 
00690                         &kz2[kbzero[jtype - 1] - 1], &lbsign[jtype - 1], &
00691                         rmagn[kbmagn[jtype - 1]], &c_b28, &rmagn[ktrian[jtype 
00692                         - 1] * kbmagn[jtype - 1]], &c__2, &iseed[1], &b[
00693                         b_offset], lda);
00694                 iadd = kadd[kbzero[jtype - 1] - 1];
00695                 if (iadd != 0 && iadd <= n) {
00696                     i__3 = iadd + iadd * b_dim1;
00697                     i__4 = kbmagn[jtype - 1];
00698                     b[i__3].r = rmagn[i__4], b[i__3].i = 0.;
00699                 }
00700 
00701                 if (kclass[jtype - 1] == 2 && n > 0) {
00702 
00703 /*                 Include rotations */
00704 
00705 /*                 Generate Q, Z as Householder transformations times */
00706 /*                 a diagonal matrix. */
00707 
00708                     i__3 = n - 1;
00709                     for (jc = 1; jc <= i__3; ++jc) {
00710                         i__4 = n;
00711                         for (jr = jc; jr <= i__4; ++jr) {
00712                             i__5 = jr + jc * q_dim1;
00713                             zlarnd_(&z__1, &c__3, &iseed[1]);
00714                             q[i__5].r = z__1.r, q[i__5].i = z__1.i;
00715                             i__5 = jr + jc * z_dim1;
00716                             zlarnd_(&z__1, &c__3, &iseed[1]);
00717                             z__[i__5].r = z__1.r, z__[i__5].i = z__1.i;
00718 /* L30: */
00719                         }
00720                         i__4 = n + 1 - jc;
00721                         zlarfg_(&i__4, &q[jc + jc * q_dim1], &q[jc + 1 + jc * 
00722                                 q_dim1], &c__1, &work[jc]);
00723                         i__4 = (n << 1) + jc;
00724                         i__5 = jc + jc * q_dim1;
00725                         d__2 = q[i__5].r;
00726                         d__1 = d_sign(&c_b28, &d__2);
00727                         work[i__4].r = d__1, work[i__4].i = 0.;
00728                         i__4 = jc + jc * q_dim1;
00729                         q[i__4].r = 1., q[i__4].i = 0.;
00730                         i__4 = n + 1 - jc;
00731                         zlarfg_(&i__4, &z__[jc + jc * z_dim1], &z__[jc + 1 + 
00732                                 jc * z_dim1], &c__1, &work[n + jc]);
00733                         i__4 = n * 3 + jc;
00734                         i__5 = jc + jc * z_dim1;
00735                         d__2 = z__[i__5].r;
00736                         d__1 = d_sign(&c_b28, &d__2);
00737                         work[i__4].r = d__1, work[i__4].i = 0.;
00738                         i__4 = jc + jc * z_dim1;
00739                         z__[i__4].r = 1., z__[i__4].i = 0.;
00740 /* L40: */
00741                     }
00742                     zlarnd_(&z__1, &c__3, &iseed[1]);
00743                     ctemp.r = z__1.r, ctemp.i = z__1.i;
00744                     i__3 = n + n * q_dim1;
00745                     q[i__3].r = 1., q[i__3].i = 0.;
00746                     i__3 = n;
00747                     work[i__3].r = 0., work[i__3].i = 0.;
00748                     i__3 = n * 3;
00749                     d__1 = z_abs(&ctemp);
00750                     z__1.r = ctemp.r / d__1, z__1.i = ctemp.i / d__1;
00751                     work[i__3].r = z__1.r, work[i__3].i = z__1.i;
00752                     zlarnd_(&z__1, &c__3, &iseed[1]);
00753                     ctemp.r = z__1.r, ctemp.i = z__1.i;
00754                     i__3 = n + n * z_dim1;
00755                     z__[i__3].r = 1., z__[i__3].i = 0.;
00756                     i__3 = n << 1;
00757                     work[i__3].r = 0., work[i__3].i = 0.;
00758                     i__3 = n << 2;
00759                     d__1 = z_abs(&ctemp);
00760                     z__1.r = ctemp.r / d__1, z__1.i = ctemp.i / d__1;
00761                     work[i__3].r = z__1.r, work[i__3].i = z__1.i;
00762 
00763 /*                 Apply the diagonal matrices */
00764 
00765                     i__3 = n;
00766                     for (jc = 1; jc <= i__3; ++jc) {
00767                         i__4 = n;
00768                         for (jr = 1; jr <= i__4; ++jr) {
00769                             i__5 = jr + jc * a_dim1;
00770                             i__6 = (n << 1) + jr;
00771                             d_cnjg(&z__3, &work[n * 3 + jc]);
00772                             z__2.r = work[i__6].r * z__3.r - work[i__6].i * 
00773                                     z__3.i, z__2.i = work[i__6].r * z__3.i + 
00774                                     work[i__6].i * z__3.r;
00775                             i__7 = jr + jc * a_dim1;
00776                             z__1.r = z__2.r * a[i__7].r - z__2.i * a[i__7].i, 
00777                                     z__1.i = z__2.r * a[i__7].i + z__2.i * a[
00778                                     i__7].r;
00779                             a[i__5].r = z__1.r, a[i__5].i = z__1.i;
00780                             i__5 = jr + jc * b_dim1;
00781                             i__6 = (n << 1) + jr;
00782                             d_cnjg(&z__3, &work[n * 3 + jc]);
00783                             z__2.r = work[i__6].r * z__3.r - work[i__6].i * 
00784                                     z__3.i, z__2.i = work[i__6].r * z__3.i + 
00785                                     work[i__6].i * z__3.r;
00786                             i__7 = jr + jc * b_dim1;
00787                             z__1.r = z__2.r * b[i__7].r - z__2.i * b[i__7].i, 
00788                                     z__1.i = z__2.r * b[i__7].i + z__2.i * b[
00789                                     i__7].r;
00790                             b[i__5].r = z__1.r, b[i__5].i = z__1.i;
00791 /* L50: */
00792                         }
00793 /* L60: */
00794                     }
00795                     i__3 = n - 1;
00796                     zunm2r_("L", "N", &n, &n, &i__3, &q[q_offset], ldq, &work[
00797                             1], &a[a_offset], lda, &work[(n << 1) + 1], &ierr);
00798                     if (ierr != 0) {
00799                         goto L90;
00800                     }
00801                     i__3 = n - 1;
00802                     zunm2r_("R", "C", &n, &n, &i__3, &z__[z_offset], ldq, &
00803                             work[n + 1], &a[a_offset], lda, &work[(n << 1) + 
00804                             1], &ierr);
00805                     if (ierr != 0) {
00806                         goto L90;
00807                     }
00808                     i__3 = n - 1;
00809                     zunm2r_("L", "N", &n, &n, &i__3, &q[q_offset], ldq, &work[
00810                             1], &b[b_offset], lda, &work[(n << 1) + 1], &ierr);
00811                     if (ierr != 0) {
00812                         goto L90;
00813                     }
00814                     i__3 = n - 1;
00815                     zunm2r_("R", "C", &n, &n, &i__3, &z__[z_offset], ldq, &
00816                             work[n + 1], &b[b_offset], lda, &work[(n << 1) + 
00817                             1], &ierr);
00818                     if (ierr != 0) {
00819                         goto L90;
00820                     }
00821                 }
00822             } else {
00823 
00824 /*              Random matrices */
00825 
00826                 i__3 = n;
00827                 for (jc = 1; jc <= i__3; ++jc) {
00828                     i__4 = n;
00829                     for (jr = 1; jr <= i__4; ++jr) {
00830                         i__5 = jr + jc * a_dim1;
00831                         i__6 = kamagn[jtype - 1];
00832                         zlarnd_(&z__2, &c__4, &iseed[1]);
00833                         z__1.r = rmagn[i__6] * z__2.r, z__1.i = rmagn[i__6] * 
00834                                 z__2.i;
00835                         a[i__5].r = z__1.r, a[i__5].i = z__1.i;
00836                         i__5 = jr + jc * b_dim1;
00837                         i__6 = kbmagn[jtype - 1];
00838                         zlarnd_(&z__2, &c__4, &iseed[1]);
00839                         z__1.r = rmagn[i__6] * z__2.r, z__1.i = rmagn[i__6] * 
00840                                 z__2.i;
00841                         b[i__5].r = z__1.r, b[i__5].i = z__1.i;
00842 /* L70: */
00843                     }
00844 /* L80: */
00845                 }
00846             }
00847 
00848 L90:
00849 
00850             if (ierr != 0) {
00851                 io___40.ciunit = *nounit;
00852                 s_wsfe(&io___40);
00853                 do_fio(&c__1, "Generator", (ftnlen)9);
00854                 do_fio(&c__1, (char *)&ierr, (ftnlen)sizeof(integer));
00855                 do_fio(&c__1, (char *)&n, (ftnlen)sizeof(integer));
00856                 do_fio(&c__1, (char *)&jtype, (ftnlen)sizeof(integer));
00857                 do_fio(&c__4, (char *)&ioldsd[0], (ftnlen)sizeof(integer));
00858                 e_wsfe();
00859                 *info = abs(ierr);
00860                 return 0;
00861             }
00862 
00863 L100:
00864 
00865             for (i__ = 1; i__ <= 7; ++i__) {
00866                 result[i__] = -1.;
00867 /* L110: */
00868             }
00869 
00870 /*           Call ZGGEV to compute eigenvalues and eigenvectors. */
00871 
00872             zlacpy_(" ", &n, &n, &a[a_offset], lda, &s[s_offset], lda);
00873             zlacpy_(" ", &n, &n, &b[b_offset], lda, &t[t_offset], lda);
00874             zggev_("V", "V", &n, &s[s_offset], lda, &t[t_offset], lda, &alpha[
00875                     1], &beta[1], &q[q_offset], ldq, &z__[z_offset], ldq, &
00876                     work[1], lwork, &rwork[1], &ierr);
00877             if (ierr != 0 && ierr != n + 1) {
00878                 result[1] = ulpinv;
00879                 io___42.ciunit = *nounit;
00880                 s_wsfe(&io___42);
00881                 do_fio(&c__1, "ZGGEV1", (ftnlen)6);
00882                 do_fio(&c__1, (char *)&ierr, (ftnlen)sizeof(integer));
00883                 do_fio(&c__1, (char *)&n, (ftnlen)sizeof(integer));
00884                 do_fio(&c__1, (char *)&jtype, (ftnlen)sizeof(integer));
00885                 do_fio(&c__4, (char *)&ioldsd[0], (ftnlen)sizeof(integer));
00886                 e_wsfe();
00887                 *info = abs(ierr);
00888                 goto L190;
00889             }
00890 
00891 /*           Do the tests (1) and (2) */
00892 
00893             zget52_(&c_true, &n, &a[a_offset], lda, &b[b_offset], lda, &q[
00894                     q_offset], ldq, &alpha[1], &beta[1], &work[1], &rwork[1], 
00895                     &result[1]);
00896             if (result[2] > *thresh) {
00897                 io___43.ciunit = *nounit;
00898                 s_wsfe(&io___43);
00899                 do_fio(&c__1, "Left", (ftnlen)4);
00900                 do_fio(&c__1, "ZGGEV1", (ftnlen)6);
00901                 do_fio(&c__1, (char *)&result[2], (ftnlen)sizeof(doublereal));
00902                 do_fio(&c__1, (char *)&n, (ftnlen)sizeof(integer));
00903                 do_fio(&c__1, (char *)&jtype, (ftnlen)sizeof(integer));
00904                 do_fio(&c__4, (char *)&ioldsd[0], (ftnlen)sizeof(integer));
00905                 e_wsfe();
00906             }
00907 
00908 /*           Do the tests (3) and (4) */
00909 
00910             zget52_(&c_false, &n, &a[a_offset], lda, &b[b_offset], lda, &z__[
00911                     z_offset], ldq, &alpha[1], &beta[1], &work[1], &rwork[1], 
00912                     &result[3]);
00913             if (result[4] > *thresh) {
00914                 io___44.ciunit = *nounit;
00915                 s_wsfe(&io___44);
00916                 do_fio(&c__1, "Right", (ftnlen)5);
00917                 do_fio(&c__1, "ZGGEV1", (ftnlen)6);
00918                 do_fio(&c__1, (char *)&result[4], (ftnlen)sizeof(doublereal));
00919                 do_fio(&c__1, (char *)&n, (ftnlen)sizeof(integer));
00920                 do_fio(&c__1, (char *)&jtype, (ftnlen)sizeof(integer));
00921                 do_fio(&c__4, (char *)&ioldsd[0], (ftnlen)sizeof(integer));
00922                 e_wsfe();
00923             }
00924 
00925 /*           Do test (5) */
00926 
00927             zlacpy_(" ", &n, &n, &a[a_offset], lda, &s[s_offset], lda);
00928             zlacpy_(" ", &n, &n, &b[b_offset], lda, &t[t_offset], lda);
00929             zggev_("N", "N", &n, &s[s_offset], lda, &t[t_offset], lda, &
00930                     alpha1[1], &beta1[1], &q[q_offset], ldq, &z__[z_offset], 
00931                     ldq, &work[1], lwork, &rwork[1], &ierr);
00932             if (ierr != 0 && ierr != n + 1) {
00933                 result[1] = ulpinv;
00934                 io___45.ciunit = *nounit;
00935                 s_wsfe(&io___45);
00936                 do_fio(&c__1, "ZGGEV2", (ftnlen)6);
00937                 do_fio(&c__1, (char *)&ierr, (ftnlen)sizeof(integer));
00938                 do_fio(&c__1, (char *)&n, (ftnlen)sizeof(integer));
00939                 do_fio(&c__1, (char *)&jtype, (ftnlen)sizeof(integer));
00940                 do_fio(&c__4, (char *)&ioldsd[0], (ftnlen)sizeof(integer));
00941                 e_wsfe();
00942                 *info = abs(ierr);
00943                 goto L190;
00944             }
00945 
00946             i__3 = n;
00947             for (j = 1; j <= i__3; ++j) {
00948                 i__4 = j;
00949                 i__5 = j;
00950                 i__6 = j;
00951                 i__7 = j;
00952                 if (alpha[i__4].r != alpha1[i__5].r || alpha[i__4].i != 
00953                         alpha1[i__5].i || (beta[i__6].r != beta1[i__7].r || 
00954                         beta[i__6].i != beta1[i__7].i)) {
00955                     result[5] = ulpinv;
00956                 }
00957 /* L120: */
00958             }
00959 
00960 /*           Do test (6): Compute eigenvalues and left eigenvectors, */
00961 /*           and test them */
00962 
00963             zlacpy_(" ", &n, &n, &a[a_offset], lda, &s[s_offset], lda);
00964             zlacpy_(" ", &n, &n, &b[b_offset], lda, &t[t_offset], lda);
00965             zggev_("V", "N", &n, &s[s_offset], lda, &t[t_offset], lda, &
00966                     alpha1[1], &beta1[1], &qe[qe_offset], ldqe, &z__[z_offset]
00967 , ldq, &work[1], lwork, &rwork[1], &ierr);
00968             if (ierr != 0 && ierr != n + 1) {
00969                 result[1] = ulpinv;
00970                 io___46.ciunit = *nounit;
00971                 s_wsfe(&io___46);
00972                 do_fio(&c__1, "ZGGEV3", (ftnlen)6);
00973                 do_fio(&c__1, (char *)&ierr, (ftnlen)sizeof(integer));
00974                 do_fio(&c__1, (char *)&n, (ftnlen)sizeof(integer));
00975                 do_fio(&c__1, (char *)&jtype, (ftnlen)sizeof(integer));
00976                 do_fio(&c__4, (char *)&ioldsd[0], (ftnlen)sizeof(integer));
00977                 e_wsfe();
00978                 *info = abs(ierr);
00979                 goto L190;
00980             }
00981 
00982             i__3 = n;
00983             for (j = 1; j <= i__3; ++j) {
00984                 i__4 = j;
00985                 i__5 = j;
00986                 i__6 = j;
00987                 i__7 = j;
00988                 if (alpha[i__4].r != alpha1[i__5].r || alpha[i__4].i != 
00989                         alpha1[i__5].i || (beta[i__6].r != beta1[i__7].r || 
00990                         beta[i__6].i != beta1[i__7].i)) {
00991                     result[6] = ulpinv;
00992                 }
00993 /* L130: */
00994             }
00995 
00996             i__3 = n;
00997             for (j = 1; j <= i__3; ++j) {
00998                 i__4 = n;
00999                 for (jc = 1; jc <= i__4; ++jc) {
01000                     i__5 = j + jc * q_dim1;
01001                     i__6 = j + jc * qe_dim1;
01002                     if (q[i__5].r != qe[i__6].r || q[i__5].i != qe[i__6].i) {
01003                         result[6] = ulpinv;
01004                     }
01005 /* L140: */
01006                 }
01007 /* L150: */
01008             }
01009 
01010 /*           Do test (7): Compute eigenvalues and right eigenvectors, */
01011 /*           and test them */
01012 
01013             zlacpy_(" ", &n, &n, &a[a_offset], lda, &s[s_offset], lda);
01014             zlacpy_(" ", &n, &n, &b[b_offset], lda, &t[t_offset], lda);
01015             zggev_("N", "V", &n, &s[s_offset], lda, &t[t_offset], lda, &
01016                     alpha1[1], &beta1[1], &q[q_offset], ldq, &qe[qe_offset], 
01017                     ldqe, &work[1], lwork, &rwork[1], &ierr);
01018             if (ierr != 0 && ierr != n + 1) {
01019                 result[1] = ulpinv;
01020                 io___47.ciunit = *nounit;
01021                 s_wsfe(&io___47);
01022                 do_fio(&c__1, "ZGGEV4", (ftnlen)6);
01023                 do_fio(&c__1, (char *)&ierr, (ftnlen)sizeof(integer));
01024                 do_fio(&c__1, (char *)&n, (ftnlen)sizeof(integer));
01025                 do_fio(&c__1, (char *)&jtype, (ftnlen)sizeof(integer));
01026                 do_fio(&c__4, (char *)&ioldsd[0], (ftnlen)sizeof(integer));
01027                 e_wsfe();
01028                 *info = abs(ierr);
01029                 goto L190;
01030             }
01031 
01032             i__3 = n;
01033             for (j = 1; j <= i__3; ++j) {
01034                 i__4 = j;
01035                 i__5 = j;
01036                 i__6 = j;
01037                 i__7 = j;
01038                 if (alpha[i__4].r != alpha1[i__5].r || alpha[i__4].i != 
01039                         alpha1[i__5].i || (beta[i__6].r != beta1[i__7].r || 
01040                         beta[i__6].i != beta1[i__7].i)) {
01041                     result[7] = ulpinv;
01042                 }
01043 /* L160: */
01044             }
01045 
01046             i__3 = n;
01047             for (j = 1; j <= i__3; ++j) {
01048                 i__4 = n;
01049                 for (jc = 1; jc <= i__4; ++jc) {
01050                     i__5 = j + jc * z_dim1;
01051                     i__6 = j + jc * qe_dim1;
01052                     if (z__[i__5].r != qe[i__6].r || z__[i__5].i != qe[i__6]
01053                             .i) {
01054                         result[7] = ulpinv;
01055                     }
01056 /* L170: */
01057                 }
01058 /* L180: */
01059             }
01060 
01061 /*           End of Loop -- Check for RESULT(j) > THRESH */
01062 
01063 L190:
01064 
01065             ntestt += 7;
01066 
01067 /*           Print out tests which fail. */
01068 
01069             for (jr = 1; jr <= 7; ++jr) {
01070                 if (result[jr] >= *thresh) {
01071 
01072 /*                 If this is the first test to fail, */
01073 /*                 print a header to the data file. */
01074 
01075                     if (nerrs == 0) {
01076                         io___48.ciunit = *nounit;
01077                         s_wsfe(&io___48);
01078                         do_fio(&c__1, "ZGV", (ftnlen)3);
01079                         e_wsfe();
01080 
01081 /*                    Matrix types */
01082 
01083                         io___49.ciunit = *nounit;
01084                         s_wsfe(&io___49);
01085                         e_wsfe();
01086                         io___50.ciunit = *nounit;
01087                         s_wsfe(&io___50);
01088                         e_wsfe();
01089                         io___51.ciunit = *nounit;
01090                         s_wsfe(&io___51);
01091                         do_fio(&c__1, "Orthogonal", (ftnlen)10);
01092                         e_wsfe();
01093 
01094 /*                    Tests performed */
01095 
01096                         io___52.ciunit = *nounit;
01097                         s_wsfe(&io___52);
01098                         e_wsfe();
01099 
01100                     }
01101                     ++nerrs;
01102                     if (result[jr] < 1e4) {
01103                         io___53.ciunit = *nounit;
01104                         s_wsfe(&io___53);
01105                         do_fio(&c__1, (char *)&n, (ftnlen)sizeof(integer));
01106                         do_fio(&c__1, (char *)&jtype, (ftnlen)sizeof(integer))
01107                                 ;
01108                         do_fio(&c__4, (char *)&ioldsd[0], (ftnlen)sizeof(
01109                                 integer));
01110                         do_fio(&c__1, (char *)&jr, (ftnlen)sizeof(integer));
01111                         do_fio(&c__1, (char *)&result[jr], (ftnlen)sizeof(
01112                                 doublereal));
01113                         e_wsfe();
01114                     } else {
01115                         io___54.ciunit = *nounit;
01116                         s_wsfe(&io___54);
01117                         do_fio(&c__1, (char *)&n, (ftnlen)sizeof(integer));
01118                         do_fio(&c__1, (char *)&jtype, (ftnlen)sizeof(integer))
01119                                 ;
01120                         do_fio(&c__4, (char *)&ioldsd[0], (ftnlen)sizeof(
01121                                 integer));
01122                         do_fio(&c__1, (char *)&jr, (ftnlen)sizeof(integer));
01123                         do_fio(&c__1, (char *)&result[jr], (ftnlen)sizeof(
01124                                 doublereal));
01125                         e_wsfe();
01126                     }
01127                 }
01128 /* L200: */
01129             }
01130 
01131 L210:
01132             ;
01133         }
01134 /* L220: */
01135     }
01136 
01137 /*     Summary */
01138 
01139     alasvm_("ZGV", nounit, &nerrs, &ntestt, &c__0);
01140 
01141     work[1].r = (doublereal) maxwrk, work[1].i = 0.;
01142 
01143     return 0;
01144 
01145 
01146 
01147 
01148 
01149 
01150 
01151 /*     End of ZDRGEV */
01152 
01153 } /* zdrgev_ */


swiftnav
Author(s):
autogenerated on Sat Jun 8 2019 18:56:21