59 INTEGER block_cyclic_2d, csrc_, ctxt_, dlen_, dt_,
60 $ lld_, mb_, m_, nb_, n_, rsrc_
61 parameter( block_cyclic_2d = 1, dlen_ = 9, dt_ = 1,
62 $ ctxt_ = 2, m_ = 3, n_ = 4, mb_ = 5, nb_ = 6,
63 $ rsrc_ = 7, csrc_ = 8, lld_ = 9 )
64 INTEGER realsz, totmem, memsiz, ntests
65 REAL padval, zero, one
66 parameter( realsz = 4, totmem = 2000000,
67 $ memsiz = totmem / realsz, ntests = 20,
68 $ padval = -9923.0e+0, zero = 0.0e+0,
75 INTEGER i, iam, iaseed, ictxt, iii, imidpad, info, ipa,
76 $ ipostpad, iprepad, ipw, ipwi, ipwr, ipz, j, k,
77 $ kfail, kpass, kskip, ktests, lda, ldz, lwork,
78 $ mycol, myrow, n, nb, ngrids, nmat, nnb, nout,
79 $ np, npcol, nprocs, nprow, nq, worksiz
81 REAL anorm, fresid, qresid, znorm
82 DOUBLE PRECISION nops, tmflops
85 INTEGER desca( dlen_ ), descz( dlen_ ), ierr( 2 ),
86 $ idum( 1 ), nbval( ntests ), nval( ntests ),
87 $ pval( ntests ), qval( ntests )
89 DOUBLE PRECISION ctime( 1 ), wtime( 1 )
92 EXTERNAL blacs_barrier, blacs_exit, blacs_get,
93 $ blacs_gridexit, blacs_gridinfo, blacs_gridinit,
105 INTRINSIC dble,
max,
min, real
108 DATA kfail, kpass, kskip, ktests / 4*0 /
114 CALL blacs_pinfo( iam, nprocs )
116 CALL psnepinfo( outfile, nout, nmat, nval, ntests, nnb, nbval,
117 $ ntests, ngrids, pval, ntests, qval, ntests,
118 $ thresh, mem, iam, nprocs )
119 check = ( thresh.GE.0.0e+0 )
124 WRITE( nout, fmt = * )
125 WRITE( nout, fmt = 9995 )
126 WRITE( nout, fmt = 9994 )
127 WRITE( nout, fmt = * )
140 IF( nprow.LT.1 )
THEN
142 $
WRITE( nout, fmt = 9999 )
'GRID',
'nprow', nprow
144 ELSE IF( npcol.LT.1 )
THEN
146 $
WRITE( nout, fmt = 9999 )
'GRID',
'npcol', npcol
148 ELSE IF( nprow*npcol.GT.nprocs )
THEN
150 $
WRITE( nout, fmt = 9998 )nprow*npcol, nprocs
154 IF( ierr( 1 ).GT.0 )
THEN
156 $
WRITE( nout, fmt = 9997 )
'grid'
163 CALL blacs_get( -1, 0, ictxt )
164 CALL blacs_gridinit( ictxt,
'Row-major', nprow, npcol )
165 CALL blacs_gridinfo( ictxt, nprow, npcol, myrow, mycol )
170 IF( myrow.GE.nprow .OR. mycol.GE.npcol )
182 $
WRITE( nout, fmt = 9999 )
'MATRIX',
'N', n
188 CALL igsum2d( ictxt,
'All',
' ', 1, 1, ierr, 1, -1, 0 )
190 IF( ierr( 1 ).GT.0 )
THEN
192 $
WRITE( nout, fmt = 9997 )
'matrix'
207 $
WRITE( nout, fmt = 9999 )
'NB',
'NB', nb
212 CALL igsum2d( ictxt,
'All',
' ', 1, 1, ierr, 1, -1, 0 )
214 IF( ierr( 1 ).GT.0 )
THEN
216 $
WRITE( nout, fmt = 9997 )
'NB'
223 np =
numroc( n, nb, myrow, 0, nprow )
224 nq =
numroc( n, nb, mycol, 0, npcol )
226 iprepad =
max( nb, np )
228 ipostpad =
max( nb, nq )
229 iprepad = iprepad + 1000
230 imidpad = imidpad + 1000
231 ipostpad = ipostpad + 1000
240 CALL descinit( desca, n, n, nb, nb, 0, 0, ictxt,
241 $
max( 1, np )+imidpad, ierr( 1 ) )
245 CALL descinit( descz, n, n, nb, nb, 0, 0, ictxt,
246 $
max( 1, np )+imidpad, ierr( 2 ) )
253 CALL igsum2d( ictxt,
'All',
' ', 2, 1, ierr, 2, -1, 0 )
255 IF( ierr( 1 ).LT.0 .OR. ierr( 2 ).LT.0 )
THEN
257 $
WRITE( nout, fmt = 9997 )
'descriptor'
266 ipz = ipa + desca( lld_ )*nq + ipostpad + iprepad
267 ipwr = ipz + descz( lld_ )*nq + ipostpad + iprepad
268 ipwi = ipwr + n + ipostpad + iprepad
269 ipw = ipwi + n + ipostpad + iprepad
273 iii = 7*iii /
ilcm( nprow, npcol )
276 lwork = 3*n +
max( 2*
max( lda, ldz )+2*nq, iii )
277 lwork = lwork +
max(2*n, (8*
ilcm(nprow,npcol)+2)**2 )
284 worksiz = lwork +
max( np*desca( nb_ ),
285 $ desca( mb_ )*nq ) + ipostpad
289 worksiz = lwork + ipostpad
296 IF( ipw+worksiz.GT.memsiz )
THEN
298 $
WRITE( nout, fmt = 9996 )
'Schur reduction',
299 $ ( ipw+worksiz )*realsz
305 CALL igsum2d( ictxt,
'All',
' ', 1, 1, ierr, 1, -1, 0 )
307 IF( ierr( 1 ).GT.0 )
THEN
309 $
WRITE( nout, fmt = 9997 )
'MEMORY'
316 CALL pslaset(
'All', n, n, zero, one, mem( ipz ), 1, 1,
321 CALL psmatgen( ictxt,
'No transpose',
'No transpose',
322 $ desca( m_ ), desca( n_ ), desca( mb_ ),
323 $ desca( nb_ ), mem( ipa ), desca( lld_ ),
324 $ desca( rsrc_ ), desca( csrc_ ), iaseed, 0,
325 $ np, 0, nq, myrow, mycol, nprow, npcol )
327 $ zero, zero, mem( ipa ),
min( n, 3 ), 1,
333 CALL psfillpad( ictxt, np, nq, mem( ipa-iprepad ),
334 $ desca( lld_ ), iprepad, ipostpad,
336 CALL psfillpad( ictxt, np, nq, mem( ipz-iprepad ),
337 $ descz( lld_ ), iprepad, ipostpad,
339 CALL psfillpad( ictxt, worksiz-ipostpad, 1,
340 $ mem( ipw-iprepad ), worksiz-ipostpad,
341 $ iprepad, ipostpad, padval )
342 anorm =
pslanhs(
'I', n, mem( ipa ), 1, 1, desca,
344 CALL pschekpad( ictxt,
'PSLANHS', np, nq,
345 $ mem( ipa-iprepad ), desca( lld_ ),
346 $ iprepad, ipostpad, padval )
347 CALL pschekpad( ictxt,
'PSLANHS', worksiz-ipostpad, 1,
348 $ mem( ipw-iprepad ), worksiz-ipostpad,
349 $ iprepad, ipostpad, padval )
351 CALL psfillpad( ictxt, n, 1, mem( ipwr-iprepad ), n,
352 $ iprepad, ipostpad, padval )
353 CALL psfillpad( ictxt, n, 1, mem( ipwi-iprepad ), n,
354 $ iprepad, ipostpad, padval )
355 CALL psfillpad( ictxt, lwork, 1, mem( ipw-iprepad ),
356 $ lwork, iprepad, ipostpad, padval )
361 CALL blacs_barrier( ictxt,
'All' )
366 CALL pslahqr( .true., .true., n, 1, n, mem( ipa ), desca,
367 $ mem( ipwr ), mem( ipwi ), 1, n, mem( ipz ),
368 $ descz, mem( ipw ), lwork, idum, 0, info )
374 $
WRITE( nout, fmt = * )
'PSLAHQR INFO=', info
383 CALL pschekpad( ictxt,
'PSLAHQR (A)', np, nq,
384 $ mem( ipa-iprepad ), desca( lld_ ),
385 $ iprepad, ipostpad, padval )
386 CALL pschekpad( ictxt,
'PSLAHQR (Z)', np, nq,
387 $ mem( ipz-iprepad ), descz( lld_ ),
388 $ iprepad, ipostpad, padval )
389 CALL pschekpad( ictxt,
'PSLAHQR (WR)', n, 1,
390 $ mem( ipwr-iprepad ), n, iprepad,
392 CALL pschekpad( ictxt,
'PSLAHQR (WI)', n, 1,
393 $ mem( ipwi-iprepad ), n, iprepad,
395 CALL pschekpad( ictxt,
'PSLAHQR (WORK)', lwork, 1,
396 $ mem( ipw-iprepad ), lwork, iprepad,
399 CALL psfillpad( ictxt, worksiz-ipostpad, 1,
400 $ mem( ipw-iprepad ), worksiz-ipostpad,
401 $ iprepad, ipostpad, padval )
405 CALL psnepfchk( n, mem( ipa ), 1, 1, desca, iaseed,
406 $ mem( ipz ), 1, 1, descz, anorm,
407 $ fresid, mem( ipw ) )
409 CALL pschekpad( ictxt,
'PSNEPFCHK (A)', np, nq,
410 $ mem( ipa-iprepad ), desca( lld_ ),
411 $ iprepad, ipostpad, padval )
412 CALL pschekpad( ictxt,
'PSNEPFCHK (Z)', np, nq,
413 $ mem( ipz-iprepad ), descz( lld_ ),
414 $ iprepad, ipostpad, padval )
415 CALL pschekpad( ictxt,
'PSNEPFCHK (WORK)',
416 $ worksiz-ipostpad, 1,
417 $ mem( ipw-iprepad ), worksiz-ipostpad,
418 $ iprepad, ipostpad, padval )
422 CALL pslaset(
'All', n, n, zero, one, mem( ipa ), 1,
424 CALL psgemm(
'Transpose',
'No transpose', n, n, n,
425 $ -one, mem( ipz ), 1, 1, descz,
426 $ mem( ipz ), 1, 1, descz, one, mem( ipa ),
428 znorm =
pslange(
'1', n, n, mem( ipa ), 1, 1, desca,
430 qresid = znorm / ( real( n )*
pslamch( ictxt,
'P' ) )
434 IF( ( fresid.LE.thresh ) .AND.
435 $ ( ( fresid-fresid ).EQ.0.0e+0 ) .AND.
436 $ ( qresid.LE.thresh ) .AND.
437 $ ( ( qresid-qresid ).EQ.0.0e+0 ) )
THEN
444 WRITE( nout, fmt = 9986 )fresid
445 WRITE( nout, fmt = 9985 )qresid
454 fresid = fresid - fresid
455 qresid = qresid - qresid
462 CALL slcombine( ictxt,
'All',
'>',
'W', 1, 1, wtime )
463 CALL slcombine( ictxt,
'All',
'>',
'C', 1, 1, ctime )
467 IF( myrow.EQ.0 .AND. mycol.EQ.0 )
THEN
471 nops = 18.0d+0*dble( n )**3
478 IF( wtime( 1 ).GT.0.0d+0 )
THEN
479 tmflops = nops / ( wtime( 1 )*1.0d+6 )
483 IF( wtime( 1 ).GE.0.0d+0 )
484 $
WRITE( nout, fmt = 9993 )
'WALL', n, nb, nprow,
485 $ npcol, wtime( 1 ), tmflops, passed
489 IF( ctime( 1 ).GT.0.0d+0 )
THEN
490 tmflops = nops / ( ctime( 1 )*1.0d+6 )
495 IF( ctime( 1 ).GE.0.0d+0 )
496 $
WRITE( nout, fmt = 9993 )
'CPU ', n, nb, nprow,
497 $ npcol, ctime( 1 ), tmflops, passed
504 CALL blacs_gridexit( ictxt )
511 ktests = kpass + kfail + kskip
512 WRITE( nout, fmt = * )
513 WRITE( nout, fmt = 9992 )ktests
515 WRITE( nout, fmt = 9991 )kpass
516 WRITE( nout, fmt = 9989 )kfail
518 WRITE( nout, fmt = 9990 )kpass
520 WRITE( nout, fmt = 9988 )kskip
521 WRITE( nout, fmt = * )
522 WRITE( nout, fmt = * )
523 WRITE( nout, fmt = 9987 )
524 IF( nout.NE.6 .AND. nout.NE.0 )
530 9999
FORMAT(
'ILLEGAL ', a6,
': ', a5,
' = ', i3,
531 $
'; It should be at least 1' )
532 9998
FORMAT(
'ILLEGAL GRID: nprow*npcol = ', i4,
'. It can be at most',
534 9997
FORMAT(
'Bad ', a6,
' parameters: going on to next test case.' )
535 9996
FORMAT(
'Unable to perform ', a,
': need TOTMEM of at least',
537 9995
FORMAT(
'TIME N NB P Q NEP Time MFLOPS CHECK' )
538 9994
FORMAT(
'---- ----- --- ---- ---- -------- -------- ------' )
539 9993
FORMAT( a4, 1x, i5, 1x, i3, 1x, i4, 1x, i4, 1x, f8.2, 1x, f8.2,
541 9992
FORMAT(
'Finished ', i6,
' tests, with the following results:' )
542 9991
FORMAT( i5,
' tests completed and passed residual checks.' )
543 9990
FORMAT( i5,
' tests completed without checking.' )
544 9989
FORMAT( i5,
' tests completed and failed residual checks.' )
545 9988
FORMAT( i5,
' tests skipped because of illegal input values.' )
546 9987
FORMAT(
'END OF TESTS.' )
547 9986
FORMAT(
'||H - Q*S*Q^T|| / (||H|| * N * eps) = ', g25.7 )
548 9985
FORMAT(
'||Q^T*Q - I|| / ( N * eps ) ', g25.7 )