ScaLAPACK 2.1  2.1
ScaLAPACK: Scalable Linear Algebra PACKage
pdgsepsubtst.f
Go to the documentation of this file.
1 *
2 *
3  SUBROUTINE pdgsepsubtst( WKNOWN, IBTYPE, JOBZ, RANGE, UPLO, N, VL,
4  $ VU, IL, IU, THRESH, ABSTOL, A, COPYA, B,
5  $ COPYB, Z, IA, JA, DESCA, WIN, WNEW,
6  $ IFAIL, ICLUSTR, GAP, IPREPAD, IPOSTPAD,
7  $ WORK, LWORK, LWORK1, IWORK, LIWORK,
8  $ RESULT, TSTNRM, QTQNRM, NOUT )
9 *
10 * -- ScaLAPACK routine (version 1.7) --
11 * University of Tennessee, Knoxville, Oak Ridge National Laboratory,
12 * and University of California, Berkeley.
13 * May 1, 1997
14 *
15 * .. Scalar Arguments ..
16  LOGICAL WKNOWN
17  CHARACTER JOBZ, RANGE, UPLO
18  INTEGER IA, IBTYPE, IL, IPOSTPAD, IPREPAD, IU, JA,
19  $ LIWORK, LWORK, LWORK1, N, NOUT, RESULT
20  DOUBLE PRECISION ABSTOL, QTQNRM, THRESH, TSTNRM, VL, VU
21 * ..
22 * .. Array Arguments ..
23  INTEGER DESCA( * ), ICLUSTR( * ), IFAIL( * ),
24  $ IWORK( * )
25  DOUBLE PRECISION A( * ), B( * ), COPYA( * ), COPYB( * ),
26  $ GAP( * ), WIN( * ), WNEW( * ), WORK( * ),
27  $ Z( * )
28 * ..
29 *
30 * Purpose
31 * =======
32 *
33 * PDGSEPSUBTST calls PDSYGVX and then tests the output of
34 * PDSYGVX
35 * If JOBZ = 'V' then the following two tests are performed:
36 * |AQ -QL| / (abstol + eps * norm(A) ) < THRESH
37 * |QT * Q - I| / eps * norm(A) < THRESH
38 * If WKNOWN then
39 * we check to make sure that the eigenvalues match expectations
40 * i.e. |WIN - WNEW(1+IPREPAD)| / (eps * |WIN|) < THRESH
41 * where WIN is the array of eigenvalues as computed by
42 * PDSYGVX when eigenvectors are requested
43 *
44 * Arguments
45 * =========
46 *
47 * NP = the number of rows local to a given process.
48 * NQ = the number of columns local to a given process.
49 *
50 * WKNOWN (global input) INTEGER
51 * .FALSE.: WIN does not contain the eigenvalues
52 * .TRUE.: WIN does contain the eigenvalues
53 *
54 * IBTYPE (global input) INTEGER
55 * Specifies the problem type to be solved:
56 * = 1: sub( A )*x = (lambda)*sub( B )*x
57 * = 2: sub( A )*sub( B )*x = (lambda)*x
58 * = 3: sub( B )*sub( A )*x = (lambda)*x
59 *
60 *
61 * JOBZ (global input) CHARACTER*1
62 * Specifies whether or not to compute the eigenvectors:
63 * = 'N': Compute eigenvalues only.
64 * = 'V': Compute eigenvalues and eigenvectors.
65 * Must be 'V' on first call to PDGSEPSUBTST
66 *
67 * RANGE (global input) CHARACTER*1
68 * = 'A': all eigenvalues will be found.
69 * = 'V': all eigenvalues in the interval [VL,VU]
70 * will be found.
71 * = 'I': the IL-th through IU-th eigenvalues will be found.
72 * Must be 'A' on first call to PDGSEPSUBTST
73 *
74 * UPLO (global input) CHARACTER*1
75 * Specifies whether the upper or lower triangular part of the
76 * symmetric matrix A is stored:
77 * = 'U': Upper triangular
78 * = 'L': Lower triangular
79 *
80 * N (global input) INTEGER
81 * Size of the matrix to be tested. (global size)
82 *
83 * VL (global input) DOUBLE PRECISION
84 * If RANGE='V', the lower bound of the interval to be searched
85 * for eigenvalues. Not referenced if RANGE = 'A' or 'I'.
86 *
87 * VU (global input) DOUBLE PRECISION
88 * If RANGE='V', the upper bound of the interval to be searched
89 * for eigenvalues. Not referenced if RANGE = 'A' or 'I'.
90 *
91 * IL (global input) INTEGER
92 * If RANGE='I', the index (from smallest to largest) of the
93 * smallest eigenvalue to be returned. IL >= 1.
94 * Not referenced if RANGE = 'A' or 'V'.
95 *
96 * IU (global input) INTEGER
97 * If RANGE='I', the index (from smallest to largest) of the
98 * largest eigenvalue to be returned. min(IL,N) <= IU <= N.
99 * Not referenced if RANGE = 'A' or 'V'.
100 *
101 * THRESH (global input) DOUBLE PRECISION
102 * A test will count as "failed" if the "error", computed as
103 * described below, exceeds THRESH. Note that the error
104 * is scaled to be O(1), so THRESH should be a reasonably
105 * small multiple of 1, e.g., 10 or 100. In particular,
106 * it should not depend on the precision (single vs. double)
107 * or the size of the matrix. It must be at least zero.
108 *
109 * ABSTOL (global input) DOUBLE PRECISION
110 * The absolute tolerance for the eigenvalues. An
111 * eigenvalue is considered to be located if it has
112 * been determined to lie in an interval whose width
113 * is "abstol" or less. If "abstol" is less than or equal
114 * to zero, then ulp*|T| will be used, where |T| is
115 * the 1-norm of the matrix. If eigenvectors are
116 * desired later by inverse iteration ("PDSTEIN"),
117 * "abstol" MUST NOT be bigger than ulp*|T|.
118 *
119 * A (local workspace) DOUBLE PRECISION array
120 * global dimension (N, N), local dimension (DESCA(DLEN_), NQ)
121 * A is distributed in a block cyclic manner over both rows
122 * and columns.
123 * See PDSYGVX for a description of block cyclic layout.
124 * The test matrix, which is then modified by PDSYGVX
125 * A has already been padded front and back, use A(1+IPREPAD)
126 *
127 * COPYA (local input) DOUBLE PRECISION array, dimension(N*N)
128 * COPYA holds a copy of the original matrix A
129 * identical in both form and content to A
130 *
131 * B (local workspace) DOUBLE PRECISION array, dim (N*N)
132 * global dimension (N, N), local dimension (LDA, NQ)
133 * A is distributed in a block cyclic manner over both rows
134 * and columns.
135 * The B test matrix, which is then modified by PDSYGVX
136 *
137 * COPYB (local input) DOUBLE PRECISION array, dim (N, N)
138 * COPYB is used to hold an identical copy of the array B
139 * identical in both form and content to B
140 *
141 * Z (local workspace) DOUBLE PRECISION array, dim (N*N)
142 * Z is distributed in the same manner as A
143 * Z contains the eigenvector matrix
144 * Z is used as workspace by the test routines
145 * PDGSEPCHK and PDSEPQTQ.
146 * Z has already been padded front and back, use Z(1+IPREPAD)
147 *
148 * IA (global input) INTEGER
149 * On entry, IA specifies the global row index of the submatrix
150 * of the global matrix A, COPYA and Z to operate on.
151 *
152 * JA (global input) INTEGER
153 * On entry, IA specifies the global column index of the submat
154 * of the global matrix A, COPYA and Z to operate on.
155 *
156 * DESCA (global/local input) INTEGER array of dimension 8
157 * The array descriptor for the matrix A, COPYA and Z.
158 *
159 * WIN (global input) DOUBLE PRECISION array, dimension (N)
160 * If .not. WKNOWN, WIN is ignored on input
161 * Otherwise, WIN() is taken as the standard by which the
162 * eigenvalues are to be compared against.
163 *
164 * WNEW (global workspace) DOUBLE PRECISION array, dimension (N)
165 * The eigenvalues as copmuted by this call to PDSYGVX
166 * If JOBZ <> 'V' or RANGE <> 'A' these eigenvalues are
167 * compared against those in WIN().
168 * WNEW has already been padded front and back,
169 * use WNEW(1+IPREPAD)
170 *
171 * IFAIL (global output) INTEGER array, dimension (N)
172 * If JOBZ = 'V', then on normal exit, the first M elements of
173 * IFAIL are zero. If INFO > 0 on exit, then IFAIL contains the
174 * indices of the eigenvectors that failed to converge.
175 * If JOBZ = 'N', then IFAIL is not referenced.
176 * IFAIL has already been padded front and back,
177 * use IFAIL(1+IPREPAD)
178 *
179 * ICLUSTR (global workspace) integer array, dimension (2*NPROW*NPCOL)
180 *
181 * GAP (global workspace) DOUBLE PRECISION array,
182 * dimension (NPROW*NPCOL)
183 *
184 * WORK (local workspace) DOUBLE PRECISION array, dimension (LWORK)
185 * WORK has already been padded front and back,
186 * use WORK(1+IPREPAD)
187 *
188 * LWORK (local input) INTEGER
189 * The actual length of the array WORK after padding.
190 *
191 *
192 * LWORK1 (local input) INTEGER
193 * The amount of real workspace to pass to PDSYGVX
194 *
195 * IWORK (local workspace) INTEGER array, dimension (LIWORK)
196 * IWORK has already been padded front and back,
197 * use IWORK(1+IPREPAD)
198 *
199 * LIWORK (local input) INTEGER
200 * The length of the array IWORK after padding.
201 *
202 * RESULT (global output) INTEGER
203 * The result of this call to PDSYGVX
204 * RESULT = -3 => This process did not participate
205 * RESULT = 0 => All tests passed
206 * RESULT = 1 => ONe or more tests failed
207 *
208 * TSTNRM (global output) DOUBLE PRECISION
209 * |AQ- QL| / |A|*N*EPS
210 *
211 * QTQNRM (global output) DOUBLE PRECISION
212 * |QTQ -I| / N*EPS
213 *
214 * .. Parameters ..
215 *
216  INTEGER BLOCK_CYCLIC_2D, DLEN_, DTYPE_, CTXT_, M_, N_,
217  $ MB_, NB_, RSRC_, CSRC_, LLD_
218  PARAMETER ( BLOCK_CYCLIC_2D = 1, dlen_ = 9, dtype_ = 1,
219  $ ctxt_ = 2, m_ = 3, n_ = 4, mb_ = 5, nb_ = 6,
220  $ rsrc_ = 7, csrc_ = 8, lld_ = 9 )
221  DOUBLE PRECISION PADVAL, FIVE, NEGONE
222  PARAMETER ( PADVAL = 13.5285d+0, five = 5.0d+0,
223  $ negone = -1.0d+0 )
224  INTEGER IPADVAL
225  PARAMETER ( IPADVAL = 927 )
226 * ..
227 * .. Local Scalars ..
228  LOGICAL MISSLARGEST, MISSSMALLEST
229  INTEGER I, IAM, INDIWRK, INFO, ISIZESUBTST, ISIZESYEVX,
230  $ isizetst, j, m, maxeigs, maxil, maxiu, maxsize,
231  $ minil, mq, mycol, myil, myrow, nclusters, np,
232  $ npcol, nprow, nq, nz, oldil, oldiu, oldnz, res,
233  $ sizechk, sizemqrleft, sizemqrright, sizeqrf,
234  $ sizeqtq, sizesubtst, sizesyevx, sizetms,
235  $ sizetst, valsize, vecsize
236  DOUBLE PRECISION EPS, ERROR, MAXERROR, MAXVU, MINERROR, MINVL,
237  $ NORMWIN, OLDVL, OLDVU, ORFAC, SAFMIN
238 * ..
239 * .. Local Arrays ..
240  INTEGER DESCZ( DLEN_ ), DSEED( 4 ), ITMP( 2 )
241 * ..
242 * .. External Functions ..
243 *
244  LOGICAL LSAME
245  INTEGER NUMROC
246  DOUBLE PRECISION PDLAMCH
247  EXTERNAL LSAME, NUMROC, PDLAMCH
248 * ..
249 * .. External Subroutines ..
250  EXTERNAL blacs_gridinfo, descinit, dgamn2d, dgamx2d,
251  $ dlacpy, igamn2d, igamx2d, pdchekpad, pdelset,
254  $ slboot, sltimer
255 * ..
256 * .. Intrinsic Functions ..
257  INTRINSIC abs, max, min, mod
258 * ..
259 * .. Executable Statements ..
260 * This is just to keep ftnchek happy
261  IF( block_cyclic_2d*csrc_*ctxt_*dlen_*dtype_*lld_*mb_*m_*nb_*n_*
262  $ rsrc_.LT.0 )RETURN
263  CALL pdlasizegsep( desca, iprepad, ipostpad, sizemqrleft,
264  $ sizemqrright, sizeqrf, sizetms, sizeqtq,
265  $ sizechk, sizesyevx, isizesyevx, sizesubtst,
266  $ isizesubtst, sizetst, isizetst )
267 *
268  tstnrm = negone
269  qtqnrm = negone
270  eps = pdlamch( desca( ctxt_ ), 'Eps' )
271  safmin = pdlamch( desca( ctxt_ ), 'Safe min' )
272 *
273  normwin = safmin / eps
274  IF( n.GE.1 )
275  $ normwin = max( abs( win( 1 ) ), abs( win( n ) ), normwin )
276 *
277 * Make sure that we aren't using information from previous calls
278 *
279  nz = -13
280  oldnz = nz
281  oldil = il
282  oldiu = iu
283  oldvl = vl
284  oldvu = vu
285 *
286  DO 10 i = 1, lwork1, 1
287  work( i+iprepad ) = 14.3d+0
288  10 CONTINUE
289  DO 20 i = 1, liwork, 1
290  iwork( i+iprepad ) = 14
291  20 CONTINUE
292 *
293  DO 30 i = 1, n
294  wnew( i+iprepad ) = 3.14159d+0
295  30 CONTINUE
296 *
297  iclustr( 1+iprepad ) = 139
298 *
299  IF( lsame( jobz, 'N' ) ) THEN
300  maxeigs = 0
301  ELSE
302  IF( lsame( range, 'A' ) ) THEN
303  maxeigs = n
304  ELSE IF( lsame( range, 'I' ) ) THEN
305  maxeigs = iu - il + 1
306  ELSE
307  minvl = vl - normwin*five*eps - abstol
308  maxvu = vu + normwin*five*eps + abstol
309  minil = 1
310  maxiu = 0
311  DO 40 i = 1, n
312  IF( win( i ).LT.minvl )
313  $ minil = minil + 1
314  IF( win( i ).LE.maxvu )
315  $ maxiu = maxiu + 1
316  40 CONTINUE
317 *
318  maxeigs = maxiu - minil + 1
319  END IF
320  END IF
321 *
322 *
323  CALL descinit( descz, desca( m_ ), desca( n_ ), desca( mb_ ),
324  $ desca( nb_ ), desca( rsrc_ ), desca( csrc_ ),
325  $ desca( ctxt_ ), desca( lld_ ), info )
326 *
327  CALL blacs_gridinfo( desca( ctxt_ ), nprow, npcol, myrow, mycol )
328  indiwrk = 1 + iprepad + nprow*npcol + 1
329 *
330  iam = 1
331  IF( myrow.EQ.0 .AND. mycol.EQ.0 )
332  $ iam = 0
333 *
334 * If this process is not involved in this test, bail out now
335 *
336  result = -3
337  IF( myrow.GE.nprow .OR. myrow.LT.0 )
338  $ GO TO 150
339  result = 0
340 *
341 *
342 * DSEED is not used in this call to PDLASIZESYEVX, the
343 * following line just makes ftnchek happy.
344 *
345  dseed( 1 ) = 1
346 *
347  CALL pdlasizesyevx( wknown, range, n, desca, vl, vu, il, iu,
348  $ dseed, win, maxsize, vecsize, valsize )
349 *
350  np = numroc( n, desca( mb_ ), myrow, 0, nprow )
351  nq = numroc( n, desca( nb_ ), mycol, 0, npcol )
352  mq = numroc( maxeigs, desca( nb_ ), mycol, 0, npcol )
353 *
354  CALL dlacpy( 'A', np, nq, copya, desca( lld_ ), a( 1+iprepad ),
355  $ desca( lld_ ) )
356 *
357  CALL dlacpy( 'A', np, nq, copyb, desca( lld_ ), b( 1+iprepad ),
358  $ desca( lld_ ) )
359 *
360  CALL pdfillpad( desca( ctxt_ ), np, nq, b, desca( lld_ ), iprepad,
361  $ ipostpad, padval+1.0d+2 )
362 *
363  CALL pdfillpad( desca( ctxt_ ), np, nq, a, desca( lld_ ), iprepad,
364  $ ipostpad, padval )
365 *
366  CALL pdfillpad( descz( ctxt_ ), np, mq, z, descz( lld_ ), iprepad,
367  $ ipostpad, padval+1.0d+0 )
368 *
369  CALL pdfillpad( desca( ctxt_ ), n, 1, wnew, n, iprepad, ipostpad,
370  $ padval+2.0d+0 )
371 *
372  CALL pdfillpad( desca( ctxt_ ), nprow*npcol, 1, gap, nprow*npcol,
373  $ iprepad, ipostpad, padval+3.0d+0 )
374 *
375  CALL pdfillpad( desca( ctxt_ ), lwork1, 1, work, lwork1, iprepad,
376  $ ipostpad, padval+4.0d+0 )
377 *
378  CALL pifillpad( desca( ctxt_ ), liwork, 1, iwork, liwork, iprepad,
379  $ ipostpad, ipadval )
380 *
381  CALL pifillpad( desca( ctxt_ ), n, 1, ifail, n, iprepad, ipostpad,
382  $ ipadval )
383 *
384  CALL pifillpad( desca( ctxt_ ), 2*nprow*npcol, 1, iclustr,
385  $ 2*nprow*npcol, iprepad, ipostpad, ipadval )
386 *
387 * Make sure that PDSYGVX does not cheat (i.e. use answers
388 * already computed.)
389 *
390  DO 60 i = 1, n, 1
391  DO 50 j = 1, maxeigs, 1
392  CALL pdelset( z( 1+iprepad ), i, j, desca, 13.0d+0 )
393  50 CONTINUE
394  60 CONTINUE
395 *
396  orfac = -1.0d+0
397 *
398  CALL slboot
399  CALL sltimer( 1 )
400  CALL sltimer( 6 )
401  CALL pdsygvx( ibtype, jobz, range, uplo, n, a( 1+iprepad ), ia,
402  $ ja, desca, b( 1+iprepad ), ia, ja, desca, vl, vu,
403  $ il, iu, abstol, m, nz, wnew( 1+iprepad ), orfac,
404  $ z( 1+iprepad ), ia, ja, desca, work( 1+iprepad ),
405  $ lwork1, iwork( 1+iprepad ), liwork,
406  $ ifail( 1+iprepad ), iclustr( 1+iprepad ),
407  $ gap( 1+iprepad ), info )
408  CALL sltimer( 6 )
409  CALL sltimer( 1 )
410 *
411  IF( thresh.LE.0 ) THEN
412  result = 0
413  ELSE
414  CALL pdchekpad( desca( ctxt_ ), 'PDSYGVX-B', np, nq, b,
415  $ desca( lld_ ), iprepad, ipostpad,
416  $ padval+1.0d+2 )
417 *
418  CALL pdchekpad( desca( ctxt_ ), 'PDSYGVX-A', np, nq, a,
419  $ desca( lld_ ), iprepad, ipostpad, padval )
420 *
421  CALL pdchekpad( descz( ctxt_ ), 'PDSYGVX-Z', np, mq, z,
422  $ descz( lld_ ), iprepad, ipostpad,
423  $ padval+1.0d+0 )
424 *
425  CALL pdchekpad( desca( ctxt_ ), 'PDSYGVX-WNEW', n, 1, wnew, n,
426  $ iprepad, ipostpad, padval+2.0d+0 )
427 *
428  CALL pdchekpad( desca( ctxt_ ), 'PDSYGVX-GAP', nprow*npcol, 1,
429  $ gap, nprow*npcol, iprepad, ipostpad,
430  $ padval+3.0d+0 )
431 *
432  CALL pdchekpad( desca( ctxt_ ), 'PDSYGVX-WORK', lwork1, 1,
433  $ work, lwork1, iprepad, ipostpad,
434  $ padval+4.0d+0 )
435 *
436  CALL pichekpad( desca( ctxt_ ), 'PDSYGVX-IWORK', liwork, 1,
437  $ iwork, liwork, iprepad, ipostpad, ipadval )
438 *
439  CALL pichekpad( desca( ctxt_ ), 'PDSYGVX-IFAIL', n, 1, ifail,
440  $ n, iprepad, ipostpad, ipadval )
441 *
442  CALL pichekpad( desca( ctxt_ ), 'PDSYGVX-ICLUSTR',
443  $ 2*nprow*npcol, 1, iclustr, 2*nprow*npcol,
444  $ iprepad, ipostpad, ipadval )
445 *
446 *
447 * Since we now know the spectrum, we can potentially reduce MAXSIZE.
448 *
449  IF( lsame( range, 'A' ) ) THEN
450  CALL pdlasizesyevx( .true., range, n, desca, vl, vu, il, iu,
451  $ dseed, wnew( 1+iprepad ), maxsize,
452  $ vecsize, valsize )
453  END IF
454 *
455 *
456 * Check INFO
457 *
458 *
459 * Make sure that all processes return the same value of INFO
460 *
461  itmp( 1 ) = info
462  itmp( 2 ) = info
463 *
464  CALL igamn2d( desca( ctxt_ ), 'a', ' ', 1, 1, itmp, 1, 1, 1,
465  $ -1, -1, 0 )
466  CALL igamx2d( desca( ctxt_ ), 'a', ' ', 1, 1, itmp( 2 ), 1, 1,
467  $ 1, -1, -1, 0 )
468 *
469 *
470  IF( itmp( 1 ).NE.itmp( 2 ) ) THEN
471  IF( iam.EQ.0 )
472  $ WRITE( nout, fmt = * )
473  $ 'Different processes return different INFO'
474  result = 1
475  ELSE IF( mod( info, 2 ).EQ.1 .OR. info.GT.7 .OR. info.LT.0 )
476  $ THEN
477  IF( iam.EQ.0 )
478  $ WRITE( nout, fmt = 9999 )info
479  result = 1
480  ELSE IF( mod( info / 2, 2 ).EQ.1 .AND. lwork1.GE.maxsize ) THEN
481  IF( iam.EQ.0 )
482  $ WRITE( nout, fmt = 9996 )info
483  result = 1
484  ELSE IF( mod( info / 4, 2 ).EQ.1 .AND. lwork1.GE.vecsize ) THEN
485  IF( iam.EQ.0 )
486  $ WRITE( nout, fmt = 9996 )info
487  result = 1
488  END IF
489 *
490 *
491  IF( lsame( jobz, 'V' ) .AND. ( iclustr( 1+iprepad ).NE.
492  $ 0 ) .AND. ( mod( info / 2, 2 ).NE.1 ) ) THEN
493  IF( iam.EQ.0 )
494  $ WRITE( nout, fmt = 9995 )
495  result = 1
496  END IF
497 *
498 * Check M
499 *
500  IF( ( m.LT.0 ) .OR. ( m.GT.n ) ) THEN
501  IF( iam.EQ.0 )
502  $ WRITE( nout, fmt = 9994 )
503  result = 1
504  ELSE IF( lsame( range, 'A' ) .AND. ( m.NE.n ) ) THEN
505  IF( iam.EQ.0 )
506  $ WRITE( nout, fmt = 9993 )
507  result = 1
508  ELSE IF( lsame( range, 'I' ) .AND. ( m.NE.iu-il+1 ) ) THEN
509  IF( iam.EQ.0 )
510  $ WRITE( nout, fmt = 9992 )
511  result = 1
512  ELSE IF( lsame( jobz, 'V' ) .AND.
513  $ ( .NOT.( lsame( range, 'V' ) ) ) .AND. ( m.NE.nz ) )
514  $ THEN
515  IF( iam.EQ.0 )
516  $ WRITE( nout, fmt = 9991 )
517  result = 1
518  END IF
519 *
520 * Check NZ
521 *
522  IF( lsame( jobz, 'V' ) ) THEN
523  IF( lsame( range, 'V' ) ) THEN
524  IF( nz.GT.m ) THEN
525  IF( iam.EQ.0 )
526  $ WRITE( nout, fmt = 9990 )
527  result = 1
528  END IF
529  IF( nz.LT.m .AND. mod( info / 4, 2 ).NE.1 ) THEN
530  IF( iam.EQ.0 )
531  $ WRITE( nout, fmt = 9989 )
532  result = 1
533  END IF
534  ELSE
535  IF( nz.NE.m ) THEN
536  IF( iam.EQ.0 )
537  $ WRITE( nout, fmt = 9988 )
538  result = 1
539  END IF
540  END IF
541  END IF
542  IF( result.EQ.0 ) THEN
543 *
544 * Make sure that all processes return the same # of eigenvalues
545 *
546  itmp( 1 ) = m
547  itmp( 2 ) = m
548 *
549  CALL igamn2d( desca( ctxt_ ), 'a', ' ', 1, 1, itmp, 1, 1, 1,
550  $ -1, -1, 0 )
551  CALL igamx2d( desca( ctxt_ ), 'a', ' ', 1, 1, itmp( 2 ), 1,
552  $ 1, 1, -1, -1, 0 )
553 *
554  IF( itmp( 1 ).NE.itmp( 2 ) ) THEN
555  IF( iam.EQ.0 )
556  $ WRITE( nout, fmt = 9987 )
557  result = 1
558  ELSE
559 *
560 * Make sure that different processes return the same eigenvalues
561 *
562  DO 70 i = 1, m
563  work( i ) = wnew( i+iprepad )
564  work( i+m ) = wnew( i+iprepad )
565  70 CONTINUE
566 *
567  CALL dgamn2d( desca( ctxt_ ), 'a', ' ', m, 1, work, m, 1,
568  $ 1, -1, -1, 0 )
569  CALL dgamx2d( desca( ctxt_ ), 'a', ' ', m, 1,
570  $ work( 1+m ), m, 1, 1, -1, -1, 0 )
571 *
572  DO 80 i = 1, m
573 *
574  IF( result.EQ.0 .AND. ( abs( work( i )-work( m+
575  $ i ) ).GT.five*eps*abs( work( i ) ) ) ) THEN
576  IF( iam.EQ.0 )
577  $ WRITE( nout, fmt = 9986 )
578  result = 1
579  END IF
580  80 CONTINUE
581  END IF
582  END IF
583 *
584 * Make sure that all processes return the same # of clusters
585 *
586  IF( lsame( jobz, 'V' ) ) THEN
587  nclusters = 0
588  DO 90 i = 0, nprow*npcol - 1
589  IF( iclustr( 1+iprepad+2*i ).EQ.0 )
590  $ GO TO 100
591  nclusters = nclusters + 1
592  90 CONTINUE
593  100 CONTINUE
594  itmp( 1 ) = nclusters
595  itmp( 2 ) = nclusters
596 *
597  CALL igamn2d( desca( ctxt_ ), 'a', ' ', 1, 1, itmp, 1, 1, 1,
598  $ -1, -1, 0 )
599  CALL igamx2d( desca( ctxt_ ), 'a', ' ', 1, 1, itmp( 2 ), 1,
600  $ 1, 1, -1, -1, 0 )
601 *
602  IF( itmp( 1 ).NE.itmp( 2 ) ) THEN
603  IF( iam.EQ.0 )
604  $ WRITE( nout, fmt = 9985 )
605  result = 1
606  ELSE
607 *
608 * Make sure that different processes return the same clusters
609 *
610  DO 110 i = 1, nclusters
611  iwork( indiwrk+i ) = iclustr( i+iprepad )
612  iwork( indiwrk+i+nclusters ) = iclustr( i+iprepad )
613  110 CONTINUE
614  CALL igamn2d( desca( ctxt_ ), 'a', ' ', nclusters*2+1, 1,
615  $ iwork( indiwrk+1 ), nclusters*2+1, 1, 1,
616  $ -1, -1, 0 )
617  CALL igamx2d( desca( ctxt_ ), 'a', ' ', nclusters*2+1, 1,
618  $ iwork( indiwrk+1+nclusters ),
619  $ nclusters*2+1, 1, 1, -1, -1, 0 )
620 *
621 *
622  DO 120 i = 1, nclusters
623  IF( result.EQ.0 .AND. iwork( indiwrk+i ).NE.
624  $ iwork( indiwrk+nclusters+i ) ) THEN
625  IF( iam.EQ.0 )
626  $ WRITE( nout, fmt = 9984 )
627  result = 1
628  END IF
629  120 CONTINUE
630 *
631  IF( iclustr( 1+iprepad+nclusters*2 ).NE.0 ) THEN
632  IF( iam.EQ.0 )
633  $ WRITE( nout, fmt = 9983 )
634  result = 1
635  END IF
636  END IF
637  END IF
638 *
639 *
640  CALL igamx2d( desca( ctxt_ ), 'a', ' ', 1, 1, result, 1, 1, 1,
641  $ -1, -1, 0 )
642  IF( result.NE.0 )
643  $ GO TO 150
644 *
645 * Note that a couple key variables get redefined in PDGSEPCHK
646 * as described by this table:
647 *
648 * PDGSEPTST name PDGSEPCHK name
649 * ------------- -------------
650 * COPYA A
651 * Z Q
652 * B B
653 * A C
654 *
655 *
656  IF( lsame( jobz, 'V' ) ) THEN
657 *
658 * Perform the residual check
659 *
660  CALL pdfillpad( desca( ctxt_ ), sizechk, 1, work, sizechk,
661  $ iprepad, ipostpad, 4.3d+0 )
662 *
663  CALL pdgsepchk( ibtype, n, nz, copya, ia, ja, desca, copyb,
664  $ ia, ja, desca, thresh, z( 1+iprepad ), ia,
665  $ ja, descz, a( 1+iprepad ), ia, ja, desca,
666  $ wnew( 1+iprepad ), work( 1+iprepad ),
667  $ sizechk, tstnrm, res )
668 *
669  CALL pdchekpad( desca( ctxt_ ), 'PDGSEPCHK-WORK', sizechk,
670  $ 1, work, sizechk, iprepad, ipostpad,
671  $ 4.3d+0 )
672 *
673  IF( res.NE.0 )
674  $ result = 1
675  END IF
676 *
677 * Check to make sure that we have the right eigenvalues
678 *
679  IF( wknown ) THEN
680 *
681 * Set up MYIL if necessary
682 *
683  myil = il
684 *
685  IF( lsame( range, 'V' ) ) THEN
686  myil = 1
687  minil = 1
688  maxil = n - m + 1
689  ELSE
690  IF( lsame( range, 'A' ) ) THEN
691  myil = 1
692  END IF
693  minil = myil
694  maxil = myil
695  END IF
696 *
697 * Find the largest difference between the computed
698 * and expected eigenvalues
699 *
700  minerror = normwin
701 *
702  DO 140 myil = minil, maxil
703  maxerror = 0
704 *
705 * Make sure that we aren't skipping any important eigenvalues
706 *
707  misssmallest = .true.
708  IF( .NOT.lsame( range, 'V' ) .OR. ( myil.EQ.1 ) )
709  $ misssmallest = .false.
710  IF( misssmallest .AND. ( win( myil-1 ).LT.vl+normwin*
711  $ five*thresh*eps ) )misssmallest = .false.
712  misslargest = .true.
713  IF( .NOT.lsame( range, 'V' ) .OR. ( myil.EQ.maxil ) )
714  $ misslargest = .false.
715  IF( misslargest .AND. ( win( myil+m ).GT.vu-normwin*five*
716  $ thresh*eps ) )misslargest = .false.
717  IF( .NOT.misssmallest ) THEN
718  IF( .NOT.misslargest ) THEN
719 *
720 * Make sure that the eigenvalues that we report are OK
721 *
722  DO 130 i = 1, m
723  error = abs( win( i+myil-1 )-wnew( i+iprepad ) )
724  maxerror = max( maxerror, error )
725  130 CONTINUE
726 *
727  minerror = min( maxerror, minerror )
728  END IF
729  END IF
730  140 CONTINUE
731 *
732 *
733 * If JOBZ = 'V' and RANGE='A', we might be comparing
734 * against our estimate of what the eigenvalues ought to
735 * be, rather than comparing against what PxSYGVX computed
736 * last time around, so we have to be more generous.
737 *
738  IF( lsame( jobz, 'V' ) .AND. lsame( range, 'A' ) ) THEN
739  IF( minerror.GT.normwin*five*five*thresh*eps ) THEN
740  IF( iam.EQ.0 )
741  $ WRITE( nout, fmt = 9997 )minerror, normwin
742  result = 1
743  END IF
744  ELSE
745  IF( minerror.GT.normwin*five*thresh*eps ) THEN
746  IF( iam.EQ.0 )
747  $ WRITE( nout, fmt = 9997 )minerror, normwin
748  result = 1
749  END IF
750  END IF
751  END IF
752 *
753 *
754 * Make sure that the IL, IU, VL and VU were not altered
755 *
756  IF( il.NE.oldil .OR. iu.NE.oldiu .OR. vl.NE.oldvl .OR. vu.NE.
757  $ oldvu ) THEN
758  IF( iam.EQ.0 )
759  $ WRITE( nout, fmt = 9982 )
760  result = 1
761  END IF
762 *
763  IF( lsame( jobz, 'N' ) .AND. ( nz.NE.oldnz ) ) THEN
764  IF( iam.EQ.0 )
765  $ WRITE( nout, fmt = 9981 )
766  result = 1
767  END IF
768 *
769  END IF
770 *
771 * All processes should report the same result
772 *
773  CALL igamx2d( desca( ctxt_ ), 'a', ' ', 1, 1, result, 1, 1, 1, -1,
774  $ -1, 0 )
775 *
776  150 CONTINUE
777 *
778 *
779  RETURN
780 *
781  9999 FORMAT( 'PDSYGVX returned INFO=', i7 )
782  9998 FORMAT( 'PDSEPQTQ returned INFO=', i7 )
783  9997 FORMAT( 'PDGSEPSUBTST minerror =', d11.2, ' normwin=', d11.2 )
784  9996 FORMAT( 'PDSYGVX returned INFO=', i7,
785  $ ' despite adequate workspace' )
786  9995 FORMAT( .NE..NE.'ICLUSTR(1)0 but mod(INFO/2,2)1' )
787  9994 FORMAT( 'M not in the range 0 to N' )
788  9993 FORMAT( 'M not equal to N' )
789  9992 FORMAT( 'M not equal to IU-IL+1' )
790  9991 FORMAT( 'M not equal to NZ' )
791  9990 FORMAT( 'NZ > M' )
792  9989 FORMAT( 'NZ < M' )
793  9988 FORMAT( 'NZ not equal to M' )
794  9987 FORMAT( 'Different processes return different values for M' )
795  9986 FORMAT( 'Different processes return different eigenvalues' )
796  9985 FORMAT( 'Different processes return ',
797  $ 'different numbers of clusters' )
798  9984 FORMAT( 'Different processes return different clusters' )
799  9983 FORMAT( 'ICLUSTR not zero terminated' )
800  9982 FORMAT( 'IL, IU, VL or VU altered by PDSYGVX' )
801  9981 FORMAT( 'NZ altered by PDSYGVX with JOBZ=N' )
802 *
803 * End of PDGSEPSUBTST
804 *
805  END
pichekpad
subroutine pichekpad(ICTXT, MESS, M, N, A, LDA, IPRE, IPOST, CHKVAL)
Definition: pichekpad.f:3
max
#define max(A, B)
Definition: pcgemr.c:180
pdlasizesyevx
subroutine pdlasizesyevx(WKNOWN, RANGE, N, DESCA, VL, VU, IL, IU, ISEED, WIN, MAXSIZE, VECSIZE, VALSIZE)
Definition: pdlasizesyevx.f:5
pdgsepsubtst
subroutine pdgsepsubtst(WKNOWN, IBTYPE, JOBZ, RANGE, UPLO, N, VL, VU, IL, IU, THRESH, ABSTOL, A, COPYA, B, COPYB, Z, IA, JA, DESCA, WIN, WNEW, IFAIL, ICLUSTR, GAP, IPREPAD, IPOSTPAD, WORK, LWORK, LWORK1, IWORK, LIWORK, RESULT, TSTNRM, QTQNRM, NOUT)
Definition: pdgsepsubtst.f:9
pdchekpad
subroutine pdchekpad(ICTXT, MESS, M, N, A, LDA, IPRE, IPOST, CHKVAL)
Definition: pdchekpad.f:3
sltimer
subroutine sltimer(I)
Definition: sltimer.f:47
pdgsepchk
subroutine pdgsepchk(IBTYPE, MS, NV, A, IA, JA, DESCA, B, IB, JB, DESCB, THRESH, Q, IQ, JQ, DESCQ, C, IC, JC, DESCC, W, WORK, LWORK, TSTNRM, RESULT)
Definition: pdgsepchk.f:6
descinit
subroutine descinit(DESC, M, N, MB, NB, IRSRC, ICSRC, ICTXT, LLD, INFO)
Definition: descinit.f:3
slboot
subroutine slboot()
Definition: sltimer.f:2
pdfillpad
subroutine pdfillpad(ICTXT, M, N, A, LDA, IPRE, IPOST, CHKVAL)
Definition: pdfillpad.f:2
pdlasizegsep
subroutine pdlasizegsep(DESCA, IPREPAD, IPOSTPAD, SIZEMQRLEFT, SIZEMQRRIGHT, SIZEQRF, SIZETMS, SIZEQTQ, SIZECHK, SIZESYEVX, ISIZESYEVX, SIZESUBTST, ISIZESUBTST, SIZETST, ISIZETST)
Definition: pdlasizegsep.f:8
pdelset
subroutine pdelset(A, IA, JA, DESCA, ALPHA)
Definition: pdelset.f:2
pdsygvx
subroutine pdsygvx(IBTYPE, JOBZ, RANGE, UPLO, N, A, IA, JA, DESCA, B, IB, JB, DESCB, VL, VU, IL, IU, ABSTOL, M, NZ, W, ORFAC, Z, IZ, JZ, DESCZ, WORK, LWORK, IWORK, LIWORK, IFAIL, ICLUSTR, GAP, INFO)
Definition: pdsygvx.f:6
pifillpad
subroutine pifillpad(ICTXT, M, N, A, LDA, IPRE, IPOST, CHKVAL)
Definition: pifillpad.f:2
min
#define min(A, B)
Definition: pcgemr.c:181