ScaLAPACK 2.1  2.1
ScaLAPACK: Scalable Linear Algebra PACKage
pspbdriver.f
Go to the documentation of this file.
1  PROGRAM pspbdriver
2 *
3 *
4 * -- ScaLAPACK routine (version 1.7) --
5 * University of Tennessee, Knoxville, Oak Ridge National Laboratory,
6 * and University of California, Berkeley.
7 * November 15, 1997
8 *
9 * Purpose
10 * =======
11 *
12 * PSPBDRIVER is a test program for the
13 * ScaLAPACK Band Cholesky routines corresponding to the options
14 * indicated by SPB. This test driver performs an
15 * A = L*L**T factorization
16 * and solves a linear system with the factors for 1 or more RHS.
17 *
18 * The program must be driven by a short data file.
19 * Here's an example file:
20 *'ScaLAPACK, Version 1.2, banded linear systems input file'
21 *'PVM.'
22 *'' output file name (if any)
23 *6 device out
24 *'L' define Lower or Upper
25 *9 number of problem sizes
26 *1 5 17 28 37 121 200 1023 2048 3073 values of N
27 *6 number of bandwidths
28 *1 2 4 10 31 64 values of BW
29 *1 number of NB's
30 *-1 3 4 5 values of NB (-1 for automatic choice)
31 *1 number of NRHS's (must be 1)
32 *8 values of NRHS
33 *1 number of NBRHS's (ignored)
34 *1 values of NBRHS (ignored)
35 *6 number of process grids
36 *1 2 3 4 5 7 8 15 26 47 64 values of "Number of Process Columns"
37 *3.0 threshold
38 *
39 * Internal Parameters
40 * ===================
41 *
42 * TOTMEM INTEGER, default = 6200000.
43 * TOTMEM is a machine-specific parameter indicating the
44 * maximum amount of available memory in bytes.
45 * The user should customize TOTMEM to his platform. Remember
46 * to leave room in memory for the operating system, the BLACS
47 * buffer, etc. For example, on a system with 8 MB of memory
48 * per process (e.g., one processor on an Intel iPSC/860), the
49 * parameters we use are TOTMEM=6200000 (leaving 1.8 MB for OS,
50 * code, BLACS buffer, etc). However, for PVM, we usually set
51 * TOTMEM = 2000000. Some experimenting with the maximum value
52 * of TOTMEM may be required.
53 *
54 * INTGSZ INTEGER, default = 4 bytes.
55 * REALSZ INTEGER, default = 4 bytes.
56 * INTGSZ and REALSZ indicate the length in bytes on the
57 * given platform for an integer and a single precision real.
58 * MEM REAL array, dimension ( TOTMEM / REALSZ )
59 * All arrays used by ScaLAPACK routines are allocated from
60 * this array and referenced by pointers. The integer IPB,
61 * for example, is a pointer to the starting element of MEM for
62 * the solution vector(s) B.
63 *
64 * =====================================================================
65 *
66 * Code Developer: Andrew J. Cleary, University of Tennessee.
67 * Current address: Lawrence Livermore National Labs.
68 * This version released: August, 2001.
69 *
70 * =====================================================================
71 *
72 * .. Parameters ..
73  INTEGER totmem
74  parameter( totmem = 3000000 )
75  INTEGER block_cyclic_2d, csrc_, ctxt_, dlen_, dtype_,
76  $ lld_, mb_, m_, nb_, n_, rsrc_
77  parameter( block_cyclic_2d = 1, dlen_ = 9, dtype_ = 1,
78  $ ctxt_ = 2, m_ = 3, n_ = 4, mb_ = 5, nb_ = 6,
79  $ rsrc_ = 7, csrc_ = 8, lld_ = 9 )
80 *
81  REAL zero
82  INTEGER memsiz, ntests, realsz
83  REAL padval
84  parameter( realsz = 4,
85  $ memsiz = totmem / realsz, ntests = 20,
86  $ padval = -9923.0e+0, zero = 0.0e+0 )
87  INTEGER int_one
88  parameter( int_one = 1 )
89 * ..
90 * .. Local Scalars ..
91  LOGICAL check
92  CHARACTER uplo
93  CHARACTER*6 passed
94  CHARACTER*80 outfile
95  INTEGER bw, bw_num, fillin_size, free_ptr, h, hh, i,
96  $ iam, iaseed, ibseed, ictxt, ictxtb, ierr_temp,
97  $ imidpad, info, ipa, ipb, ipostpad, iprepad,
98  $ ipw, ipw_size, ipw_solve, ipw_solve_size,
99  $ ip_driver_w, ip_fillin, j, k, kfail, kpass,
100  $ kskip, ktests, mycol, myrhs_size, myrow, n, nb,
101  $ nbw, ngrids, nmat, nnb, nnbr, nnr, nout, np,
102  $ npcol, nprocs, nprocs_real, nprow, nq, nrhs,
103  $ n_first, n_last, worksiz
104  REAL anorm, sresid, thresh
105  DOUBLE PRECISION nops, nops2, tmflops, tmflops2
106 * ..
107 * .. Local Arrays ..
108  INTEGER bwval( ntests ), desca( 7 ), desca2d( dlen_ ),
109  $ descb( 7 ), descb2d( dlen_ ), ierr( 1 ),
110  $ nbrval( ntests ), nbval( ntests ),
111  $ nrval( ntests ), nval( ntests ),
112  $ pval( ntests ), qval( ntests )
113  REAL mem( memsiz )
114  DOUBLE PRECISION ctime( 2 ), wtime( 2 )
115 * ..
116 * .. External Subroutines ..
117  EXTERNAL blacs_barrier, blacs_exit, blacs_get,
118  $ blacs_gridexit, blacs_gridinfo, blacs_gridinit,
119  $ blacs_pinfo, descinit, igsum2d, psbmatgen,
122  $ slcombine, sltimer
123 * ..
124 * .. External Functions ..
125  INTEGER numroc
126  LOGICAL lsame
127  REAL pslange
128  EXTERNAL lsame, numroc, pslange
129 * ..
130 * .. Intrinsic Functions ..
131  INTRINSIC dble, max, min, mod
132 * ..
133 * .. Data Statements ..
134  DATA kfail, kpass, kskip, ktests / 4*0 /
135 * ..
136 *
137 *
138 *
139 * .. Executable Statements ..
140 *
141 * Get starting information
142 *
143  CALL blacs_pinfo( iam, nprocs )
144  iaseed = 100
145  ibseed = 200
146 *
147  CALL pspbinfo( outfile, nout, uplo, nmat, nval, ntests, nbw,
148  $ bwval, ntests, nnb, nbval, ntests, nnr, nrval,
149  $ ntests, nnbr, nbrval, ntests, ngrids, pval, ntests,
150  $ qval, ntests, thresh, mem, iam, nprocs )
151 *
152  check = ( thresh.GE.0.0e+0 )
153 *
154 * Print headings
155 *
156  IF( iam.EQ.0 ) THEN
157  WRITE( nout, fmt = * )
158  WRITE( nout, fmt = 9995 )
159  WRITE( nout, fmt = 9994 )
160  WRITE( nout, fmt = * )
161  END IF
162 *
163 * Loop over different process grids
164 *
165  DO 60 i = 1, ngrids
166 *
167  nprow = pval( i )
168  npcol = qval( i )
169 *
170 * Make sure grid information is correct
171 *
172  ierr( 1 ) = 0
173  IF( nprow.LT.1 ) THEN
174  IF( iam.EQ.0 )
175  $ WRITE( nout, fmt = 9999 ) 'GRID', 'nprow', nprow
176  ierr( 1 ) = 1
177  ELSE IF( npcol.LT.1 ) THEN
178  IF( iam.EQ.0 )
179  $ WRITE( nout, fmt = 9999 ) 'GRID', 'npcol', npcol
180  ierr( 1 ) = 1
181  ELSE IF( nprow*npcol.GT.nprocs ) THEN
182  IF( iam.EQ.0 )
183  $ WRITE( nout, fmt = 9998 ) nprow*npcol, nprocs
184  ierr( 1 ) = 1
185  END IF
186 *
187  IF( ierr( 1 ).GT.0 ) THEN
188  IF( iam.EQ.0 )
189  $ WRITE( nout, fmt = 9997 ) 'grid'
190  kskip = kskip + 1
191  GO TO 50
192  END IF
193 *
194 * Define process grid
195 *
196  CALL blacs_get( -1, 0, ictxt )
197  CALL blacs_gridinit( ictxt, 'Row-major', nprow, npcol )
198 *
199 *
200 * Define transpose process grid
201 *
202  CALL blacs_get( -1, 0, ictxtb )
203  CALL blacs_gridinit( ictxtb, 'Column-major', npcol, nprow )
204 *
205 * Go to bottom of process grid loop if this case doesn't use my
206 * process
207 *
208  CALL blacs_gridinfo( ictxt, nprow, npcol, myrow, mycol )
209 *
210  IF( myrow.LT.0 .OR. mycol.LT.0 ) THEN
211  GO TO 50
212  ENDIF
213 *
214  DO 40 j = 1, nmat
215 *
216  ierr( 1 ) = 0
217 *
218  n = nval( j )
219 *
220 * Make sure matrix information is correct
221 *
222  IF( n.LT.1 ) THEN
223  IF( iam.EQ.0 )
224  $ WRITE( nout, fmt = 9999 ) 'MATRIX', 'N', n
225  ierr( 1 ) = 1
226  END IF
227 *
228 * Check all processes for an error
229 *
230  CALL igsum2d( ictxt, 'All', ' ', 1, 1, ierr, 1,
231  $ -1, 0 )
232 *
233  IF( ierr( 1 ).GT.0 ) THEN
234  IF( iam.EQ.0 )
235  $ WRITE( nout, fmt = 9997 ) 'size'
236  kskip = kskip + 1
237  GO TO 40
238  END IF
239 *
240 *
241  DO 45 bw_num = 1, nbw
242 *
243  ierr( 1 ) = 0
244 *
245  bw = bwval( bw_num )
246  IF( bw.LT.0 ) THEN
247  IF( iam.EQ.0 )
248  $ WRITE( nout, fmt = 9999 ) 'Band', 'bw', bw
249  ierr( 1 ) = 1
250  END IF
251 *
252  IF( bw.GT.n-1 ) THEN
253  ierr( 1 ) = 1
254  END IF
255 *
256 * Check all processes for an error
257 *
258  CALL igsum2d( ictxt, 'All', ' ', 1, 1, ierr, 1,
259  $ -1, 0 )
260 *
261  IF( ierr( 1 ).GT.0 ) THEN
262  kskip = kskip + 1
263  GO TO 45
264  END IF
265 *
266  DO 30 k = 1, nnb
267 *
268  ierr( 1 ) = 0
269 *
270  nb = nbval( k )
271  IF( nb.LT.0 ) THEN
272  nb =( (n-(npcol-1)*bw-1)/npcol + 1 )
273  $ + bw
274  nb = max( nb, 2*bw )
275  nb = min( n, nb )
276  END IF
277 *
278 * Make sure NB is legal
279 *
280  ierr( 1 ) = 0
281  IF( nb.LT.min( 2*bw, n ) ) THEN
282  ierr( 1 ) = 1
283  ENDIF
284 *
285 * Check all processes for an error
286 *
287  CALL igsum2d( ictxt, 'All', ' ', 1, 1, ierr, 1,
288  $ -1, 0 )
289 *
290  IF( ierr( 1 ).GT.0 ) THEN
291  kskip = kskip + 1
292  GO TO 30
293  END IF
294 *
295 * Padding constants
296 *
297  np = numroc( (bw+1), (bw+1),
298  $ myrow, 0, nprow )
299  nq = numroc( n, nb, mycol, 0, npcol )
300 *
301  IF( check ) THEN
302  iprepad = ((bw+1)+10)
303  imidpad = 10
304  ipostpad = ((bw+1)+10)
305  ELSE
306  iprepad = 0
307  imidpad = 0
308  ipostpad = 0
309  END IF
310 *
311 * Initialize the array descriptor for the matrix A
312 *
313  CALL descinit( desca2d, (bw+1), n,
314  $ (bw+1), nb, 0, 0,
315  $ ictxt,((bw+1)+10), ierr( 1 ) )
316 *
317 * Convert this to 1D descriptor
318 *
319  desca( 1 ) = 501
320  desca( 3 ) = n
321  desca( 4 ) = nb
322  desca( 5 ) = 0
323  desca( 2 ) = ictxt
324  desca( 6 ) = ((bw+1)+10)
325  desca( 7 ) = 0
326 *
327  ierr_temp = ierr( 1 )
328  ierr( 1 ) = 0
329  ierr( 1 ) = min( ierr( 1 ), ierr_temp )
330 *
331 * Check all processes for an error
332 *
333  CALL igsum2d( ictxt, 'All', ' ', 1, 1, ierr, 1, -1, 0 )
334 *
335  IF( ierr( 1 ).LT.0 ) THEN
336  IF( iam.EQ.0 )
337  $ WRITE( nout, fmt = 9997 ) 'descriptor'
338  kskip = kskip + 1
339  GO TO 30
340  END IF
341 *
342 * Assign pointers into MEM for SCALAPACK arrays, A is
343 * allocated starting at position MEM( IPREPAD+1 )
344 *
345  free_ptr = 1
346  ipb = 0
347 *
348 * Save room for prepadding
349  free_ptr = free_ptr + iprepad
350 *
351  ipa = free_ptr
352  free_ptr = free_ptr + desca2d( lld_ )*
353  $ desca2d( nb_ )
354  $ + ipostpad
355 *
356 * Add memory for fillin
357 * Fillin space needs to store:
358 * Fillin spike:
359 * Contribution to previous proc's diagonal block of
360 * reduced system:
361 * Off-diagonal block of reduced system:
362 * Diagonal block of reduced system:
363 *
364  fillin_size =
365  $ (nb+2*bw)*bw
366 *
367 * Claim memory for fillin
368 *
369  free_ptr = free_ptr + iprepad
370  ip_fillin = free_ptr
371  free_ptr = free_ptr + fillin_size
372 *
373 * Workspace needed by computational routines:
374 *
375  ipw_size = 0
376 *
377 * factorization:
378 *
379  ipw_size = bw*bw
380 *
381 * Claim memory for IPW
382 *
383  ipw = free_ptr
384  free_ptr = free_ptr + ipw_size
385 *
386 * Check for adequate memory for problem size
387 *
388  ierr( 1 ) = 0
389  IF( free_ptr.GT.memsiz ) THEN
390  IF( iam.EQ.0 )
391  $ WRITE( nout, fmt = 9996 )
392  $ 'divide and conquer factorization',
393  $ (free_ptr )*realsz
394  ierr( 1 ) = 1
395  END IF
396 *
397 * Check all processes for an error
398 *
399  CALL igsum2d( ictxt, 'All', ' ', 1, 1, ierr,
400  $ 1, -1, 0 )
401 *
402  IF( ierr( 1 ).GT.0 ) THEN
403  IF( iam.EQ.0 )
404  $ WRITE( nout, fmt = 9997 ) 'MEMORY'
405  kskip = kskip + 1
406  GO TO 30
407  END IF
408 *
409 * Worksize needed for LAPRNT
410  worksiz = max( ((bw+1)+10), nb )
411 *
412  IF( check ) THEN
413 *
414 * Calculate the amount of workspace required by
415 * the checking routines.
416 *
417 * PSLANGE
418  worksiz = max( worksiz, desca2d( nb_ ) )
419 *
420 * PSPBLASCHK
421  worksiz = max( worksiz,
422  $ max(5,max(bw*(bw+2),nb))+2*nb )
423  END IF
424 *
425  free_ptr = free_ptr + iprepad
426  ip_driver_w = free_ptr
427  free_ptr = free_ptr + worksiz + ipostpad
428 *
429 *
430 * Check for adequate memory for problem size
431 *
432  ierr( 1 ) = 0
433  IF( free_ptr.GT.memsiz ) THEN
434  IF( iam.EQ.0 )
435  $ WRITE( nout, fmt = 9996 ) 'factorization',
436  $ ( free_ptr )*realsz
437  ierr( 1 ) = 1
438  END IF
439 *
440 * Check all processes for an error
441 *
442  CALL igsum2d( ictxt, 'All', ' ', 1, 1, ierr,
443  $ 1, -1, 0 )
444 *
445  IF( ierr( 1 ).GT.0 ) THEN
446  IF( iam.EQ.0 )
447  $ WRITE( nout, fmt = 9997 ) 'MEMORY'
448  kskip = kskip + 1
449  GO TO 30
450  END IF
451 *
452  CALL psbmatgen( ictxt, uplo, 'B', bw, bw, n, (bw+1), nb,
453  $ mem( ipa ), ((bw+1)+10), 0, 0, iaseed,
454  $ myrow, mycol, nprow, npcol )
455 *
456  CALL psfillpad( ictxt, np, nq, mem( ipa-iprepad ),
457  $ ((bw+1)+10), iprepad, ipostpad,
458  $ padval )
459 *
460  CALL psfillpad( ictxt, worksiz, 1,
461  $ mem( ip_driver_w-iprepad ), worksiz,
462  $ iprepad, ipostpad, padval )
463 *
464 * Calculate norm of A for residual error-checking
465 *
466  IF( check ) THEN
467 *
468  anorm = pslange( '1', (bw+1),
469  $ n, mem( ipa ), 1, 1,
470  $ desca2d, mem( ip_driver_w ) )
471  CALL pschekpad( ictxt, 'PSLANGE', np, nq,
472  $ mem( ipa-iprepad ), ((bw+1)+10),
473  $ iprepad, ipostpad, padval )
474  CALL pschekpad( ictxt, 'PSLANGE',
475  $ worksiz, 1,
476  $ mem( ip_driver_w-iprepad ), worksiz,
477  $ iprepad, ipostpad, padval )
478  END IF
479 *
480 *
481  CALL slboot()
482  CALL blacs_barrier( ictxt, 'All' )
483 *
484 * Perform factorization
485 *
486  CALL sltimer( 1 )
487 *
488  CALL pspbtrf( uplo, n, bw, mem( ipa ), 1, desca,
489  $ mem( ip_fillin ), fillin_size, mem( ipw ),
490  $ ipw_size, info )
491 *
492  CALL sltimer( 1 )
493 *
494  IF( info.NE.0 ) THEN
495  IF( iam.EQ.0 ) THEN
496  WRITE( nout, fmt = * ) 'PSPBTRF INFO=', info
497  ENDIF
498  kfail = kfail + 1
499  GO TO 30
500  END IF
501 *
502  IF( check ) THEN
503 *
504 * Check for memory overwrite in factorization
505 *
506  CALL pschekpad( ictxt, 'PSPBTRF', np,
507  $ nq, mem( ipa-iprepad ), ((bw+1)+10),
508  $ iprepad, ipostpad, padval )
509  END IF
510 *
511 *
512 * Loop over the different values for NRHS
513 *
514  DO 20 hh = 1, nnr
515 *
516  ierr( 1 ) = 0
517 *
518  nrhs = nrval( hh )
519 *
520 * Initialize Array Descriptor for rhs
521 *
522  CALL descinit( descb2d, n, nrhs, nb, 1, 0, 0,
523  $ ictxtb, nb+10, ierr( 1 ) )
524 *
525 * Convert this to 1D descriptor
526 *
527  descb( 1 ) = 502
528  descb( 3 ) = n
529  descb( 4 ) = nb
530  descb( 5 ) = 0
531  descb( 2 ) = ictxt
532  descb( 6 ) = descb2d( lld_ )
533  descb( 7 ) = 0
534 *
535 * reset free_ptr to reuse space for right hand sides
536 *
537  IF( ipb .GT. 0 ) THEN
538  free_ptr = ipb
539  ENDIF
540 *
541  free_ptr = free_ptr + iprepad
542  ipb = free_ptr
543  free_ptr = free_ptr + nrhs*descb2d( lld_ )
544  $ + ipostpad
545 *
546 * Allocate workspace for workspace in TRS routine:
547 *
548  ipw_solve_size = (bw*nrhs)
549 *
550  ipw_solve = free_ptr
551  free_ptr = free_ptr + ipw_solve_size
552 *
553  ierr( 1 ) = 0
554  IF( free_ptr.GT.memsiz ) THEN
555  IF( iam.EQ.0 )
556  $ WRITE( nout, fmt = 9996 )'solve',
557  $ ( free_ptr )*realsz
558  ierr( 1 ) = 1
559  END IF
560 *
561 * Check all processes for an error
562 *
563  CALL igsum2d( ictxt, 'All', ' ', 1, 1,
564  $ ierr, 1, -1, 0 )
565 *
566  IF( ierr( 1 ).GT.0 ) THEN
567  IF( iam.EQ.0 )
568  $ WRITE( nout, fmt = 9997 ) 'MEMORY'
569  kskip = kskip + 1
570  GO TO 15
571  END IF
572 *
573  myrhs_size = numroc( n, nb, mycol, 0, npcol )
574 *
575 * Generate RHS
576 *
577  CALL psmatgen(ictxtb, 'No', 'No',
578  $ descb2d( m_ ), descb2d( n_ ),
579  $ descb2d( mb_ ), descb2d( nb_ ),
580  $ mem( ipb ),
581  $ descb2d( lld_ ), descb2d( rsrc_ ),
582  $ descb2d( csrc_ ),
583  $ ibseed, 0, myrhs_size, 0, nrhs, mycol,
584  $ myrow, npcol, nprow )
585 *
586  IF( check ) THEN
587  CALL psfillpad( ictxtb, nb, nrhs,
588  $ mem( ipb-iprepad ),
589  $ descb2d( lld_ ),
590  $ iprepad, ipostpad,
591  $ padval )
592  CALL psfillpad( ictxt, worksiz, 1,
593  $ mem( ip_driver_w-iprepad ),
594  $ worksiz, iprepad,
595  $ ipostpad, padval )
596  END IF
597 *
598 *
599  CALL blacs_barrier( ictxt, 'All')
600  CALL sltimer( 2 )
601 *
602 * Solve linear system via factorization
603 *
604  CALL pspbtrs( uplo, n, bw, nrhs, mem( ipa ), 1,
605  $ desca, mem( ipb ), 1, descb,
606  $ mem( ip_fillin ), fillin_size,
607  $ mem( ipw_solve ), ipw_solve_size,
608  $ info )
609 *
610  CALL sltimer( 2 )
611 *
612  IF( info.NE.0 ) THEN
613  IF( iam.EQ.0 )
614  $ WRITE( nout, fmt = * ) 'PSPBTRS INFO=', info
615  kfail = kfail + 1
616  passed = 'FAILED'
617  GO TO 20
618  END IF
619 *
620  IF( check ) THEN
621 *
622 * check for memory overwrite
623 *
624  CALL pschekpad( ictxt, 'PSPBTRS-work',
625  $ worksiz, 1,
626  $ mem( ip_driver_w-iprepad ),
627  $ worksiz, iprepad,
628  $ ipostpad, padval )
629 *
630 * check the solution to rhs
631 *
632  sresid = zero
633 *
634  CALL pspblaschk( 'S', uplo, n, bw, bw, nrhs,
635  $ mem( ipb ), 1, 1, descb2d,
636  $ iaseed, mem( ipa ), 1, 1, desca2d,
637  $ ibseed, anorm, sresid,
638  $ mem( ip_driver_w ), worksiz )
639 *
640  IF( iam.EQ.0 ) THEN
641  IF( sresid.GT.thresh )
642  $ WRITE( nout, fmt = 9985 ) sresid
643  END IF
644 *
645 * The second test is a NaN trap
646 *
647  IF( ( sresid.LE.thresh ).AND.
648  $ ( (sresid-sresid).EQ.0.0e+0 ) ) THEN
649  kpass = kpass + 1
650  passed = 'PASSED'
651  ELSE
652  kfail = kfail + 1
653  passed = 'FAILED'
654  END IF
655 *
656  END IF
657 *
658  15 CONTINUE
659 * Skipped tests jump to here to print out "SKIPPED"
660 *
661 * Gather maximum of all CPU and WALL clock timings
662 *
663  CALL slcombine( ictxt, 'All', '>', 'W', 2, 1,
664  $ wtime )
665  CALL slcombine( ictxt, 'All', '>', 'C', 2, 1,
666  $ ctime )
667 *
668 * Print results
669 *
670  IF( myrow.EQ.0 .AND. mycol.EQ.0 ) THEN
671 *
672  nops = 0
673  nops2 = 0
674 *
675  n_first = nb
676  nprocs_real = ( n-1 )/nb + 1
677  n_last = mod( n-1, nb ) + 1
678 *
679 *
680  nops = nops + dble(bw)*( -2.d0 / 3.d0+dble(bw)*
681  $ ( -1.d0+dble(bw)*( -1.d0 / 3.d0 ) ) ) +
682  $ dble(n)*( 1.d0+dble(bw)*( 3.d0 /
683  $ 2.d0+dble(bw)*( 1.d0 / 2.d0 ) ) )
684  nops = nops + dble(bw)*( -1.d0 / 6.d0+dble(bw)
685  $ *( -1.d0 /2.d0+dble(bw)
686  $ *( -1.d0 / 3.d0 ) ) ) +
687  $ dble(n)*( dble(bw) /
688  $ 2.d0*( 1.d0+dble(bw) ) )
689 *
690  nops = nops +
691  $ dble(nrhs)*( ( 2*dble(n)-dble(bw) )*
692  $ ( dble(bw)+1.d0 ) )+ dble(nrhs)*
693  $ ( dble(bw)*( 2*dble(n)-
694  $ ( dble(bw)+1.d0 ) ) )
695 *
696 *
697 * Second calc to represent actual hardware speed
698 *
699 * NB bw^2 flops for LLt factorization in 1st proc
700 *
701  nops2 = ( (dble(n_first))* dble(bw)**2 )
702 *
703  IF ( nprocs_real .GT. 1) THEN
704 * 4 NB bw^2 flops for LLt factorization and
705 * spike calc in last processor
706 *
707  nops2 = nops2 +
708  $ 4*( (dble(n_last)*dble(bw)**2) )
709  ENDIF
710 *
711  IF ( nprocs_real .GT. 2) THEN
712 * 4 NB bw^2 flops for LLt factorization and
713 * spike calc in other processors
714 *
715  nops2 = nops2 + (nprocs_real-2)*
716  $ 4*( (dble(nb)*dble(bw)**2) )
717  ENDIF
718 *
719 * Reduced system
720 *
721  nops2 = nops2 +
722  $ ( nprocs_real-1 ) * ( bw*bw*bw/3 )
723  IF( nprocs_real .GT. 1 ) THEN
724  nops2 = nops2 +
725  $ ( nprocs_real-2 ) * ( 2 * bw*bw*bw )
726  ENDIF
727 *
728 *
729 * nrhs * 4 n_first*bw flops for LLt solve in proc 1.
730 *
731  nops2 = nops2 +
732  $ ( 4.0d+0*(dble(n_first)*dble(bw))*dble(nrhs) )
733 *
734  IF ( nprocs_real .GT. 1 ) THEN
735 *
736 * 2*nrhs*4 n_last*bw flops for LLt solve in last.
737 *
738  nops2 = nops2 +
739  $ 2*( 4.0d+0*(dble(n_last)*dble(bw))*dble(nrhs) )
740  ENDIF
741 *
742  IF ( nprocs_real .GT. 2 ) THEN
743 *
744 * 2 * nrhs * 4 NB*bw flops for LLt solve in others.
745 *
746  nops2 = nops2 +
747  $ ( nprocs_real-2)*2*
748  $ ( 4.0d+0*(dble(nb)*dble(bw))*dble(nrhs) )
749  ENDIF
750 *
751 * Reduced system
752 *
753  nops2 = nops2 +
754  $ nrhs*( nprocs_real-1 ) * ( bw*bw )
755  IF( nprocs_real .GT. 1 ) THEN
756  nops2 = nops2 +
757  $ nrhs*( nprocs_real-2 ) * ( 3 * bw*bw )
758  ENDIF
759 *
760 *
761 * Calculate total megaflops - factorization and/or
762 * solve -- for WALL and CPU time, and print output
763 *
764 * Print WALL time if machine supports it
765 *
766  IF( wtime( 1 ) + wtime( 2 ) .GT. 0.0d+0 ) THEN
767  tmflops = nops /
768  $ ( ( wtime( 1 )+wtime( 2 ) ) * 1.0d+6 )
769  ELSE
770  tmflops = 0.0d+0
771  END IF
772 *
773  IF( wtime( 1 )+wtime( 2 ).GT.0.0d+0 ) THEN
774  tmflops2 = nops2 /
775  $ ( ( wtime( 1 )+wtime( 2 ) ) * 1.0d+6 )
776  ELSE
777  tmflops2 = 0.0d+0
778  END IF
779 *
780  IF( wtime( 2 ).GE.0.0d+0 )
781  $ WRITE( nout, fmt = 9993 ) 'WALL', uplo,
782  $ n,
783  $ bw,
784  $ nb, nrhs, nprow, npcol,
785  $ wtime( 1 ), wtime( 2 ), tmflops,
786  $ tmflops2, passed
787 *
788 * Print CPU time if machine supports it
789 *
790  IF( ctime( 1 )+ctime( 2 ).GT.0.0d+0 ) THEN
791  tmflops = nops /
792  $ ( ( ctime( 1 )+ctime( 2 ) ) * 1.0d+6 )
793  ELSE
794  tmflops = 0.0d+0
795  END IF
796 *
797  IF( ctime( 1 )+ctime( 2 ).GT.0.0d+0 ) THEN
798  tmflops2 = nops2 /
799  $ ( ( ctime( 1 )+ctime( 2 ) ) * 1.0d+6 )
800  ELSE
801  tmflops2 = 0.0d+0
802  END IF
803 *
804  IF( ctime( 2 ).GE.0.0d+0 )
805  $ WRITE( nout, fmt = 9993 ) 'CPU ', uplo,
806  $ n,
807  $ bw,
808  $ nb, nrhs, nprow, npcol,
809  $ ctime( 1 ), ctime( 2 ), tmflops,
810  $ tmflops2, passed
811 *
812  END IF
813  20 CONTINUE
814 *
815 *
816  30 CONTINUE
817 * NNB loop
818 *
819  45 CONTINUE
820 * BW[] loop
821 *
822  40 CONTINUE
823 * NMAT loop
824 *
825  CALL blacs_gridexit( ictxt )
826  CALL blacs_gridexit( ictxtb )
827 *
828  50 CONTINUE
829 * NGRIDS DROPOUT
830  60 CONTINUE
831 * NGRIDS loop
832 *
833 * Print ending messages and close output file
834 *
835  IF( iam.EQ.0 ) THEN
836  ktests = kpass + kfail + kskip
837  WRITE( nout, fmt = * )
838  WRITE( nout, fmt = 9992 ) ktests
839  IF( check ) THEN
840  WRITE( nout, fmt = 9991 ) kpass
841  WRITE( nout, fmt = 9989 ) kfail
842  ELSE
843  WRITE( nout, fmt = 9990 ) kpass
844  END IF
845  WRITE( nout, fmt = 9988 ) kskip
846  WRITE( nout, fmt = * )
847  WRITE( nout, fmt = * )
848  WRITE( nout, fmt = 9987 )
849  IF( nout.NE.6 .AND. nout.NE.0 )
850  $ CLOSE ( nout )
851  END IF
852 *
853  CALL blacs_exit( 0 )
854 *
855  9999 FORMAT( 'ILLEGAL ', a6, ': ', a5, ' = ', i3,
856  $ '; It should be at least 1' )
857  9998 FORMAT( 'ILLEGAL GRID: nprow*npcol = ', i4, '. It can be at most',
858  $ i4 )
859  9997 FORMAT( 'Bad ', a6, ' parameters: going on to next test case.' )
860  9996 FORMAT( 'Unable to perform ', a, ': need TOTMEM of at least',
861  $ i11 )
862  9995 FORMAT( 'TIME UL N BW NB NRHS P Q L*U Time ',
863  $ 'Slv Time MFLOPS MFLOP2 CHECK' )
864  9994 FORMAT( '---- -- ------ --- ---- ----- -- ---- -------- ',
865  $ '-------- ------ ------ ------' )
866  9993 FORMAT( a4, 2x, a1, 1x, i6, 1x, i3, 1x, i4, 1x,
867  $ i5, 1x, i2, 1x,
868  $ i4, 1x, f8.3, f9.4, f9.2, f9.2, 1x, a6 )
869  9992 FORMAT( 'Finished ', i6, ' tests, with the following results:' )
870  9991 FORMAT( i5, ' tests completed and passed residual checks.' )
871  9990 FORMAT( i5, ' tests completed without checking.' )
872  9989 FORMAT( i5, ' tests completed and failed residual checks.' )
873  9988 FORMAT( i5, ' tests skipped because of illegal input values.' )
874  9987 FORMAT( 'END OF TESTS.' )
875  9986 FORMAT( '||A - ', a4, '|| / (||A|| * N * eps) = ', g25.7 )
876  9985 FORMAT( '||Ax-b||/(||x||*||A||*eps*N) ', f25.7 )
877 *
878  stop
879 *
880 * End of PSPBTRS_DRIVER
881 *
882  END
883 *
pspbtrs
subroutine pspbtrs(UPLO, N, BW, NRHS, A, JA, DESCA, B, IB, DESCB, AF, LAF, WORK, LWORK, INFO)
Definition: pspbtrs.f:3
max
#define max(A, B)
Definition: pcgemr.c:180
pspbdriver
program pspbdriver
Definition: pspbdriver.f:1
pslange
real function pslange(NORM, M, N, A, IA, JA, DESCA, WORK)
Definition: pslange.f:3
sltimer
subroutine sltimer(I)
Definition: sltimer.f:47
lsame
logical function lsame(CA, CB)
Definition: tools.f:1724
pschekpad
subroutine pschekpad(ICTXT, MESS, M, N, A, LDA, IPRE, IPOST, CHKVAL)
Definition: pschekpad.f:3
pspbinfo
subroutine pspbinfo(SUMMRY, NOUT, UPLO, NMAT, NVAL, LDNVAL, NBW, BWVAL, LDBWVAL, NNB, NBVAL, LDNBVAL, NNR, NRVAL, LDNRVAL, NNBR, NBRVAL, LDNBRVAL, NGRIDS, PVAL, LDPVAL, QVAL, LDQVAL, THRESH, WORK, IAM, NPROCS)
Definition: pspbinfo.f:6
pspbtrf
subroutine pspbtrf(UPLO, N, BW, A, JA, DESCA, AF, LAF, WORK, LWORK, INFO)
Definition: pspbtrf.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
pspblaschk
subroutine pspblaschk(SYMM, UPLO, N, BWL, BWU, NRHS, X, IX, JX, DESCX, IASEED, A, IA, JA, DESCA, IBSEED, ANORM, RESID, WORK, WORKSIZ)
Definition: pspblaschk.f:4
psmatgen
subroutine psmatgen(ICTXT, AFORM, DIAG, M, N, MB, NB, A, LDA, IAROW, IACOL, ISEED, IROFF, IRNUM, ICOFF, ICNUM, MYROW, MYCOL, NPROW, NPCOL)
Definition: psmatgen.f:4
numroc
integer function numroc(N, NB, IPROC, ISRCPROC, NPROCS)
Definition: numroc.f:2
psbmatgen
subroutine psbmatgen(ICTXT, AFORM, AFORM2, BWL, BWU, N, MB, NB, A, LDA, IAROW, IACOL, ISEED, MYROW, MYCOL, NPROW, NPCOL)
Definition: psbmatgen.f:5
psfillpad
subroutine psfillpad(ICTXT, M, N, A, LDA, IPRE, IPOST, CHKVAL)
Definition: psfillpad.f:2
slcombine
subroutine slcombine(ICTXT, SCOPE, OP, TIMETYPE, N, IBEG, TIMES)
Definition: sltimer.f:267
min
#define min(A, B)
Definition: pcgemr.c:181