ScaLAPACK 2.1  2.1
ScaLAPACK: Scalable Linear Algebra PACKage
pdtrddriver.f
Go to the documentation of this file.
1  PROGRAM pdtrddriver
2 *
3 * -- ScaLAPACK testing driver (version 1.7) --
4 * University of Tennessee, Knoxville, Oak Ridge National Laboratory,
5 * and University of California, Berkeley.
6 * October 15, 1999
7 *
8 * Purpose
9 * ========
10 *
11 * PDTRDDRIVER is the main test program for the DOUBLE PRECISION
12 * SCALAPACK TRD (symmetric tridiagonal reduction) routines.
13 *
14 * The program must be driven by a short data file. An annotated
15 * example of a data file can be obtained by deleting the first 3
16 * characters from the following 13 lines:
17 * 'ScaLAPACK TRD computation input file'
18 * 'PVM machine'
19 * 'TRD.out' output file name
20 * 6 device out
21 * 'L' define Lower or Upper
22 * 3 number of problems sizes
23 * 5 31 201 values of N
24 * 3 number of NB's
25 * 2 3 5 values of NB
26 * 7 number of process grids (ordered pairs of P & Q)
27 * 1 2 1 4 2 3 8 values of P
28 * 1 2 4 1 3 2 1 values of Q
29 * 1.0 threshold
30 *
31 * Internal Parameters
32 * ===================
33 *
34 * TOTMEM INTEGER, default = 2000000
35 * TOTMEM is a machine-specific parameter indicating the
36 * maximum amount of available memory in bytes.
37 * The user should customize TOTMEM to his platform. Remember
38 * to leave room in memory for the operating system, the BLACS
39 * buffer, etc. For example, on a system with 8 MB of memory
40 * per process (e.g., one processor on an Intel iPSC/860), the
41 * parameters we use are TOTMEM=6200000 (leaving 1.8 MB for OS,
42 * code, BLACS buffer, etc). However, for PVM, we usually set
43 * TOTMEM = 2000000. Some experimenting with the maximum value
44 * of TOTMEM may be required.
45 *
46 * INTGSZ INTEGER, default = 4 bytes.
47 * DBLESZ INTEGER, default = 8 bytes.
48 * INTGSZ and DBLESZ indicate the length in bytes on the
49 * given platform for an integer and a double precision real.
50 * MEM DOUBLE PRECISION array, dimension ( TOTMEM / DBLESZ )
51 *
52 * All arrays used by SCALAPACK routines are allocated from
53 * this array and referenced by pointers. The integer IPA,
54 * for example, is a pointer to the starting element of MEM for
55 * the matrix A.
56 *
57 * =====================================================================
58 *
59 * .. Parameters ..
60  INTEGER block_cyclic_2d, csrc_, ctxt_, dlen_, dtype_,
61  $ lld_, mb_, m_, nb_, n_, rsrc_
62  parameter( block_cyclic_2d = 1, dlen_ = 9, dtype_ = 1,
63  $ ctxt_ = 2, m_ = 3, n_ = 4, mb_ = 5, nb_ = 6,
64  $ rsrc_ = 7, csrc_ = 8, lld_ = 9 )
65  INTEGER dblesz, totmem, memsiz, ntests
66  DOUBLE PRECISION padval
67  parameter( dblesz = 8, totmem = 2000000,
68  $ memsiz = totmem / dblesz, ntests = 20,
69  $ padval = -9923.0d+0 )
70 * ..
71 * .. Local Scalars ..
72  LOGICAL check
73  CHARACTER uplo
74  CHARACTER*6 passed
75  CHARACTER*80 outfile
76  INTEGER i, iam, iaseed, ictxt, imidpad, info, ipa, ipd,
77  $ ipe, ipostpad, iprepad, ipt, ipw, itemp, j, k,
78  $ kfail, kpass, kskip, ktests, lcm, lwork, mycol,
79  $ myrow, n, nb, ndiag, ngrids, nmat, nnb, noffd,
80  $ nout, np, npcol, nprocs, nprow, nq, worksiz,
81  $ worktrd
82  REAL thresh
83  DOUBLE PRECISION anorm, fresid, nops, tmflops
84 * ..
85 * .. Local Arrays ..
86  INTEGER desca( dlen_ ), ierr( 1 ), nbval( ntests ),
87  $ nval( ntests ), pval( ntests ), qval( ntests )
88  DOUBLE PRECISION ctime( 1 ), mem( memsiz ), wtime( 1 )
89 * ..
90 * .. External Subroutines ..
91  EXTERNAL blacs_barrier, blacs_exit, blacs_get,
92  $ blacs_gridexit, blacs_gridinfo, blacs_gridinit,
93  $ blacs_pinfo, descinit, igsum2d, pdchekpad,
97 * ..
98 * .. External Functions ..
99  LOGICAL lsame
100  INTEGER iceil, ilcm, numroc
101  DOUBLE PRECISION pdlansy
102  EXTERNAL lsame, iceil, ilcm, numroc, pdlansy
103 * ..
104 * .. Intrinsic Functions ..
105  INTRINSIC dble, max
106 * ..
107 * .. Data statements ..
108  DATA ktests, kpass, kfail, kskip / 4*0 /
109 * ..
110 * .. Executable Statements ..
111 *
112  IF( block_cyclic_2d*csrc_*ctxt_*dlen_*dtype_*lld_*mb_*m_*nb_*n_*
113  $ rsrc_.LT.0 )stop
114 * Get starting information
115 *
116  CALL blacs_pinfo( iam, nprocs )
117  iaseed = 100
118  CALL pdtrdinfo( outfile, nout, uplo, nmat, nval, ntests, nnb,
119  $ nbval, ntests, ngrids, pval, ntests, qval, ntests,
120  $ thresh, mem, iam, nprocs )
121  check = ( thresh.GE.0.0e+0 )
122 *
123 * Print headings
124 *
125  IF( iam.EQ.0 ) THEN
126  WRITE( nout, fmt = * )
127  WRITE( nout, fmt = 9995 )
128  WRITE( nout, fmt = 9994 )
129  WRITE( nout, fmt = * )
130  END IF
131 *
132 * Loop over different process grids
133 *
134  DO 30 i = 1, ngrids
135 *
136  nprow = pval( i )
137  npcol = qval( i )
138 *
139 * Make sure grid information is correct
140 *
141  ierr( 1 ) = 0
142  IF( nprow.LT.1 ) THEN
143  IF( iam.EQ.0 )
144  $ WRITE( nout, fmt = 9999 )'GRID', 'nprow', nprow
145  ierr( 1 ) = 1
146  ELSE IF( npcol.LT.1 ) THEN
147  IF( iam.EQ.0 )
148  $ WRITE( nout, fmt = 9999 )'GRID', 'npcol', npcol
149  ierr( 1 ) = 1
150  ELSE IF( nprow*npcol.GT.nprocs ) THEN
151  IF( iam.EQ.0 )
152  $ WRITE( nout, fmt = 9998 )nprow*npcol, nprocs
153  ierr( 1 ) = 1
154  END IF
155 *
156  IF( ierr( 1 ).GT.0 ) THEN
157  IF( iam.EQ.0 )
158  $ WRITE( nout, fmt = 9997 )'grid'
159  kskip = kskip + 1
160  GO TO 30
161  END IF
162 *
163 * Define process grid
164 *
165  CALL blacs_get( -1, 0, ictxt )
166  CALL blacs_gridinit( ictxt, 'Row-major', nprow, npcol )
167  CALL blacs_gridinfo( ictxt, nprow, npcol, myrow, mycol )
168 *
169 * Go to bottom of loop if this case doesn't use my process
170 *
171  IF( myrow.GE.nprow .OR. mycol.GE.npcol )
172  $ GO TO 30
173 *
174  DO 20 j = 1, nmat
175 *
176  n = nval( j )
177 *
178 * Make sure matrix information is correct
179 *
180  ierr( 1 ) = 0
181  IF( n.LT.1 ) THEN
182  IF( iam.EQ.0 )
183  $ WRITE( nout, fmt = 9999 )'MATRIX', 'N', n
184  ierr( 1 ) = 1
185  END IF
186 *
187 * Make sure no one had error
188 *
189  CALL igsum2d( ictxt, 'All', ' ', 1, 1, ierr, 1, -1, 0 )
190 *
191  IF( ierr( 1 ).GT.0 ) THEN
192  IF( iam.EQ.0 )
193  $ WRITE( nout, fmt = 9997 )'matrix'
194  kskip = kskip + 1
195  GO TO 20
196  END IF
197 *
198 * Loop over different blocking sizes
199 *
200  DO 10 k = 1, nnb
201 *
202  nb = nbval( k )
203 *
204 * Make sure nb is legal
205 *
206  ierr( 1 ) = 0
207  IF( nb.LT.1 ) THEN
208  ierr( 1 ) = 1
209  IF( iam.EQ.0 )
210  $ WRITE( nout, fmt = 9999 )'NB', 'NB', nb
211  END IF
212 *
213 * Check all processes for an error
214 *
215  CALL igsum2d( ictxt, 'All', ' ', 1, 1, ierr, 1, -1, 0 )
216 *
217  IF( ierr( 1 ).GT.0 ) THEN
218  IF( iam.EQ.0 )
219  $ WRITE( nout, fmt = 9997 )'NB'
220  kskip = kskip + 1
221  GO TO 10
222  END IF
223 *
224 * Padding constants
225 *
226  np = numroc( n, nb, myrow, 0, nprow )
227  nq = numroc( n, nb, mycol, 0, npcol )
228  IF( check ) THEN
229  iprepad = max( nb, np )
230  imidpad = nb
231  ipostpad = max( nb, nq )
232  ELSE
233  iprepad = 0
234  imidpad = 0
235  ipostpad = 0
236  END IF
237 *
238 * Initialize the array descriptor for the matrix A
239 *
240  CALL descinit( desca, n, n, nb, nb, 0, 0, ictxt,
241  $ max( 1, np )+imidpad, ierr( 1 ) )
242 *
243 * Check all processes for an error
244 *
245  CALL igsum2d( ictxt, 'All', ' ', 1, 1, ierr, 1, -1, 0 )
246 *
247  IF( ierr( 1 ).LT.0 ) THEN
248  IF( iam.EQ.0 )
249  $ WRITE( nout, fmt = 9997 )'descriptor'
250  kskip = kskip + 1
251  GO TO 10
252  END IF
253 *
254 * Assign pointers into MEM for SCALAPACK arrays, A is
255 * allocated starting at position MEM( IPREPAD+1 )
256 *
257  ndiag = nq
258  IF( lsame( uplo, 'U' ) ) THEN
259  noffd = nq
260  ELSE
261  noffd = numroc( n-1, nb, mycol, 0, npcol )
262  END IF
263 *
264  ipa = iprepad + 1
265  ipd = ipa + desca( lld_ )*nq + ipostpad + iprepad
266  ipe = ipd + ndiag + ipostpad + iprepad
267  ipt = ipe + noffd + ipostpad + iprepad
268  ipw = ipt + nq + ipostpad + iprepad
269 *
270 * Calculate the amount of workspace required for the
271 * reduction
272 *
273  lwork = max( nb*( np+1 ), 3*nb )
274  worktrd = lwork + ipostpad
275  worksiz = worktrd
276 *
277 * Figure the amount of workspace required by the check
278 *
279  IF( check ) THEN
280  itemp = 2*nq + np
281  IF( nprow.NE.npcol ) THEN
282  lcm = ilcm( nprow, npcol )
283  itemp = nb*iceil( iceil( np, nb ), lcm / nprow ) +
284  $ itemp
285  END IF
286  itemp = max( itemp, 2*( nb+np )*nb )
287  worksiz = max( lwork, itemp ) + ipostpad
288  END IF
289 *
290 * Check for adequate memory for problem size
291 *
292  ierr( 1 ) = 0
293  IF( ipw+worksiz.GT.memsiz ) THEN
294  IF( iam.EQ.0 )
295  $ WRITE( nout, fmt = 9996 )'Tridiagonal reduction',
296  $ ( ipw+worksiz )*dblesz
297  ierr( 1 ) = 1
298  END IF
299 *
300 * Check all processes for an error
301 *
302  CALL igsum2d( ictxt, 'All', ' ', 1, 1, ierr, 1, -1, 0 )
303 *
304  IF( ierr( 1 ).GT.0 ) THEN
305  IF( iam.EQ.0 )
306  $ WRITE( nout, fmt = 9997 )'MEMORY'
307  kskip = kskip + 1
308  GO TO 10
309  END IF
310 *
311 * Generate the matrix A
312 *
313  CALL pdmatgen( ictxt, 'Symm', 'N', desca( m_ ),
314  $ desca( n_ ), desca( mb_ ), desca( nb_ ),
315  $ mem( ipa ), desca( lld_ ), desca( rsrc_ ),
316  $ desca( csrc_ ), iaseed, 0, np, 0, nq,
317  $ myrow, mycol, nprow, npcol )
318 *
319 * Need Infinity-norm of A for checking
320 *
321  IF( check ) THEN
322  CALL pdfillpad( ictxt, np, nq, mem( ipa-iprepad ),
323  $ desca( lld_ ), iprepad, ipostpad,
324  $ padval )
325  CALL pdfillpad( ictxt, ndiag, 1, mem( ipd-iprepad ),
326  $ ndiag, iprepad, ipostpad, padval )
327  CALL pdfillpad( ictxt, noffd, 1, mem( ipe-iprepad ),
328  $ noffd, iprepad, ipostpad, padval )
329  CALL pdfillpad( ictxt, nq, 1, mem( ipt-iprepad ), nq,
330  $ iprepad, ipostpad, padval )
331  CALL pdfillpad( ictxt, worksiz-ipostpad, 1,
332  $ mem( ipw-iprepad ), worksiz-ipostpad,
333  $ iprepad, ipostpad, padval )
334  anorm = pdlansy( 'I', uplo, n, mem( ipa ), 1, 1,
335  $ desca, mem( ipw ) )
336  CALL pdchekpad( ictxt, 'PDLANSY', np, nq,
337  $ mem( ipa-iprepad ), desca( lld_ ),
338  $ iprepad, ipostpad, padval )
339  CALL pdchekpad( ictxt, 'PDLANSY', worksiz-ipostpad, 1,
340  $ mem( ipw-iprepad ), worksiz-ipostpad,
341  $ iprepad, ipostpad, padval )
342  CALL pdfillpad( ictxt, worktrd-ipostpad, 1,
343  $ mem( ipw-iprepad ), worktrd-ipostpad,
344  $ iprepad, ipostpad, padval )
345  END IF
346 *
347  CALL slboot
348  CALL blacs_barrier( ictxt, 'All' )
349  CALL sltimer( 1 )
350 *
351 * Reduce to symmetric tridiagonal form
352 *
353  CALL pdsytrd( uplo, n, mem( ipa ), 1, 1, desca,
354  $ mem( ipd ), mem( ipe ), mem( ipt ),
355  $ mem( ipw ), lwork, info )
356 *
357  CALL sltimer( 1 )
358 *
359  IF( check ) THEN
360 *
361 * Check for memory overwrite
362 *
363  CALL pdchekpad( ictxt, 'PDSYTRD', np, nq,
364  $ mem( ipa-iprepad ), desca( lld_ ),
365  $ iprepad, ipostpad, padval )
366  CALL pdchekpad( ictxt, 'PDSYTRD', ndiag, 1,
367  $ mem( ipd-iprepad ), ndiag, iprepad,
368  $ ipostpad, padval )
369  CALL pdchekpad( ictxt, 'PDSYTRD', noffd, 1,
370  $ mem( ipe-iprepad ), noffd, iprepad,
371  $ ipostpad, padval )
372  CALL pdchekpad( ictxt, 'PDSYTRD', nq, 1,
373  $ mem( ipt-iprepad ), nq, iprepad,
374  $ ipostpad, padval )
375  CALL pdchekpad( ictxt, 'PDSYTRD', worktrd-ipostpad, 1,
376  $ mem( ipw-iprepad ), worktrd-ipostpad,
377  $ iprepad, ipostpad, padval )
378  CALL pdfillpad( ictxt, worksiz-ipostpad, 1,
379  $ mem( ipw-iprepad ), worksiz-ipostpad,
380  $ iprepad, ipostpad, padval )
381 *
382 * Compute fctres = ||A - QTQ'|| / (||A|| * N * eps)
383 *
384  CALL pdsytdrv( uplo, n, mem( ipa ), 1, 1, desca,
385  $ mem( ipd ), mem( ipe ), mem( ipt ),
386  $ mem( ipw ), ierr( 1 ) )
387  CALL pdlafchk( 'Symm', 'No', n, n, mem( ipa ), 1, 1,
388  $ desca, iaseed, anorm, fresid,
389  $ mem( ipw ) )
390 *
391 * Check for memory overwrite
392 *
393  CALL pdchekpad( ictxt, 'PDSYTDRV', np, nq,
394  $ mem( ipa-iprepad ), desca( lld_ ),
395  $ iprepad, ipostpad, padval )
396  CALL pdchekpad( ictxt, 'PDSYTDRV', ndiag, 1,
397  $ mem( ipd-iprepad ), ndiag, iprepad,
398  $ ipostpad, padval )
399  CALL pdchekpad( ictxt, 'PDSYTDRV', noffd, 1,
400  $ mem( ipe-iprepad ), noffd, iprepad,
401  $ ipostpad, padval )
402  CALL pdchekpad( ictxt, 'PDSYTDRV', worksiz-ipostpad,
403  $ 1, mem( ipw-iprepad ),
404  $ worksiz-ipostpad, iprepad, ipostpad,
405  $ padval )
406 *
407 * Test residual and detect NaN result
408 *
409  IF( fresid.LE.thresh .AND. fresid-fresid.EQ.
410  $ 0.0d+0 .AND. ierr( 1 ).EQ.0 ) THEN
411  kpass = kpass + 1
412  passed = 'PASSED'
413  ELSE
414  IF( myrow.EQ.0 .AND. mycol.EQ.0 )
415  $ WRITE( nout, fmt = 9986 )fresid
416  kfail = kfail + 1
417  passed = 'FAILED'
418  END IF
419 *
420  IF( myrow.EQ.0 .AND. mycol.EQ.0 .AND. ierr( 1 ).NE.0 )
421  $ WRITE( nout, fmt = * )'D or E copies incorrect ...'
422  ELSE
423 *
424 * Don't perform the checking, only the timing operation
425 *
426  kpass = kpass + 1
427  fresid = fresid - fresid
428  passed = 'BYPASS'
429  END IF
430 *
431 * Gather maximum of all CPU and WALL clock timings
432 *
433  CALL slcombine( ictxt, 'All', '>', 'W', 1, 1, wtime )
434  CALL slcombine( ictxt, 'All', '>', 'C', 1, 1, ctime )
435 *
436 * Print results
437 *
438  IF( myrow.EQ.0 .AND. mycol.EQ.0 ) THEN
439 *
440 * TRD requires 4/3 N^3 floating point operations
441 *
442  nops = dble( n )
443 *
444  nops = ( 4.0d+0 / 3.0d+0 )*nops**3
445  nops = nops / 1.0d+6
446 *
447 * Print WALL time
448 *
449  IF( wtime( 1 ).GT.0.0d+0 ) THEN
450  tmflops = nops / wtime( 1 )
451  ELSE
452  tmflops = 0.0d+0
453  END IF
454  IF( wtime( 1 ).GE.0.0d+0 )
455  $ WRITE( nout, fmt = 9993 )'WALL', uplo, n, nb,
456  $ nprow, npcol, wtime( 1 ), tmflops, fresid, passed
457 *
458 * Print CPU time
459 *
460  IF( ctime( 1 ).GT.0.0d+0 ) THEN
461  tmflops = nops / ctime( 1 )
462  ELSE
463  tmflops = 0.0d+0
464  END IF
465  IF( ctime( 1 ).GE.0.0d+0 )
466  $ WRITE( nout, fmt = 9993 )'CPU ', uplo, n, nb,
467  $ nprow, npcol, ctime( 1 ), tmflops, fresid, passed
468  END IF
469  10 CONTINUE
470  20 CONTINUE
471 *
472  CALL blacs_gridexit( ictxt )
473  30 CONTINUE
474 *
475  CALL pdttrdtester( iam, nprocs, check, nout, thresh, nval, nmat,
476  $ mem, totmem, kpass, kfail, kskip )
477 *
478 * Print ending messages and close output file
479 *
480  IF( iam.EQ.0 ) THEN
481  ktests = kpass + kfail + kskip
482  WRITE( nout, fmt = * )
483  WRITE( nout, fmt = 9992 )ktests
484  IF( check ) THEN
485  WRITE( nout, fmt = 9991 )kpass
486  WRITE( nout, fmt = 9989 )kfail
487  ELSE
488  WRITE( nout, fmt = 9990 )kpass
489  END IF
490  WRITE( nout, fmt = 9988 )kskip
491  WRITE( nout, fmt = * )
492  WRITE( nout, fmt = * )
493  WRITE( nout, fmt = 9987 )
494  IF( nout.NE.6 .AND. nout.NE.0 )
495  $ CLOSE ( nout )
496  END IF
497 *
498  CALL blacs_exit( 0 )
499 *
500  9999 FORMAT( 'ILLEGAL ', a6, ': ', a5, ' = ', i3,
501  $ '; It should be at least 1' )
502  9998 FORMAT( 'ILLEGAL GRID: nprow*npcol = ', i4, '. It can be at most',
503  $ i4 )
504  9997 FORMAT( 'Bad ', a6, ' parameters: going on to next test case.' )
505  9996 FORMAT( 'Unable to perform ', a, ': need TOTMEM of at least',
506  $ i11 )
507  9995 FORMAT( 'TIME UPLO N NB P Q TRD Time ',
508  $ ' MFLOPS Residual CHECK' )
509  9994 FORMAT( '---- ---- ------ --- ----- ----- --------- ',
510  $ '----------- -------- ------' )
511  9993 FORMAT( a4, 1x, a4, 1x, i6, 1x, i3, 1x, i5, 1x, i5, 1x, f9.2, 1x,
512  $ f11.2, 1x, f8.2, 1x, a6 )
513  9992 FORMAT( 'Finished', i4, ' tests, with the following results:' )
514  9991 FORMAT( i5, ' tests completed and passed residual checks.' )
515  9990 FORMAT( i5, ' tests completed without checking.' )
516  9989 FORMAT( i5, ' tests completed and failed residual checks.' )
517  9988 FORMAT( i5, ' tests skipped because of illegal input values.' )
518  9987 FORMAT( 'END OF TESTS.' )
519  9986 FORMAT( '||A - Q*T*Q''|| / (||A|| * N * eps) = ', g25.7 )
520 *
521  stop
522 *
523 * End of PDTRDDRIVER
524 *
525  END
max
#define max(A, B)
Definition: pcgemr.c:180
pdlansy
double precision function pdlansy(NORM, UPLO, N, A, IA, JA, DESCA, WORK)
Definition: pdlansy.f:3
pdchekpad
subroutine pdchekpad(ICTXT, MESS, M, N, A, LDA, IPRE, IPOST, CHKVAL)
Definition: pdchekpad.f:3
ilcm
integer function ilcm(M, N)
Definition: ilcm.f:2
sltimer
subroutine sltimer(I)
Definition: sltimer.f:47
lsame
logical function lsame(CA, CB)
Definition: tools.f:1724
pdsytdrv
subroutine pdsytdrv(UPLO, N, A, IA, JA, DESCA, D, E, TAU, WORK, INFO)
Definition: pdsytdrv.f:3
descinit
subroutine descinit(DESC, M, N, MB, NB, IRSRC, ICSRC, ICTXT, LLD, INFO)
Definition: descinit.f:3
slboot
subroutine slboot()
Definition: sltimer.f:2
pdmatgen
subroutine pdmatgen(ICTXT, AFORM, DIAG, M, N, MB, NB, A, LDA, IAROW, IACOL, ISEED, IROFF, IRNUM, ICOFF, ICNUM, MYROW, MYCOL, NPROW, NPCOL)
Definition: pdmatgen.f:4
pdsytrd
subroutine pdsytrd(UPLO, N, A, IA, JA, DESCA, D, E, TAU, WORK, LWORK, INFO)
Definition: pdsytrd.f:3
pdlafchk
subroutine pdlafchk(AFORM, DIAG, M, N, A, IA, JA, DESCA, IASEED, ANORM, FRESID, WORK)
Definition: pdlafchk.f:3
numroc
integer function numroc(N, NB, IPROC, ISRCPROC, NPROCS)
Definition: numroc.f:2
pdfillpad
subroutine pdfillpad(ICTXT, M, N, A, LDA, IPRE, IPOST, CHKVAL)
Definition: pdfillpad.f:2
pdtrddriver
program pdtrddriver
Definition: pdtrddriver.f:1
slcombine
subroutine slcombine(ICTXT, SCOPE, OP, TIMETYPE, N, IBEG, TIMES)
Definition: sltimer.f:267
pdttrdtester
subroutine pdttrdtester(IAM, NPROCS, CHECK, NOUT, THRESH, NVAL, NMAT, MEM, TOTMEM, KPASS, KFAIL, KSKIP)
Definition: pdttrdtester.f:3
pdtrdinfo
subroutine pdtrdinfo(SUMMRY, NOUT, UPLO, NMAT, NVAL, LDNVAL, NNB, NBVAL, LDNBVAL, NGRIDS, PVAL, LDPVAL, QVAL, LDQVAL, THRESH, WORK, IAM, NPROCS)
Definition: pdtrdinfo.f:4
iceil
integer function iceil(INUM, IDENOM)
Definition: iceil.f:2