68 INTEGER block_cyclic_2d, csrc_, ctxt_, dlen_, dtype_,
69 $ lld_, mb_, m_, nb_, n_, rsrc_
70 parameter( block_cyclic_2d = 1, dlen_ = 9, dtype_ = 1,
71 $ ctxt_ = 2, m_ = 3, n_ = 4, mb_ = 5, nb_ = 6,
72 $ rsrc_ = 7, csrc_ = 8, lld_ = 9 )
73 INTEGER cplxsz, intgsz, memsiz, ntests, realsz, totmem
75 parameter( cplxsz = 8, intgsz = 4, realsz = 4,
76 $ totmem = 2000000, memsiz = totmem / cplxsz,
78 $ padval = ( -9923.0e+0, -9923.0e+0 ) )
87 INTEGER i, iam, iaseed, ictxt, imidpad, info, ipa,
88 $ ipostpad, ippiv, iprepad, iptau, iprw, ipw, j,
89 $ k, kfail, kpass, kskip, ktests, l, lipiv,
90 $ lrwork, ltau, lwork, m, maxmn, mb, minmn, mnp,
91 $ mnq, mp, mycol, myrow, n, nb, nfact, ngrids,
92 $ nmat, nnb, nout, npcol, nprocs, nprow, nq,
93 $ workfct, workrfct, worksiz
94 REAL anorm, fresid, thresh
95 DOUBLE PRECISION nops, tmflops
98 CHARACTER*2 factor( ntests )
99 INTEGER desca( dlen_ ), ierr( 1 ), mbval( ntests ),
100 $ mval( ntests ), nbval( ntests ),
101 $ nval( ntests ), pval( ntests ), qval( ntests )
102 DOUBLE PRECISION ctime( 1 ), wtime( 1 )
103 COMPLEX mem( memsiz )
106 EXTERNAL blacs_barrier, blacs_exit, blacs_get,
107 $ blacs_gridexit, blacs_gridinfo, blacs_gridinit,
126 DATA ktests, kpass, kfail, kskip /4*0/
132 CALL blacs_pinfo( iam, nprocs )
134 CALL pcqrinfo( outfile, nout, nfact, factor, ntests, nmat, mval,
135 $ ntests, nval, ntests, nnb, mbval, ntests, nbval,
136 $ ntests, ngrids, pval, ntests, qval, ntests,
137 $ thresh, mem, iam, nprocs )
138 check = ( thresh.GE.0.0e+0 )
149 WRITE( nout, fmt = * )
150 IF(
lsamen( 2, fact,
'QR' ) )
THEN
153 WRITE( nout, fmt = 9986 )
154 $
'QR factorization tests.'
155 ELSE IF(
lsamen( 2, fact,
'QL' ) )
THEN
158 WRITE( nout, fmt = 9986 )
159 $
'QL factorization tests.'
160 ELSE IF(
lsamen( 2, fact,
'LQ' ) )
THEN
163 WRITE( nout, fmt = 9986 )
164 $
'LQ factorization tests.'
165 ELSE IF(
lsamen( 2, fact,
'RQ' ) )
THEN
168 WRITE( nout, fmt = 9986 )
169 $
'RQ factorization tests.'
170 ELSE IF(
lsamen( 2, fact,
'QP' ) )
THEN
173 WRITE( nout, fmt = 9986 )
174 $
'QR factorization with column pivoting tests.'
175 ELSE IF(
lsamen( 2, fact,
'TZ' ) )
THEN
178 WRITE( nout, fmt = 9986 )
179 $
'Complete unitary factorization tests.'
181 WRITE( nout, fmt = * )
182 WRITE( nout, fmt = 9995 )
183 WRITE( nout, fmt = 9994 )
184 WRITE( nout, fmt = * )
197 IF( nprow.LT.1 )
THEN
199 $
WRITE( nout, fmt = 9999 )
'GRID',
'nprow', nprow
201 ELSE IF( npcol.LT.1 )
THEN
203 $
WRITE( nout, fmt = 9999 )
'GRID',
'npcol', npcol
205 ELSE IF( nprow*npcol.GT.nprocs )
THEN
207 $
WRITE( nout, fmt = 9998 ) nprow*npcol, nprocs
211 IF( ierr( 1 ).GT.0 )
THEN
213 $
WRITE( nout, fmt = 9997 )
'grid'
220 CALL blacs_get( -1, 0, ictxt )
221 CALL blacs_gridinit( ictxt,
'Row-major', nprow, npcol )
222 CALL blacs_gridinfo( ictxt, nprow, npcol, myrow, mycol )
226 IF( myrow.GE.nprow .OR. mycol.GE.npcol )
239 $
WRITE( nout, fmt = 9999 )
'MATRIX',
'M', m
241 ELSE IF( n.LT.1 )
THEN
243 $
WRITE( nout, fmt = 9999 )
'MATRIX',
'N', n
249 CALL igsum2d( ictxt,
'All',
' ', 1, 1, ierr, 1, -1, 0 )
251 IF( ierr( 1 ).GT.0 )
THEN
253 $
WRITE( nout, fmt = 9997 )
'matrix'
271 $
WRITE( nout, fmt = 9999 )
'MB',
'MB', mb
276 CALL igsum2d( ictxt,
'All',
' ', 1, 1, ierr, 1, -1,
279 IF( ierr( 1 ).GT.0 )
THEN
281 $
WRITE( nout, fmt = 9997 )
'MB'
292 $
WRITE( nout, fmt = 9999 )
'NB',
'NB', nb
297 CALL igsum2d( ictxt,
'All',
' ', 1, 1, ierr, 1, -1,
300 IF( ierr( 1 ).GT.0 )
THEN
302 $
WRITE( nout, fmt = 9997 )
'NB'
309 mp =
numroc( m, mb, myrow, 0, nprow )
310 nq =
numroc( n, nb, mycol, 0, npcol )
311 mnp =
numroc(
min( m, n ), mb, myrow, 0, nprow )
312 mnq =
numroc(
min( m, n ), nb, mycol, 0, npcol )
314 iprepad =
max( mb, mp )
316 ipostpad =
max( nb, nq )
325 CALL descinit( desca, m, n, mb, nb, 0, 0, ictxt,
326 $
max( 1, mp ) + imidpad, ierr( 1 ) )
330 CALL igsum2d( ictxt,
'All',
' ', 1, 1, ierr, 1, -1,
333 IF( ierr( 1 ).LT.0 )
THEN
335 $
WRITE( nout, fmt = 9997 )
'descriptor'
344 iptau = ipa + desca( lld_ ) * nq + ipostpad + iprepad
346 IF(
lsamen( 2, fact,
'QR' ) )
THEN
349 ipw = iptau + ltau + ipostpad + iprepad
354 lwork = desca( nb_ ) * ( mp + nq + desca( nb_ ) )
355 workfct = lwork + ipostpad
364 worksiz = lwork + mp*desca( nb_ ) + ipostpad
368 ELSE IF(
lsamen( 2, fact,
'QL' ) )
THEN
371 ipw = iptau + ltau + ipostpad + iprepad
376 lwork = desca( nb_ ) * ( mp + nq + desca( nb_ ) )
377 workfct = lwork + ipostpad
386 worksiz = lwork + mp*desca( nb_ ) + ipostpad
390 ELSE IF(
lsamen( 2, fact,
'LQ' ) )
THEN
393 ipw = iptau + ltau + ipostpad + iprepad
398 lwork = desca( mb_ ) * ( mp + nq + desca( mb_ ) )
399 workfct = lwork + ipostpad
409 $
max( mp*desca( nb_ ), nq*desca( mb_ )
414 ELSE IF(
lsamen( 2, fact,
'RQ' ) )
THEN
417 ipw = iptau + ltau + ipostpad + iprepad
422 lwork = desca( mb_ ) * ( mp + nq + desca( mb_ ) )
423 workfct = lwork + ipostpad
433 $
max( mp*desca( nb_ ), nq*desca( mb_ )
438 ELSE IF(
lsamen( 2, fact,
'QP' ) )
THEN
441 ippiv = iptau + ltau + ipostpad + iprepad
442 lipiv =
iceil( intgsz*nq, cplxsz )
443 ipw = ippiv + lipiv + ipostpad + iprepad
448 lwork =
max( 3, mp +
max( 1, nq ) )
449 workfct = lwork + ipostpad
450 lrwork =
max( 1, 2 * nq )
451 workrfct =
iceil( lrwork*realsz, cplxsz ) +
453 iprw = ipw + workfct + iprepad
454 worksiz = workfct + iprepad + workrfct
462 worksiz =
max( worksiz - ipostpad,
463 $ desca( nb_ )*( 2*mp + nq + desca( nb_ ) ) ) +
467 ELSE IF(
lsamen( 2, fact,
'TZ' ) )
THEN
470 ipw = iptau + ltau + ipostpad + iprepad
475 lwork = desca( mb_ ) * ( mp + nq + desca( mb_ ) )
476 workfct = lwork + ipostpad
486 $
max( mp*desca( nb_ ), nq*desca( mb_ )
496 IF( ipw+worksiz.GT.memsiz )
THEN
498 $
WRITE( nout, fmt = 9996 )
499 $ fact //
' factorization',
500 $ ( ipw+worksiz )*cplxsz
506 CALL igsum2d( ictxt,
'All',
' ', 1, 1, ierr, 1, -1,
509 IF( ierr( 1 ).GT.0 )
THEN
511 $
WRITE( nout, fmt = 9997 )
'MEMORY'
518 CALL pcmatgen( ictxt,
'N',
'N', desca( m_ ),
519 $ desca( n_ ), desca( mb_ ),
520 $ desca( nb_ ), mem( ipa ),
521 $ desca( lld_ ), desca( rsrc_ ),
522 $ desca( csrc_ ), iaseed, 0, mp, 0, nq,
523 $ myrow, mycol, nprow, npcol )
528 CALL pcfillpad( ictxt, mp, nq, mem( ipa-iprepad ),
529 $ desca( lld_ ), iprepad, ipostpad,
531 IF(
lsamen( 2, fact,
'QP' ) )
THEN
533 $ mem( ippiv-iprepad ), lipiv,
534 $ iprepad, ipostpad, padval )
537 $ mem( iptau-iprepad ), ltau,
538 $ iprepad, ipostpad, padval )
539 CALL pcfillpad( ictxt, worksiz-ipostpad, 1,
540 $ mem( ipw-iprepad ),
542 $ iprepad, ipostpad, padval )
543 anorm =
pclange(
'I', m, n, mem( ipa ), 1, 1,
544 $ desca, mem( ipw ) )
545 CALL pcchekpad( ictxt,
'PCLANGE', mp, nq,
546 $ mem( ipa-iprepad ), desca( lld_ ),
547 $ iprepad, ipostpad, padval )
549 $ worksiz-ipostpad, 1,
550 $ mem( ipw-iprepad ),
551 $ worksiz-ipostpad, iprepad,
553 IF(
lsamen( 2, fact,
'QP' ) )
THEN
554 CALL pcfillpad( ictxt, workrfct-ipostpad, 1,
555 $ mem( iprw-iprepad ),
557 $ iprepad, ipostpad, padval )
559 CALL pcfillpad( ictxt, workfct-ipostpad, 1,
560 $ mem( ipw-iprepad ),
562 $ iprepad, ipostpad, padval )
566 CALL blacs_barrier( ictxt,
'All' )
570 IF(
lsamen( 2, fact,
'QR' ) )
THEN
572 CALL pcgeqrf( m, n, mem( ipa ), 1, 1, desca,
573 $ mem( iptau ), mem( ipw ), lwork,
576 ELSE IF(
lsamen( 2, fact,
'QL' ) )
THEN
578 CALL pcgeqlf( m, n, mem( ipa ), 1, 1, desca,
579 $ mem( iptau ), mem( ipw ), lwork,
582 ELSE IF(
lsamen( 2, fact,
'LQ' ) )
THEN
584 CALL pcgelqf( m, n, mem( ipa ), 1, 1, desca,
585 $ mem( iptau ), mem( ipw ), lwork,
588 ELSE IF(
lsamen( 2, fact,
'RQ' ) )
THEN
590 CALL pcgerqf( m, n, mem( ipa ), 1, 1, desca,
591 $ mem( iptau ), mem( ipw ), lwork,
594 ELSE IF(
lsamen( 2, fact,
'QP' ) )
THEN
596 CALL pcgeqpf( m, n, mem( ipa ), 1, 1, desca,
597 $ mem( ippiv ), mem( iptau ),
598 $ mem( ipw ), lwork, mem( iprw ),
601 ELSE IF(
lsamen( 2, fact,
'TZ' ) )
THEN
604 $
CALL pctzrzf( m, n, mem( ipa ), 1, 1, desca,
605 $ mem( iptau ), mem( ipw ), lwork,
615 $ mem( ipa-iprepad ), desca( lld_ ),
616 $ iprepad, ipostpad, padval )
618 $ mem( iptau-iprepad ), ltau,
619 $ iprepad, ipostpad, padval )
620 IF(
lsamen( 2, fact,
'QP' ) )
THEN
622 $ mem( ippiv-iprepad ), lipiv,
623 $ iprepad, ipostpad, padval )
624 CALL pcchekpad( ictxt, rout, workrfct-ipostpad,
625 $ 1, mem( iprw-iprepad ),
627 $ iprepad, ipostpad, padval )
629 CALL pcchekpad( ictxt, rout, workfct-ipostpad, 1,
630 $ mem( ipw-iprepad ),
631 $ workfct-ipostpad, iprepad,
633 CALL pcfillpad( ictxt, worksiz-ipostpad, 1,
634 $ mem( ipw-iprepad ),
636 $ iprepad, ipostpad, padval )
638 IF(
lsamen( 2, fact,
'QR' ) )
THEN
642 CALL pcgeqrrv( m, n, mem( ipa ), 1, 1, desca,
643 $ mem( iptau ), mem( ipw ) )
644 CALL pclafchk(
'No',
'No', m, n, mem( ipa ), 1,
645 $ 1, desca, iaseed, anorm, fresid,
647 ELSE IF(
lsamen( 2, fact,
'QL' ) )
THEN
651 CALL pcgeqlrv( m, n, mem( ipa ), 1, 1, desca,
652 $ mem( iptau ), mem( ipw ) )
653 CALL pclafchk(
'No',
'No', m, n, mem( ipa ), 1,
654 $ 1, desca, iaseed, anorm, fresid,
656 ELSE IF(
lsamen( 2, fact,
'LQ' ) )
THEN
660 CALL pcgelqrv( m, n, mem( ipa ), 1, 1, desca,
661 $ mem( iptau ), mem( ipw ) )
662 CALL pclafchk(
'No',
'No', m, n, mem( ipa ), 1,
663 $ 1, desca, iaseed, anorm, fresid,
665 ELSE IF(
lsamen( 2, fact,
'RQ' ) )
THEN
669 CALL pcgerqrv( m, n, mem( ipa ), 1, 1, desca,
670 $ mem( iptau ), mem( ipw ) )
671 CALL pclafchk(
'No',
'No', m, n, mem( ipa ), 1,
672 $ 1, desca, iaseed, anorm, fresid,
674 ELSE IF(
lsamen( 2, fact,
'QP' ) )
THEN
678 CALL pcgeqrrv( m, n, mem( ipa ), 1, 1, desca,
679 $ mem( iptau ), mem( ipw ) )
680 ELSE IF(
lsamen( 2, fact,
'TZ' ) )
THEN
685 CALL pctzrzrv( m, n, mem( ipa ), 1, 1, desca,
686 $ mem( iptau ), mem( ipw ) )
688 CALL pclafchk(
'No',
'No', m, n, mem( ipa ), 1,
689 $ 1, desca, iaseed, anorm, fresid,
696 $ mem( ipa-iprepad ), desca( lld_ ),
697 $ iprepad, ipostpad, padval )
699 $ mem( iptau-iprepad ), ltau,
700 $ iprepad, ipostpad, padval )
701 CALL pcchekpad( ictxt, routchk, worksiz-ipostpad,
702 $ 1, mem( ipw-iprepad ),
703 $ worksiz-ipostpad, iprepad,
706 IF(
lsamen( 2, fact,
'QP' ) )
THEN
708 CALL pcqppiv( m, n, mem( ipa ), 1, 1, desca,
713 CALL pcchekpad( ictxt,
'PCQPPIV', mp, nq,
714 $ mem( ipa-iprepad ),
716 $ iprepad, ipostpad, padval )
717 CALL pcchekpad( ictxt,
'PCQPPIV', lipiv, 1,
718 $ mem( ippiv-iprepad ), lipiv,
719 $ iprepad, ipostpad, padval )
721 CALL pclafchk(
'No',
'No', m, n, mem( ipa ), 1,
722 $ 1, desca, iaseed, anorm, fresid,
727 CALL pcchekpad( ictxt,
'PCLAFCHK', mp, nq,
728 $ mem( ipa-iprepad ),
730 $ iprepad, ipostpad, padval )
732 $ worksiz-ipostpad, 1,
733 $ mem( ipw-iprepad ),
734 $ worksiz-ipostpad, iprepad,
740 IF(
lsamen( 2, fact,
'TZ' ) .AND. n.LT.m )
THEN
744 IF( fresid.LE.thresh .AND.
745 $ (fresid-fresid).EQ.0.0e+0 )
THEN
759 fresid = fresid - fresid
766 CALL slcombine( ictxt,
'All',
'>',
'W', 1, 1, wtime )
767 CALL slcombine( ictxt,
'All',
'>',
'C', 1, 1, ctime )
771 IF( myrow.EQ.0 .AND. mycol.EQ.0 )
THEN
776 IF(
lsamen( 2, fact,
'TZ' ) )
THEN
785 $ dble( n )*( dble( m )**2 ) -
787 $ 13.0d+0*dble( n )*dble( m ) -
796 nops = 8.0d+0 * ( dble( minmn )**2 ) *
797 $ ( dble( maxmn )-dble( minmn ) / 3.0d+0 ) +
798 $ ( 6.0d+0 * dble( maxmn ) +
799 $ 8.0d+0 * dble( minmn ) ) *
805 IF( wtime( 1 ).GT.0.0d+0 )
THEN
806 tmflops = nops / ( wtime( 1 ) * 1.0d+6 )
810 IF( wtime( 1 ).GE.0.0d+0 )
811 $
WRITE( nout, fmt = 9993 )
'WALL', m, n, mb, nb,
812 $ nprow, npcol, wtime( 1 ), tmflops,
817 IF( ctime( 1 ).GT.0.0d+0 )
THEN
818 tmflops = nops / ( ctime( 1 ) * 1.0d+6 )
822 IF( ctime( 1 ).GE.0.0d+0 )
823 $
WRITE( nout, fmt = 9993 )
'CPU ', m, n, mb, nb,
824 $ nprow, npcol, ctime( 1 ), tmflops,
833 CALL blacs_gridexit( ictxt )
842 ktests = kpass + kfail + kskip
843 WRITE( nout, fmt = * )
844 WRITE( nout, fmt = 9992 ) ktests
846 WRITE( nout, fmt = 9991 ) kpass
847 WRITE( nout, fmt = 9989 ) kfail
849 WRITE( nout, fmt = 9990 ) kpass
851 WRITE( nout, fmt = 9988 ) kskip
852 WRITE( nout, fmt = * )
853 WRITE( nout, fmt = * )
854 WRITE( nout, fmt = 9987 )
855 IF( nout.NE.6 .AND. nout.NE.0 )
861 9999
FORMAT(
'ILLEGAL ', a6,
': ', a5,
' = ', i3,
862 $
'; It should be at least 1' )
863 9998
FORMAT(
'ILLEGAL GRID: nprow*npcol = ', i4,
'. It can be at most',
865 9997
FORMAT(
'Bad ', a6,
' parameters: going on to next test case.' )
866 9996
FORMAT(
'Unable to perform ', a,
': need TOTMEM of at least',
868 9995
FORMAT(
'TIME M N MB NB P Q Fact Time ',
869 $
' MFLOPS CHECK Residual' )
870 9994
FORMAT(
'---- ------ ------ --- --- ----- ----- --------- ',
871 $
'----------- ------ --------' )
872 9993
FORMAT( a4, 1x, i6, 1x, i6, 1x, i3, 1x, i3, 1x, i5, 1x, i5, 1x,
873 $ f9.2, 1x, f11.2, 1x, a6, 2x, g8.1 )
874 9992
FORMAT(
'Finished ', i6,
' tests, with the following results:' )
875 9991
FORMAT( i5,
' tests completed and passed residual checks.' )
876 9990
FORMAT( i5,
' tests completed without checking.' )
877 9989
FORMAT( i5,
' tests completed and failed residual checks.' )
878 9988
FORMAT( i5,
' tests skipped because of illegal input values.' )
879 9987
FORMAT(
'END OF TESTS.' )
888 SUBROUTINE pcqppiv( M, N, A, IA, JA, DESCA, IPIV )
899 INTEGER DESCA( * ), IPIV( * )
999 INTEGER BLOCK_CYCLIC_2D, CSRC_, CTXT_, DLEN_, DTYPE_,
1000 $ LLD_, MB_, M_, NB_, N_, RSRC_
1001 parameter( block_cyclic_2d = 1, dlen_ = 9, dtype_ = 1,
1002 $ ctxt_ = 2, m_ = 3, n_ = 4, mb_ = 5, nb_ = 6,
1003 $ rsrc_ = 7, csrc_ = 8, lld_ = 9 )
1006 INTEGER IACOL, ICOFFA, ICTXT, IITMP, IPVT, IPCOL,
1007 $ IPROW, ITMP, J, JJ, JJA, KK, MYCOL, MYROW,
1011 EXTERNAL blacs_gridinfo, igebr2d, igebs2d, igerv2d,
1012 $ igesd2d, igamn2d,
infog1l, pcswap
1015 INTEGER INDXL2G, NUMROC
1016 EXTERNAL indxl2g, numroc
1025 ictxt = desca( ctxt_ )
1026 CALL blacs_gridinfo( ictxt, nprow, npcol, myrow, mycol )
1027 CALL infog1l( ja, desca( nb_ ), npcol, mycol, desca( csrc_ ), jja,
1029 icoffa = mod( ja-1, desca( nb_ ) )
1030 nq = numroc( n+icoffa, desca( nb_ ), mycol, iacol, npcol )
1031 IF( mycol.EQ.iacol )
1034 DO 20 j = ja, ja+n-2
1041 CALL infog1l( j, desca( nb_ ), npcol, mycol, desca( csrc_ ),
1043 DO 10 kk = jj, jja+nq-1
1044 IF( ipiv( kk ).LT.ipvt )
THEN
1052 CALL igamn2d( ictxt,
'Rowwise',
' ', 1, 1, ipvt, 1, iprow,
1053 $ ipcol, 1, -1, mycol )
1057 IF( mycol.EQ.ipcol )
THEN
1058 itmp = indxl2g( iitmp, desca( nb_ ), mycol, desca( csrc_ ),
1060 CALL igebs2d( ictxt,
'Rowwise',
' ', 1, 1, itmp, 1 )
1061 IF( ipcol.NE.iacol )
THEN
1062 CALL igerv2d( ictxt, 1, 1, ipiv( iitmp ), 1, myrow,
1065 IF( mycol.EQ.iacol )
1066 $ ipiv( iitmp ) = ipiv( jj )
1069 CALL igebr2d( ictxt,
'Rowwise',
' ', 1, 1, itmp, 1, myrow,
1071 IF( mycol.EQ.iacol .AND. ipcol.NE.iacol )
1072 $
CALL igesd2d( ictxt, 1, 1, ipiv( jj ), 1, myrow, ipcol )
1077 CALL pcswap( m, a, ia, itmp, desca, 1, a, ia, j, desca, 1 )