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

◆ pzmmch2()

subroutine pzmmch2 ( integer  ictxt,
character*1  uplo,
character*1  trans,
integer  n,
integer  k,
complex*16  alpha,
complex*16, dimension( * )  a,
integer  ia,
integer  ja,
integer, dimension( * )  desca,
complex*16, dimension( * )  b,
integer  ib,
integer  jb,
integer, dimension( * )  descb,
complex*16  beta,
complex*16, dimension( * )  c,
complex*16, dimension( * )  pc,
integer  ic,
integer  jc,
integer, dimension( * )  descc,
complex*16, dimension( * )  ct,
double precision, dimension( * )  g,
double precision  err,
integer  info 
)

Definition at line 6166 of file pzblastst.f.

6169*
6170* -- PBLAS test routine (version 2.0) --
6171* University of Tennessee, Knoxville, Oak Ridge National Laboratory,
6172* and University of California, Berkeley.
6173* April 1, 1998
6174*
6175* .. Scalar Arguments ..
6176 CHARACTER*1 TRANS, UPLO
6177 INTEGER IA, IB, IC, ICTXT, INFO, JA, JB, JC, K, N
6178 DOUBLE PRECISION ERR
6179 COMPLEX*16 ALPHA, BETA
6180* ..
6181* .. Array Arguments ..
6182 INTEGER DESCA( * ), DESCB( * ), DESCC( * )
6183 DOUBLE PRECISION G( * )
6184 COMPLEX*16 A( * ), B( * ), C( * ), CT( * ),
6185 $ PC( * )
6186* ..
6187*
6188* Purpose
6189* =======
6190*
6191* PZMMCH2 checks the results of the computational tests.
6192*
6193* Notes
6194* =====
6195*
6196* A description vector is associated with each 2D block-cyclicly dis-
6197* tributed matrix. This vector stores the information required to
6198* establish the mapping between a matrix entry and its corresponding
6199* process and memory location.
6200*
6201* In the following comments, the character _ should be read as
6202* "of the distributed matrix". Let A be a generic term for any 2D
6203* block cyclicly distributed matrix. Its description vector is DESCA:
6204*
6205* NOTATION STORED IN EXPLANATION
6206* ---------------- --------------- ------------------------------------
6207* DTYPE_A (global) DESCA( DTYPE_ ) The descriptor type.
6208* CTXT_A (global) DESCA( CTXT_ ) The BLACS context handle, indicating
6209* the NPROW x NPCOL BLACS process grid
6210* A is distributed over. The context
6211* itself is global, but the handle
6212* (the integer value) may vary.
6213* M_A (global) DESCA( M_ ) The number of rows in the distribu-
6214* ted matrix A, M_A >= 0.
6215* N_A (global) DESCA( N_ ) The number of columns in the distri-
6216* buted matrix A, N_A >= 0.
6217* IMB_A (global) DESCA( IMB_ ) The number of rows of the upper left
6218* block of the matrix A, IMB_A > 0.
6219* INB_A (global) DESCA( INB_ ) The number of columns of the upper
6220* left block of the matrix A,
6221* INB_A > 0.
6222* MB_A (global) DESCA( MB_ ) The blocking factor used to distri-
6223* bute the last M_A-IMB_A rows of A,
6224* MB_A > 0.
6225* NB_A (global) DESCA( NB_ ) The blocking factor used to distri-
6226* bute the last N_A-INB_A columns of
6227* A, NB_A > 0.
6228* RSRC_A (global) DESCA( RSRC_ ) The process row over which the first
6229* row of the matrix A is distributed,
6230* NPROW > RSRC_A >= 0.
6231* CSRC_A (global) DESCA( CSRC_ ) The process column over which the
6232* first column of A is distributed.
6233* NPCOL > CSRC_A >= 0.
6234* LLD_A (local) DESCA( LLD_ ) The leading dimension of the local
6235* array storing the local blocks of
6236* the distributed matrix A,
6237* IF( Lc( 1, N_A ) > 0 )
6238* LLD_A >= MAX( 1, Lr( 1, M_A ) )
6239* ELSE
6240* LLD_A >= 1.
6241*
6242* Let K be the number of rows of a matrix A starting at the global in-
6243* dex IA,i.e, A( IA:IA+K-1, : ). Lr( IA, K ) denotes the number of rows
6244* that the process of row coordinate MYROW ( 0 <= MYROW < NPROW ) would
6245* receive if these K rows were distributed over NPROW processes. If K
6246* is the number of columns of a matrix A starting at the global index
6247* JA, i.e, A( :, JA:JA+K-1, : ), Lc( JA, K ) denotes the number of co-
6248* lumns that the process MYCOL ( 0 <= MYCOL < NPCOL ) would receive if
6249* these K columns were distributed over NPCOL processes.
6250*
6251* The values of Lr() and Lc() may be determined via a call to the func-
6252* tion PB_NUMROC:
6253* Lr( IA, K ) = PB_NUMROC( K, IA, IMB_A, MB_A, MYROW, RSRC_A, NPROW )
6254* Lc( JA, K ) = PB_NUMROC( K, JA, INB_A, NB_A, MYCOL, CSRC_A, NPCOL )
6255*
6256* Arguments
6257* =========
6258*
6259* ICTXT (local input) INTEGER
6260* On entry, ICTXT specifies the BLACS context handle, indica-
6261* ting the global context of the operation. The context itself
6262* is global, but the value of ICTXT is local.
6263*
6264* UPLO (global input) CHARACTER*1
6265* On entry, UPLO specifies which part of C should contain the
6266* result.
6267*
6268* TRANS (global input) CHARACTER*1
6269* On entry, TRANS specifies whether the matrices A and B have
6270* to be transposed or not before computing the matrix-matrix
6271* product.
6272*
6273* N (global input) INTEGER
6274* On entry, N specifies the order the submatrix operand C. N
6275* must be at least zero.
6276*
6277* K (global input) INTEGER
6278* On entry, K specifies the number of columns (resp. rows) of A
6279* and B when TRANS = 'N' (resp. TRANS <> 'N'). K must be at
6280* least zero.
6281*
6282* ALPHA (global input) COMPLEX*16
6283* On entry, ALPHA specifies the scalar alpha.
6284*
6285* A (local input) COMPLEX*16 array
6286* On entry, A is an array of dimension (DESCA( M_ ),*). This
6287* array contains a local copy of the initial entire matrix PA.
6288*
6289* IA (global input) INTEGER
6290* On entry, IA specifies A's global row index, which points to
6291* the beginning of the submatrix sub( A ).
6292*
6293* JA (global input) INTEGER
6294* On entry, JA specifies A's global column index, which points
6295* to the beginning of the submatrix sub( A ).
6296*
6297* DESCA (global and local input) INTEGER array
6298* On entry, DESCA is an integer array of dimension DLEN_. This
6299* is the array descriptor for the matrix A.
6300*
6301* B (local input) COMPLEX*16 array
6302* On entry, B is an array of dimension (DESCB( M_ ),*). This
6303* array contains a local copy of the initial entire matrix PB.
6304*
6305* IB (global input) INTEGER
6306* On entry, IB specifies B's global row index, which points to
6307* the beginning of the submatrix sub( B ).
6308*
6309* JB (global input) INTEGER
6310* On entry, JB specifies B's global column index, which points
6311* to the beginning of the submatrix sub( B ).
6312*
6313* DESCB (global and local input) INTEGER array
6314* On entry, DESCB is an integer array of dimension DLEN_. This
6315* is the array descriptor for the matrix B.
6316*
6317* BETA (global input) COMPLEX*16
6318* On entry, BETA specifies the scalar beta.
6319*
6320* C (local input/local output) COMPLEX*16 array
6321* On entry, C is an array of dimension (DESCC( M_ ),*). This
6322* array contains a local copy of the initial entire matrix PC.
6323*
6324* PC (local input) COMPLEX*16 array
6325* On entry, PC is an array of dimension (DESCC( LLD_ ),*). This
6326* array contains the local pieces of the matrix PC.
6327*
6328* IC (global input) INTEGER
6329* On entry, IC specifies C's global row index, which points to
6330* the beginning of the submatrix sub( C ).
6331*
6332* JC (global input) INTEGER
6333* On entry, JC specifies C's global column index, which points
6334* to the beginning of the submatrix sub( C ).
6335*
6336* DESCC (global and local input) INTEGER array
6337* On entry, DESCC is an integer array of dimension DLEN_. This
6338* is the array descriptor for the matrix C.
6339*
6340* CT (workspace) COMPLEX*16 array
6341* On entry, CT is an array of dimension at least MAX(M,N,K). CT
6342* holds a copy of the current column of C.
6343*
6344* G (workspace) DOUBLE PRECISION array
6345* On entry, G is an array of dimension at least MAX(M,N,K). G
6346* is used to compute the gauges.
6347*
6348* ERR (global output) DOUBLE PRECISION
6349* On exit, ERR specifies the largest error in absolute value.
6350*
6351* INFO (global output) INTEGER
6352* On exit, if INFO <> 0, the result is less than half accurate.
6353*
6354* -- Written on April 1, 1998 by
6355* Antoine Petitet, University of Tennessee, Knoxville 37996, USA.
6356*
6357* =====================================================================
6358*
6359* .. Parameters ..
6360 INTEGER BLOCK_CYCLIC_2D_INB, CSRC_, CTXT_, DLEN_,
6361 $ DTYPE_, IMB_, INB_, LLD_, MB_, M_, NB_, N_,
6362 $ RSRC_
6363 parameter( block_cyclic_2d_inb = 2, dlen_ = 11,
6364 $ dtype_ = 1, ctxt_ = 2, m_ = 3, n_ = 4,
6365 $ imb_ = 5, inb_ = 6, mb_ = 7, nb_ = 8,
6366 $ rsrc_ = 9, csrc_ = 10, lld_ = 11 )
6367 DOUBLE PRECISION RZERO, RONE
6368 parameter( rzero = 0.0d+0, rone = 1.0d+0 )
6369 COMPLEX*16 ZERO
6370 parameter( zero = ( 0.0d+0, 0.0d+0 ) )
6371* ..
6372* .. Local Scalars ..
6373 LOGICAL COLREP, HTRAN, NOTRAN, ROWREP, TRAN, UPPER
6374 INTEGER I, IBB, IBEG, ICCOL, ICROW, ICURROW, IEND, IIC,
6375 $ IN, IOFFAK, IOFFAN, IOFFBK, IOFFBN, IOFFC, J,
6376 $ JJC, KK, LDA, LDB, LDC, LDPC, MYCOL, MYROW,
6377 $ NPCOL, NPROW
6378 DOUBLE PRECISION EPS, ERRI
6379 COMPLEX*16 Z
6380* ..
6381* .. External Subroutines ..
6382 EXTERNAL blacs_gridinfo, dgamx2d, igsum2d, pb_infog2l
6383* ..
6384* .. External Functions ..
6385 LOGICAL LSAME
6386 DOUBLE PRECISION PDLAMCH
6387 EXTERNAL lsame, pdlamch
6388* ..
6389* .. Intrinsic Functions ..
6390 INTRINSIC abs, dble, dconjg, dimag, max, min, mod, sqrt
6391* ..
6392* .. Statement Functions ..
6393 DOUBLE PRECISION ABS1
6394 abs1( z ) = abs( dble( z ) ) + abs( dimag( z ) )
6395* ..
6396* .. Executable Statements ..
6397*
6398 CALL blacs_gridinfo( ictxt, nprow, npcol, myrow, mycol )
6399*
6400 eps = pdlamch( ictxt, 'eps' )
6401*
6402 upper = lsame( uplo, 'U' )
6403 htran = lsame( trans, 'H' )
6404 notran = lsame( trans, 'N' )
6405 tran = lsame( trans, 'T' )
6406*
6407 lda = max( 1, desca( m_ ) )
6408 ldb = max( 1, descb( m_ ) )
6409 ldc = max( 1, descc( m_ ) )
6410*
6411* Compute expected result in C using data in A, B and C.
6412* Compute gauges in G. This part of the computation is performed
6413* by every process in the grid.
6414*
6415 DO 140 j = 1, n
6416*
6417 IF( upper ) THEN
6418 ibeg = 1
6419 iend = j
6420 ELSE
6421 ibeg = j
6422 iend = n
6423 END IF
6424*
6425 DO 10 i = 1, n
6426 ct( i ) = zero
6427 g( i ) = rzero
6428 10 CONTINUE
6429*
6430 IF( notran ) THEN
6431 DO 30 kk = 1, k
6432 ioffak = ia + j - 1 + ( ja + kk - 2 ) * lda
6433 ioffbk = ib + j - 1 + ( jb + kk - 2 ) * ldb
6434 DO 20 i = ibeg, iend
6435 ioffan = ia + i - 1 + ( ja + kk - 2 ) * lda
6436 ioffbn = ib + i - 1 + ( jb + kk - 2 ) * ldb
6437 ct( i ) = ct( i ) + alpha * (
6438 $ a( ioffan ) * b( ioffbk ) +
6439 $ b( ioffbn ) * a( ioffak ) )
6440 g( i ) = g( i ) + abs( alpha ) * (
6441 $ abs1( a( ioffan ) ) * abs1( b( ioffbk ) ) +
6442 $ abs1( b( ioffbn ) ) * abs1( a( ioffak ) ) )
6443 20 CONTINUE
6444 30 CONTINUE
6445 ELSE IF( tran ) THEN
6446 DO 50 kk = 1, k
6447 ioffak = ia + kk - 1 + ( ja + j - 2 ) * lda
6448 ioffbk = ib + kk - 1 + ( jb + j - 2 ) * ldb
6449 DO 40 i = ibeg, iend
6450 ioffan = ia + kk - 1 + ( ja + i - 2 ) * lda
6451 ioffbn = ib + kk - 1 + ( jb + i - 2 ) * ldb
6452 ct( i ) = ct( i ) + alpha * (
6453 $ a( ioffan ) * b( ioffbk ) +
6454 $ b( ioffbn ) * a( ioffak ) )
6455 g( i ) = g( i ) + abs( alpha ) * (
6456 $ abs1( a( ioffan ) ) * abs1( b( ioffbk ) ) +
6457 $ abs1( b( ioffbn ) ) * abs1( a( ioffak ) ) )
6458 40 CONTINUE
6459 50 CONTINUE
6460 ELSE IF( htran ) THEN
6461 DO 70 kk = 1, k
6462 ioffak = ia + j - 1 + ( ja + kk - 2 ) * lda
6463 ioffbk = ib + j - 1 + ( jb + kk - 2 ) * ldb
6464 DO 60 i = ibeg, iend
6465 ioffan = ia + i - 1 + ( ja + kk - 2 ) * lda
6466 ioffbn = ib + i - 1 + ( jb + kk - 2 ) * ldb
6467 ct( i ) = ct( i ) +
6468 $ alpha * a( ioffan ) * dconjg( b( ioffbk ) ) +
6469 $ b( ioffbn ) * dconjg( alpha * a( ioffak ) )
6470 g( i ) = g( i ) + abs1( alpha ) * (
6471 $ abs1( a( ioffan ) ) * abs1( b( ioffbk ) ) +
6472 $ abs1( b( ioffbn ) ) * abs1( a( ioffak ) ) )
6473 60 CONTINUE
6474 70 CONTINUE
6475 ELSE
6476 DO 90 kk = 1, k
6477 ioffak = ia + kk - 1 + ( ja + j - 2 ) * lda
6478 ioffbk = ib + kk - 1 + ( jb + j - 2 ) * ldb
6479 DO 80 i = ibeg, iend
6480 ioffan = ia + kk - 1 + ( ja + i - 2 ) * lda
6481 ioffbn = ib + kk - 1 + ( jb + i - 2 ) * ldb
6482 ct( i ) = ct( i ) +
6483 $ alpha * dconjg( a( ioffan ) ) * b( ioffbk ) +
6484 $ dconjg( alpha * b( ioffbn ) ) * a( ioffak )
6485 g( i ) = g( i ) + abs1( alpha ) * (
6486 $ abs1( dconjg( a( ioffan ) ) * b( ioffbk ) ) +
6487 $ abs1( dconjg( b( ioffbn ) ) * a( ioffak ) ) )
6488 80 CONTINUE
6489 90 CONTINUE
6490 END IF
6491*
6492 ioffc = ic + ibeg - 1 + ( jc + j - 2 ) * ldc
6493*
6494 DO 100 i = ibeg, iend
6495 ct( i ) = ct( i ) + beta * c( ioffc )
6496 g( i ) = g( i ) + abs1( beta )*abs1( c( ioffc ) )
6497 c( ioffc ) = ct( i )
6498 ioffc = ioffc + 1
6499 100 CONTINUE
6500*
6501* Compute the error ratio for this result.
6502*
6503 err = rzero
6504 info = 0
6505 ldpc = descc( lld_ )
6506 ioffc = ic + ( jc + j - 2 ) * ldc
6507 CALL pb_infog2l( ic, jc+j-1, descc, nprow, npcol, myrow, mycol,
6508 $ iic, jjc, icrow, iccol )
6509 icurrow = icrow
6510 rowrep = ( icrow.EQ.-1 )
6511 colrep = ( iccol.EQ.-1 )
6512*
6513 IF( mycol.EQ.iccol .OR. colrep ) THEN
6514*
6515 ibb = descc( imb_ ) - ic + 1
6516 IF( ibb.LE.0 )
6517 $ ibb = ( ( -ibb ) / descc( mb_ ) + 1 )*descc( mb_ ) + ibb
6518 ibb = min( ibb, n )
6519 in = ic + ibb - 1
6520*
6521 DO 110 i = ic, in
6522*
6523 IF( myrow.EQ.icurrow .OR. rowrep ) THEN
6524 erri = abs( pc( iic+(jjc-1)*ldpc ) -
6525 $ c( ioffc ) ) / eps
6526 IF( g( i-ic+1 ).NE.rzero )
6527 $ erri = erri / g( i-ic+1 )
6528 err = max( err, erri )
6529 IF( err*sqrt( eps ).GE.rone )
6530 $ info = 1
6531 iic = iic + 1
6532 END IF
6533*
6534 ioffc = ioffc + 1
6535*
6536 110 CONTINUE
6537*
6538 icurrow = mod( icurrow+1, nprow )
6539*
6540 DO 130 i = in+1, ic+n-1, descc( mb_ )
6541 ibb = min( ic+n-i, descc( mb_ ) )
6542*
6543 DO 120 kk = 0, ibb-1
6544*
6545 IF( myrow.EQ.icurrow .OR. rowrep ) THEN
6546 erri = abs( pc( iic+(jjc-1)*ldpc ) -
6547 $ c( ioffc ) )/eps
6548 IF( g( i+kk-ic+1 ).NE.rzero )
6549 $ erri = erri / g( i+kk-ic+1 )
6550 err = max( err, erri )
6551 IF( err*sqrt( eps ).GE.rone )
6552 $ info = 1
6553 iic = iic + 1
6554 END IF
6555*
6556 ioffc = ioffc + 1
6557*
6558 120 CONTINUE
6559*
6560 icurrow = mod( icurrow+1, nprow )
6561*
6562 130 CONTINUE
6563*
6564 END IF
6565*
6566* If INFO = 0, all results are at least half accurate.
6567*
6568 CALL igsum2d( ictxt, 'All', ' ', 1, 1, info, 1, -1, mycol )
6569 CALL dgamx2d( ictxt, 'All', ' ', 1, 1, err, 1, i, j, -1, -1,
6570 $ mycol )
6571 IF( info.NE.0 )
6572 $ GO TO 150
6573*
6574 140 CONTINUE
6575*
6576 150 CONTINUE
6577*
6578 RETURN
6579*
6580* End of PZMMCH2
6581*
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: