SCALAPACK 2.2.2
LAPACK: Linear Algebra PACKage
Loading...
Searching...
No Matches

◆ pzsepsubtst()

subroutine pzsepsubtst ( logical  wknown,
character  jobz,
character  range,
character  uplo,
integer  n,
double precision  vl,
double precision  vu,
integer  il,
integer  iu,
double precision  thresh,
double precision  abstol,
complex*16, dimension( * )  a,
complex*16, dimension( * )  copya,
complex*16, dimension( * )  z,
integer  ia,
integer  ja,
integer, dimension( * )  desca,
double precision, dimension( * )  win,
double precision, dimension( * )  wnew,
integer, dimension( * )  ifail,
integer, dimension( * )  iclustr,
double precision, dimension( * )  gap,
integer  iprepad,
integer  ipostpad,
complex*16, dimension( * )  work,
integer  lwork,
double precision, dimension( * )  rwork,
integer  lrwork,
integer  lwork1,
integer, dimension( * )  iwork,
integer  liwork,
integer  result,
double precision  tstnrm,
double precision  qtqnrm,
integer  nout 
)

Definition at line 3 of file pzsepsubtst.f.

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