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

◆ pdmmch()

subroutine pdmmch ( integer  ictxt,
character*1  transa,
character*1  transb,
integer  m,
integer  n,
integer  k,
double precision  alpha,
double precision, dimension( * )  a,
integer  ia,
integer  ja,
integer, dimension( * )  desca,
double precision, dimension( * )  b,
integer  ib,
integer  jb,
integer, dimension( * )  descb,
double precision  beta,
double precision, dimension( * )  c,
double precision, dimension( * )  pc,
integer  ic,
integer  jc,
integer, dimension( * )  descc,
double precision, dimension( * )  ct,
double precision, dimension( * )  g,
double precision  err,
integer  info 
)

Definition at line 5269 of file pdblastst.f.

5272*
5273* -- PBLAS test routine (version 2.0) --
5274* University of Tennessee, Knoxville, Oak Ridge National Laboratory,
5275* and University of California, Berkeley.
5276* April 1, 1998
5277*
5278* .. Scalar Arguments ..
5279 CHARACTER*1 TRANSA, TRANSB
5280 INTEGER IA, IB, IC, ICTXT, INFO, JA, JB, JC, K, M, N
5281 DOUBLE PRECISION ALPHA, BETA, ERR
5282* ..
5283* .. Array Arguments ..
5284 INTEGER DESCA( * ), DESCB( * ), DESCC( * )
5285 DOUBLE PRECISION A( * ), B( * ), C( * ), CT( * ), G( * ),
5286 $ PC( * )
5287* ..
5288*
5289* Purpose
5290* =======
5291*
5292* PDMMCH checks the results of the computational tests.
5293*
5294* Notes
5295* =====
5296*
5297* A description vector is associated with each 2D block-cyclicly dis-
5298* tributed matrix. This vector stores the information required to
5299* establish the mapping between a matrix entry and its corresponding
5300* process and memory location.
5301*
5302* In the following comments, the character _ should be read as
5303* "of the distributed matrix". Let A be a generic term for any 2D
5304* block cyclicly distributed matrix. Its description vector is DESCA:
5305*
5306* NOTATION STORED IN EXPLANATION
5307* ---------------- --------------- ------------------------------------
5308* DTYPE_A (global) DESCA( DTYPE_ ) The descriptor type.
5309* CTXT_A (global) DESCA( CTXT_ ) The BLACS context handle, indicating
5310* the NPROW x NPCOL BLACS process grid
5311* A is distributed over. The context
5312* itself is global, but the handle
5313* (the integer value) may vary.
5314* M_A (global) DESCA( M_ ) The number of rows in the distribu-
5315* ted matrix A, M_A >= 0.
5316* N_A (global) DESCA( N_ ) The number of columns in the distri-
5317* buted matrix A, N_A >= 0.
5318* IMB_A (global) DESCA( IMB_ ) The number of rows of the upper left
5319* block of the matrix A, IMB_A > 0.
5320* INB_A (global) DESCA( INB_ ) The number of columns of the upper
5321* left block of the matrix A,
5322* INB_A > 0.
5323* MB_A (global) DESCA( MB_ ) The blocking factor used to distri-
5324* bute the last M_A-IMB_A rows of A,
5325* MB_A > 0.
5326* NB_A (global) DESCA( NB_ ) The blocking factor used to distri-
5327* bute the last N_A-INB_A columns of
5328* A, NB_A > 0.
5329* RSRC_A (global) DESCA( RSRC_ ) The process row over which the first
5330* row of the matrix A is distributed,
5331* NPROW > RSRC_A >= 0.
5332* CSRC_A (global) DESCA( CSRC_ ) The process column over which the
5333* first column of A is distributed.
5334* NPCOL > CSRC_A >= 0.
5335* LLD_A (local) DESCA( LLD_ ) The leading dimension of the local
5336* array storing the local blocks of
5337* the distributed matrix A,
5338* IF( Lc( 1, N_A ) > 0 )
5339* LLD_A >= MAX( 1, Lr( 1, M_A ) )
5340* ELSE
5341* LLD_A >= 1.
5342*
5343* Let K be the number of rows of a matrix A starting at the global in-
5344* dex IA,i.e, A( IA:IA+K-1, : ). Lr( IA, K ) denotes the number of rows
5345* that the process of row coordinate MYROW ( 0 <= MYROW < NPROW ) would
5346* receive if these K rows were distributed over NPROW processes. If K
5347* is the number of columns of a matrix A starting at the global index
5348* JA, i.e, A( :, JA:JA+K-1, : ), Lc( JA, K ) denotes the number of co-
5349* lumns that the process MYCOL ( 0 <= MYCOL < NPCOL ) would receive if
5350* these K columns were distributed over NPCOL processes.
5351*
5352* The values of Lr() and Lc() may be determined via a call to the func-
5353* tion PB_NUMROC:
5354* Lr( IA, K ) = PB_NUMROC( K, IA, IMB_A, MB_A, MYROW, RSRC_A, NPROW )
5355* Lc( JA, K ) = PB_NUMROC( K, JA, INB_A, NB_A, MYCOL, CSRC_A, NPCOL )
5356*
5357* Arguments
5358* =========
5359*
5360* ICTXT (local input) INTEGER
5361* On entry, ICTXT specifies the BLACS context handle, indica-
5362* ting the global context of the operation. The context itself
5363* is global, but the value of ICTXT is local.
5364*
5365* TRANSA (global input) CHARACTER*1
5366* On entry, TRANSA specifies if the matrix operand A is to be
5367* transposed.
5368*
5369* TRANSB (global input) CHARACTER*1
5370* On entry, TRANSB specifies if the matrix operand B is to be
5371* transposed.
5372*
5373* M (global input) INTEGER
5374* On entry, M specifies the number of rows of C.
5375*
5376* N (global input) INTEGER
5377* On entry, N specifies the number of columns of C.
5378*
5379* K (global input) INTEGER
5380* On entry, K specifies the number of columns (resp. rows) of A
5381* when TRANSA = 'N' (resp. TRANSA <> 'N') in PxGEMM, PxSYRK,
5382* PxSYR2K, PxHERK and PxHER2K.
5383*
5384* ALPHA (global input) DOUBLE PRECISION
5385* On entry, ALPHA specifies the scalar alpha.
5386*
5387* A (local input) DOUBLE PRECISION array
5388* On entry, A is an array of dimension (DESCA( M_ ),*). This
5389* array contains a local copy of the initial entire matrix PA.
5390*
5391* IA (global input) INTEGER
5392* On entry, IA specifies A's global row index, which points to
5393* the beginning of the submatrix sub( A ).
5394*
5395* JA (global input) INTEGER
5396* On entry, JA specifies A's global column index, which points
5397* to the beginning of the submatrix sub( A ).
5398*
5399* DESCA (global and local input) INTEGER array
5400* On entry, DESCA is an integer array of dimension DLEN_. This
5401* is the array descriptor for the matrix A.
5402*
5403* B (local input) DOUBLE PRECISION array
5404* On entry, B is an array of dimension (DESCB( M_ ),*). This
5405* array contains a local copy of the initial entire matrix PB.
5406*
5407* IB (global input) INTEGER
5408* On entry, IB specifies B's global row index, which points to
5409* the beginning of the submatrix sub( B ).
5410*
5411* JB (global input) INTEGER
5412* On entry, JB specifies B's global column index, which points
5413* to the beginning of the submatrix sub( B ).
5414*
5415* DESCB (global and local input) INTEGER array
5416* On entry, DESCB is an integer array of dimension DLEN_. This
5417* is the array descriptor for the matrix B.
5418*
5419* BETA (global input) DOUBLE PRECISION
5420* On entry, BETA specifies the scalar beta.
5421*
5422* C (local input/local output) DOUBLE PRECISION array
5423* On entry, C is an array of dimension (DESCC( M_ ),*). This
5424* array contains a local copy of the initial entire matrix PC.
5425*
5426* PC (local input) DOUBLE PRECISION array
5427* On entry, PC is an array of dimension (DESCC( LLD_ ),*). This
5428* array contains the local pieces of the matrix PC.
5429*
5430* IC (global input) INTEGER
5431* On entry, IC specifies C's global row index, which points to
5432* the beginning of the submatrix sub( C ).
5433*
5434* JC (global input) INTEGER
5435* On entry, JC specifies C's global column index, which points
5436* to the beginning of the submatrix sub( C ).
5437*
5438* DESCC (global and local input) INTEGER array
5439* On entry, DESCC is an integer array of dimension DLEN_. This
5440* is the array descriptor for the matrix C.
5441*
5442* CT (workspace) DOUBLE PRECISION array
5443* On entry, CT is an array of dimension at least MAX(M,N,K). CT
5444* holds a copy of the current column of C.
5445*
5446* G (workspace) DOUBLE PRECISION array
5447* On entry, G is an array of dimension at least MAX(M,N,K). G
5448* is used to compute the gauges.
5449*
5450* ERR (global output) DOUBLE PRECISION
5451* On exit, ERR specifies the largest error in absolute value.
5452*
5453* INFO (global output) INTEGER
5454* On exit, if INFO <> 0, the result is less than half accurate.
5455*
5456* -- Written on April 1, 1998 by
5457* Antoine Petitet, University of Tennessee, Knoxville 37996, USA.
5458*
5459* =====================================================================
5460*
5461* .. Parameters ..
5462 INTEGER BLOCK_CYCLIC_2D_INB, CSRC_, CTXT_, DLEN_,
5463 $ DTYPE_, IMB_, INB_, LLD_, MB_, M_, NB_, N_,
5464 $ RSRC_
5465 parameter( block_cyclic_2d_inb = 2, dlen_ = 11,
5466 $ dtype_ = 1, ctxt_ = 2, m_ = 3, n_ = 4,
5467 $ imb_ = 5, inb_ = 6, mb_ = 7, nb_ = 8,
5468 $ rsrc_ = 9, csrc_ = 10, lld_ = 11 )
5469 DOUBLE PRECISION ZERO, ONE
5470 parameter( zero = 0.0d+0, one = 1.0d+0 )
5471* ..
5472* .. Local Scalars ..
5473 LOGICAL COLREP, ROWREP, TRANA, TRANB
5474 INTEGER I, IBB, ICCOL, ICROW, ICURROW, IIC, IN, IOFFA,
5475 $ IOFFB, IOFFC, J, JJC, KK, LDA, LDB, LDC, LDPC,
5476 $ MYCOL, MYROW, NPCOL, NPROW
5477 DOUBLE PRECISION EPS, ERRI
5478* ..
5479* .. External Subroutines ..
5480 EXTERNAL blacs_gridinfo, dgamx2d, igsum2d, pb_infog2l
5481* ..
5482* .. External Functions ..
5483 LOGICAL LSAME
5484 DOUBLE PRECISION PDLAMCH
5485 EXTERNAL lsame, pdlamch
5486* ..
5487* .. Intrinsic Functions ..
5488 INTRINSIC abs, max, min, mod, sqrt
5489* ..
5490* .. Executable Statements ..
5491*
5492 CALL blacs_gridinfo( ictxt, nprow, npcol, myrow, mycol )
5493*
5494 eps = pdlamch( ictxt, 'eps' )
5495*
5496 trana = lsame( transa, 'T' ).OR.lsame( transa, 'C' )
5497 tranb = lsame( transb, 'T' ).OR.lsame( transb, 'C' )
5498*
5499 lda = max( 1, desca( m_ ) )
5500 ldb = max( 1, descb( m_ ) )
5501 ldc = max( 1, descc( m_ ) )
5502*
5503* Compute expected result in C using data in A, B and C.
5504* Compute gauges in G. This part of the computation is performed
5505* by every process in the grid.
5506*
5507 DO 240 j = 1, n
5508*
5509 ioffc = ic + ( jc + j - 2 ) * ldc
5510 DO 10 i = 1, m
5511 ct( i ) = zero
5512 g( i ) = zero
5513 10 CONTINUE
5514*
5515 IF( .NOT.trana .AND. .NOT.tranb ) THEN
5516 DO 30 kk = 1, k
5517 ioffb = ib + kk - 1 + ( jb + j - 2 ) * ldb
5518 DO 20 i = 1, m
5519 ioffa = ia + i - 1 + ( ja + kk - 2 ) * lda
5520 ct( i ) = ct( i ) + a( ioffa ) * b( ioffb )
5521 g( i ) = g( i ) + abs( a( ioffa ) ) *
5522 $ abs( b( ioffb ) )
5523 20 CONTINUE
5524 30 CONTINUE
5525 ELSE IF( trana .AND. .NOT.tranb ) THEN
5526 DO 50 kk = 1, k
5527 ioffb = ib + kk - 1 + ( jb + j - 2 ) * ldb
5528 DO 40 i = 1, m
5529 ioffa = ia + kk - 1 + ( ja + i - 2 ) * lda
5530 ct( i ) = ct( i ) + a( ioffa ) * b( ioffb )
5531 g( i ) = g( i ) + abs( a( ioffa ) ) *
5532 $ abs( b( ioffb ) )
5533 40 CONTINUE
5534 50 CONTINUE
5535 ELSE IF( .NOT.trana .AND. tranb ) THEN
5536 DO 70 kk = 1, k
5537 ioffb = ib + j - 1 + ( jb + kk - 2 ) * ldb
5538 DO 60 i = 1, m
5539 ioffa = ia + i - 1 + ( ja + kk - 2 ) * lda
5540 ct( i ) = ct( i ) + a( ioffa ) * b( ioffb )
5541 g( i ) = g( i ) + abs( a( ioffa ) ) *
5542 $ abs( b( ioffb ) )
5543 60 CONTINUE
5544 70 CONTINUE
5545 ELSE IF( trana .AND. tranb ) THEN
5546 DO 90 kk = 1, k
5547 ioffb = ib + j - 1 + ( jb + kk - 2 ) * ldb
5548 DO 80 i = 1, m
5549 ioffa = ia + kk - 1 + ( ja + i - 2 ) * lda
5550 ct( i ) = ct( i ) + a( ioffa ) * b( ioffb )
5551 g( i ) = g( i ) + abs( a( ioffa ) ) *
5552 $ abs( b( ioffb ) )
5553 80 CONTINUE
5554 90 CONTINUE
5555 END IF
5556*
5557 DO 200 i = 1, m
5558 ct( i ) = alpha*ct( i ) + beta * c( ioffc )
5559 g( i ) = abs( alpha )*g( i ) + abs( beta )*abs( c( ioffc ) )
5560 c( ioffc ) = ct( i )
5561 ioffc = ioffc + 1
5562 200 CONTINUE
5563*
5564* Compute the error ratio for this result.
5565*
5566 err = zero
5567 info = 0
5568 ldpc = descc( lld_ )
5569 ioffc = ic + ( jc + j - 2 ) * ldc
5570 CALL pb_infog2l( ic, jc+j-1, descc, nprow, npcol, myrow, mycol,
5571 $ iic, jjc, icrow, iccol )
5572 icurrow = icrow
5573 rowrep = ( icrow.EQ.-1 )
5574 colrep = ( iccol.EQ.-1 )
5575*
5576 IF( mycol.EQ.iccol .OR. colrep ) THEN
5577*
5578 ibb = descc( imb_ ) - ic + 1
5579 IF( ibb.LE.0 )
5580 $ ibb = ( ( -ibb ) / descc( mb_ ) + 1 )*descc( mb_ ) + ibb
5581 ibb = min( ibb, m )
5582 in = ic + ibb - 1
5583*
5584 DO 210 i = ic, in
5585*
5586 IF( myrow.EQ.icurrow .OR. rowrep ) THEN
5587 erri = abs( pc( iic+(jjc-1)*ldpc ) -
5588 $ c( ioffc ) ) / eps
5589 IF( g( i-ic+1 ).NE.zero )
5590 $ erri = erri / g( i-ic+1 )
5591 err = max( err, erri )
5592 IF( err*sqrt( eps ).GE.one )
5593 $ info = 1
5594 iic = iic + 1
5595 END IF
5596*
5597 ioffc = ioffc + 1
5598*
5599 210 CONTINUE
5600*
5601 icurrow = mod( icurrow+1, nprow )
5602*
5603 DO 230 i = in+1, ic+m-1, descc( mb_ )
5604 ibb = min( ic+m-i, descc( mb_ ) )
5605*
5606 DO 220 kk = 0, ibb-1
5607*
5608 IF( myrow.EQ.icurrow .OR. rowrep ) THEN
5609 erri = abs( pc( iic+(jjc-1)*ldpc ) -
5610 $ c( ioffc ) )/eps
5611 IF( g( i+kk-ic+1 ).NE.zero )
5612 $ erri = erri / g( i+kk-ic+1 )
5613 err = max( err, erri )
5614 IF( err*sqrt( eps ).GE.one )
5615 $ info = 1
5616 iic = iic + 1
5617 END IF
5618*
5619 ioffc = ioffc + 1
5620*
5621 220 CONTINUE
5622*
5623 icurrow = mod( icurrow+1, nprow )
5624*
5625 230 CONTINUE
5626*
5627 END IF
5628*
5629* If INFO = 0, all results are at least half accurate.
5630*
5631 CALL igsum2d( ictxt, 'All', ' ', 1, 1, info, 1, -1, mycol )
5632 CALL dgamx2d( ictxt, 'All', ' ', 1, 1, err, 1, i, j, -1, -1,
5633 $ mycol )
5634 IF( info.NE.0 )
5635 $ GO TO 250
5636*
5637 240 CONTINUE
5638*
5639 250 CONTINUE
5640*
5641 RETURN
5642*
5643* End of PDMMCH
5644*
subroutine pb_infog2l(i, j, desc, nprow, npcol, myrow, mycol, ii, jj, prow, pcol)
Definition pblastst.f:1673
#define max(A, B)
Definition pcgemr.c:180
#define min(A, B)
Definition pcgemr.c:181
double precision function pdlamch(ictxt, cmach)
Definition pdblastst.f:6769
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: