pbstran()

subroutine pbstran	(	integer	icontxt,
		character*1	adist,
		character*1	trans,
		integer	m,
		integer	n,
		integer	nb,
		real, dimension( lda, * )	a,
		integer	lda,
		real	beta,
		real, dimension( ldc, * )	c,
		integer	ldc,
		integer	iarow,
		integer	iacol,
		integer	icrow,
		integer	iccol,
		real, dimension( * )	work
	)

Definition at line 1 of file pbstran.f.

*
*  -- PB-BLAS routine (version 2.1) --
*     University of Tennessee, Knoxville, Oak Ridge National Laboratory.
*     April 28, 1996
*
*     Jaeyoung Choi, Oak Ridge National Laboratory
*     Jack Dongarra, University of Tennessee and Oak Ridge National Lab.
*     David Walker,  Oak Ridge National Laboratory
*
*     .. Scalar Arguments ..
      CHARACTER*1        ADIST, TRANS
      INTEGER            IACOL, IAROW, ICCOL, ICONTXT, ICROW, LDA, LDC,
     $                   M, N, NB
      REAL               BETA
*     ..
*     .. Array Arguments ..
      REAL               A( LDA, * ), C( LDC, * ), WORK( * )
*     ..
*
*  Purpose
*  =======
*
*  PBSTRAN  transposes  a column block to row block, or a row block to
*  column block by reallocating data distribution.
*
*     C := A^T + beta*C, or C := A^C + beta*C
*
*  where A is an M-by-N matrix  and C is an N-by-M matrix, and the size
*  of M or N is limited to its block size NB.
*
*  The first elements  of the matrices A, and C  should  be  located  at
*  the beginnings of their first blocks. (not the middle of the blocks.)
*
*  Parameters
*  ==========
*
*  ICONTXT (input) INTEGER
*          ICONTXT is the BLACS mechanism for partitioning communication
*          space.  A defining property of a context is that a message in
*          a context cannot be sent or received in another context.  The
*          BLACS context includes the definition of a grid, and each
*          process' coordinates in it.
*
*  ADIST  - (input) CHARACTER*1
*           ADIST specifies whether A is a column block or a row block.
*
*              ADIST = 'C',  A is a column block
*              ADIST = 'R',  A is a row block
*
*  TRANS  - (input) CHARACTER*1
*           TRANS specifies whether the transposed format is transpose
*           or conjugate transpose.  If the matrices A and C are real,
*           the argument is ignored.
*
*              TRANS = 'T',  transpose
*              TRANS = 'C',  conjugate transpose
*
*  M      - (input) INTEGER
*           M specifies the (global) number of rows of the matrix (block
*           column or block row) A and of columns of the matrix C.
*           M >= 0.
*
*  N      - (input) INTEGER
*           N specifies the (global) number of columns of the matrix
*           (block column or block row) A  and of columns of the matrix
*           C.  N >= 0.
*
*  NB     - (input) INTEGER
*           NB specifies  the column block size of the matrix A and the
*           row block size of the matrix C when ADIST = 'C'.  Otherwise,
*           it specifies  the row block size of the matrix A and the
*           column block size of the matrix C. NB >= 1.
*
*  A       (input) REAL array of DIMENSION ( LDA, Lx ),
*          where Lx is N  when ADIST = 'C', or Nq when ADIST = 'R'.
*          Before entry with  ADIST = 'C',  the leading Mp by N part of
*          the array A must contain the matrix A, otherwise the leading
*          M by Nq part of the array A  must contain the matrix A.  See
*          parameter details for the values of Mp and Nq.
*
*  LDA     (input) INTEGER
*          LDA specifies the leading dimension of (local) A as declared
*          in the calling (sub) program.  LDA >= MAX(1,Mp) when
*          ADIST = 'C', or LDA >= MAX(1,M) otherwise.
*
*  BETA    (input) REAL
*          BETA specifies scaler beta.
*
*  C       (input/output) REAL array of DIMENSION ( LDC, Lx ),
*          where Lx is Mq when ADIST = 'C', or N when ADIST = 'R'.
*          If ADIST = 'C', the leading N-by-Mq part of the array C
*          contains the (local) matrix C, otherwise the leading
*          Np-by-M part of the array C must contain the (local) matrix
*          C.  C will not be referenced if beta is zero.
*
*  LDC     (input) INTEGER
*          LDC specifies the leading dimension of (local) C as declared
*          in the calling (sub) program. LDC >= MAX(1,N) when ADIST='C',
*          or LDC >= MAX(1,Np) otherwise.
*
*  IAROW   (input) INTEGER
*          IAROW specifies  a row  of the process  template,
*          which holds the first block  of the matrix  A. If A is a row
*          of blocks (ADIST = 'R') and all rows of processes have a copy
*          of A, then set IAROW = -1.
*
*  IACOL   (input) INTEGER
*          IACOL specifies  a column of the process template,
*          which holds  the first block  of the matrix A.  If  A is  a
*          column of blocks (ADIST = 'C') and all columns of processes
*          have a copy of A, then set IACOL = -1.
*
*  ICROW   (input) INTEGER
*          ICROW specifies the current row process which holds
*          the first block  of the matrix C,  which is transposed of A.
*          If C is a row of blocks (ADIST = 'C') and the transposed
*          row block C is distributed all rows of processes, set
*          ICROW = -1.
*
*  ICCOL   (input) INTEGER
*          ICCOL specifies  the current column process which holds
*          the first block of the matrix C,  which is transposed of A.
*          If C is a column of blocks (ADIST = 'R') and the transposed
*          column block C is distributed all columns of processes,
*          set ICCOL = -1.
*
*  WORK    (workspace) REAL array of dimension Size(WORK).
*          It needs extra working space of A'.
*
*  Parameters Details
*  ==================
*
*  Lx      It is  a local portion  of L  owned  by  a process,  (L is
*          replaced by M, or N,  and x is replaced by either p (=NPROW)
*          or q (=NPCOL)).  The value is  determined by  L, LB, x,  and
*          MI, where  LB is  a block size  and  MI is a  row  or column
*          position  in a process template.  Lx is  equal to  or less
*          than Lx0 = CEIL( L, LB*x ) * LB.
*
*  Communication Scheme
*  ====================
*
*  The communication scheme of the routine is set to '1-tree', which is
*  fan-out.  (For details, see BLACS user's guide.)
*
*  Memory Requirement of WORK
*  ==========================
*
*  Mqb  = CEIL( M, NB*NPCOL )
*  Npb  = CEIL( N, NB*NPROW )
*  LCMQ = LCM / NPCOL
*  LCMP = LCM / NPROW
*
*  (1) ADIST = 'C'
*   (a) IACOL != -1
*       Size(WORK) = N * CEIL(Mqb,LCMQ)*NB
*   (b) IACOL = -1
*       Size(WORK) = N * CEIL(Mqb,LCMQ)*NB * MIN(LCMQ,CEIL(M,NB))
*
*  (2) ADIST = 'R'
*   (a) IAROW != -1
*       Size(WORK) = M * CEIL(Npb,LCMP)*NB
*   (b) IAROW = -1
*       Size(WORK) = M * CEIL(Npb,LCMP)*NB * MIN(LCMP,CEIL(N,NB))
*
*  Notes
*  -----
*  More precise space can be computed as
*
*  CEIL(Mqb,LCMQ)*NB => NUMROC( NUMROC(M,NB,0,0,NPCOL), NB, 0, 0, LCMQ )
*  CEIL(Npb,LCMP)*NB => NUMROC( NUMROC(N,NB,0,0,NPROW), NB, 0, 0, LCMP )
*
*  =====================================================================
*
*     ..
*     .. Parameters ..
      REAL               ONE, ZERO
      parameter( one = 1.0e+0, zero = 0.0e+0 )
*     ..
*     .. Local Scalars ..
      LOGICAL            COLFORM, ROWFORM
      INTEGER            I, IDEX, IGD, INFO, JCCOL, JCROW, JDEX, LCM,
     $                   LCMP, LCMQ, MCCOL, MCROW, ML, MP, MQ, MQ0,
     $                   MRCOL, MRROW, MYCOL, MYROW, NP, NP0, NPCOL,
     $                   NPROW, NQ
      REAL               TBETA
*     ..
*     .. External Functions ..
      LOGICAL            LSAME
      INTEGER            ILCM, ICEIL, NUMROC
      EXTERNAL           ilcm, iceil, lsame, numroc
*     ..
*     .. External Subroutines ..
      EXTERNAL           blacs_gridinfo, pbsmatadd, pbstr2af, pbstr2at,
     $                   pbstr2bt, pbstrget, pbstrsrt, pxerbla, sgebr2d,
     $                   sgebs2d, sgerv2d, sgesd2d
*     ..
*     .. Intrinsic Functions ..
      INTRINSIC          max, min, mod
*     ..
*     .. Executable Statements ..
*
*     Quick return if possible.
*
      IF( m.EQ.0 .OR. n.EQ.0 ) RETURN
*
      CALL blacs_gridinfo( icontxt, nprow, npcol, myrow, mycol )
*
      colform = lsame( adist, 'C' )
      rowform = lsame( adist, 'R' )
*
*     Test the input parameters.
*
      info = 0
      IF( ( .NOT.colform ) .AND. ( .NOT.rowform ) ) THEN
         info = 2
      ELSE IF( m .LT.0                            ) THEN
         info = 4
      ELSE IF( n .LT.0                            ) THEN
         info = 5
      ELSE IF( nb.LT.1                            ) THEN
         info = 6
      ELSE IF( iarow.LT.-1 .OR. iarow.GE.nprow .OR.
     $       ( iarow.EQ.-1 .AND. colform )        ) THEN
         info = 12
      ELSE IF( iacol.LT.-1 .OR. iacol.GE.npcol .OR.
     $       ( iacol.EQ.-1 .AND. rowform )        ) THEN
         info = 13
      ELSE IF( icrow.LT.-1 .OR. icrow.GE.nprow .OR.
     $       ( icrow.EQ.-1 .AND. rowform )        ) THEN
         info = 14
      ELSE IF( iccol.LT.-1 .OR. iccol.GE.npcol .OR.
     $       ( iccol.EQ.-1 .AND. colform )        ) THEN
         info = 15
      END IF
*
   10 CONTINUE
      IF( info .NE. 0 ) THEN
         CALL pxerbla( icontxt, 'PBSTRAN ', info )
         RETURN
      END IF
*
*     Start the operations.
*
*     LCM : the least common multiple of NPROW and NPCOL
*
      lcm  = ilcm( nprow, npcol )
      lcmp = lcm   / nprow
      lcmq = lcm   / npcol
      igd  = npcol / lcmp
*
*     When A is a column block
*
      IF( colform ) THEN
*
*       Form  C <== A'  ( A is a column block )
*                                         _
*                                        | |
*                                        | |
*            _____________               | |
*           |______C______|     <==      |A|
*                                        | |
*                                        | |
*                                        |_|
*
*       MRROW : row relative position in template from IAROW
*       MRCOL : column relative position in template from ICCOL
*
        mrrow = mod( nprow+myrow-iarow, nprow )
        mrcol = mod( npcol+mycol-iccol, npcol )
        jcrow = icrow
        IF( icrow.EQ.-1 ) jcrow = iarow
*
        mp  = numroc( m, nb, myrow, iarow, nprow )
        mq  = numroc( m, nb, mycol, iccol, npcol )
        mq0 = numroc( numroc(m, nb, 0, 0, npcol), nb, 0, 0, lcmq )
*
        IF( lda.LT.mp .AND.
     $         ( iacol.EQ.mycol .OR. iacol.EQ.-1 ) ) THEN
           info = 8
        ELSE IF( ldc.LT.n .AND.
     $         ( icrow.EQ.myrow .OR. icrow.EQ.-1 ) ) THEN
           info = 11
        END IF
        IF( info.NE.0 ) GO TO 10
*
*       When a column process of IACOL has a column block A,
*
        IF( iacol.GE.0 ) THEN
          tbeta = zero
          IF( myrow.EQ.jcrow ) tbeta = beta
*
          DO 20 i = 0, min( lcm, iceil(m,nb) ) - 1
            mcrow = mod( mod(i, nprow) + iarow, nprow )
            mccol = mod( mod(i, npcol) + iccol, npcol )
            IF( lcmq.EQ.1 )  mq0 = numroc( m, nb, i, 0, npcol )
            jdex = (i/npcol) * nb
*
*           A source node copies the blocks to WORK, and send it
*
            IF( myrow.EQ.mcrow .AND. mycol.EQ.iacol ) THEN
*
*             The source node is a destination node
*
              idex = (i/nprow) * nb
              IF( myrow.EQ.jcrow .AND. mycol.EQ.mccol ) THEN
                CALL pbstr2at( icontxt, 'Col', trans, mp-idex, n, nb,
     $                         a(idex+1,1), lda, tbeta, c(1,jdex+1),
     $                         ldc, lcmp, lcmq )
*
*             The source node sends blocks to a destination node
*
              ELSE
                CALL pbstr2bt( icontxt, 'Col', trans, mp-idex, n, nb,
     $                         a(idex+1,1), lda, zero, work, n,
     $                         lcmp*nb )
                CALL sgesd2d( icontxt, n, mq0, work, n, jcrow, mccol )
              END IF
*
*           A destination node receives the copied blocks
*
            ELSE IF( myrow.EQ.jcrow .AND. mycol.EQ.mccol ) THEN
              IF( lcmq.EQ.1 .AND. tbeta.EQ.zero ) THEN
                CALL sgerv2d( icontxt, n, mq0, c, ldc, mcrow, iacol )
              ELSE
                CALL sgerv2d( icontxt, n, mq0, work, n, mcrow, iacol )
                CALL pbstr2af( icontxt, 'Row', n, mq-jdex, nb, work, n,
     $                         tbeta, c(1,jdex+1), ldc, lcmp, lcmq,
     $                         mq0 )
              END IF
            END IF
   20     CONTINUE
*
*         Broadcast a row block of C in each column of template
*
          IF( icrow.EQ.-1 ) THEN
            IF( myrow.EQ.jcrow ) THEN
              CALL sgebs2d( icontxt, 'Col', '1-tree', n, mq, c, ldc )
            ELSE
              CALL sgebr2d( icontxt, 'Col', '1-tree', n, mq, c, ldc,
     $                      jcrow, mycol )
            END IF
          END IF
*
*       When all column procesors have a copy of the column block A,
*
        ELSE
          IF( lcmq.EQ.1 ) mq0 = mq
*
*         Processors, which have diagonal blocks of A, copy them to
*         WORK array in transposed form
*
          DO 30 i = 0, lcmp-1
            IF( mrcol.EQ.mod( nprow*i+mrrow, npcol ) ) THEN
              IF( lcmq.EQ.1.AND.(icrow.EQ.-1.OR.icrow.EQ.myrow) ) THEN
                 CALL pbstr2bt( icontxt, 'Col', trans, mp-i*nb, n, nb,
     $                          a(i*nb+1,1), lda, beta, c, ldc,
     $                          lcmp*nb )
              ELSE
                 CALL pbstr2bt( icontxt, 'Col', trans, mp-i*nb, n, nb,
     $                          a(i*nb+1,1), lda, zero, work, n,
     $                          lcmp*nb )
              END IF
            END IF
   30     CONTINUE
*
*         Get diagonal blocks of A for each column of the template
*
          mcrow = mod( mod(mrcol,nprow)+iarow, nprow )
          IF( lcmq.GT.1 ) THEN
            mccol = mod( npcol+mycol-iccol, npcol )
            CALL pbstrget( icontxt, 'Row', n, mq0, iceil(m,nb), work, n,
     $                     mcrow,  mccol, igd, myrow, mycol, nprow,
     $                     npcol )
          END IF
*
*         Broadcast a row block of WORK in every row of template
*
          IF( icrow.EQ.-1 ) THEN
            IF( myrow.EQ.mcrow ) THEN
              IF( lcmq.GT.1 )
     $          CALL pbstrsrt( icontxt, 'Row', n, mq, nb, work, n, beta,
     $                         c, ldc, lcmp, lcmq, mq0 )
              CALL sgebs2d( icontxt, 'Col', '1-tree', n, mq, c, ldc )
            ELSE
              CALL sgebr2d( icontxt, 'Col', '1-tree', n, mq, c, ldc,
     $                      mcrow, mycol )
            END IF
*
*         Send a row block of WORK to the destination row
*
          ELSE
            IF( lcmq.EQ.1 ) THEN
              IF( myrow.EQ.mcrow ) THEN
                IF( myrow.NE.icrow )
     $            CALL sgesd2d( icontxt, n, mq, work, n, icrow, mycol )
              ELSE IF( myrow.EQ.icrow ) THEN
                IF( beta.EQ.zero ) THEN
                  CALL sgerv2d( icontxt, n, mq, c, ldc, mcrow, mycol )
                ELSE
                  CALL sgerv2d( icontxt, n, mq, work, n, mcrow, mycol )
                  CALL pbsmatadd( icontxt, 'G', n, mq, one, work, n,
     $                            beta, c, ldc )
                END IF
              END IF
*
            ELSE
              ml = mq0 * min( lcmq, max(0,iceil(m,nb)-mccol) )
              IF( myrow.EQ.mcrow ) THEN
                IF( myrow.NE.icrow )
     $            CALL sgesd2d( icontxt, n, ml, work, n, icrow, mycol )
              ELSE IF( myrow.EQ.icrow ) THEN
                CALL sgerv2d( icontxt, n, ml, work, n, mcrow, mycol )
              END IF
*
              IF( myrow.EQ.icrow )
     $          CALL pbstrsrt( icontxt, 'Row', n, mq, nb, work, n, beta,
     $                         c, ldc, lcmp, lcmq, mq0 )
            END IF
          END IF
*
        END IF
*
*     When A is a row block
*
      ELSE
*
*        Form  C <== A'  ( A is a row block )
*            _
*           | |
*           | |
*           | |                _____________
*           |C|      <==      |______A______|
*           | |
*           | |
*           |_|
*
*        MRROW : row relative position in template from ICROW
*        MRCOL : column relative position in template from IACOL
*
         mrrow = mod( nprow+myrow-icrow, nprow )
         mrcol = mod( npcol+mycol-iacol, npcol )
         jccol = iccol
         IF( iccol.EQ.-1 ) jccol = iacol
*
         np  = numroc( n, nb, myrow, icrow, nprow )
         nq  = numroc( n, nb, mycol, iacol, npcol )
         np0 = numroc( numroc(n, nb, 0, 0, nprow), nb, 0, 0, lcmp )
*
         IF( lda.LT.m .AND.
     $          ( iarow.EQ.myrow .OR. iarow.EQ.-1 ) ) THEN
            info = 8
         ELSE IF( ldc.LT.np .AND.
     $          ( iccol.EQ.mycol .OR. iccol.EQ.-1 ) ) THEN
            info = 11
         END IF
         IF( info.NE.0 ) GO TO 10
*
*        When a row process of IAROW has a row block A,
*
         IF( iarow.GE.0 ) THEN
           tbeta = zero
           IF( mycol.EQ.jccol ) tbeta = beta
*
           DO 40 i = 0, min( lcm, iceil(n,nb) ) - 1
             mcrow = mod( mod(i, nprow) + icrow, nprow )
             mccol = mod( mod(i, npcol) + iacol, npcol )
             IF( lcmp.EQ.1 )  np0 = numroc( n, nb, i, 0, nprow )
             idex = (i/nprow) * nb
*
*            A source node copies the blocks to WORK, and send it
*
             IF( myrow.EQ.iarow .AND. mycol.EQ.mccol ) THEN
*
*              The source node is a destination node
*
               jdex = (i/npcol) * nb
               IF( myrow.EQ.mcrow .AND. mycol.EQ.jccol ) THEN
                 CALL pbstr2at( icontxt, 'Row', trans, m, nq-jdex, nb,
     $                          a(1,jdex+1), lda, tbeta, c(idex+1,1),
     $                          ldc, lcmp, lcmq )
*
*              The source node sends blocks to a destination node
*
               ELSE
                 CALL pbstr2bt( icontxt, 'Row', trans, m, nq-jdex, nb,
     $                          a(1,jdex+1), lda, zero, work, np0,
     $                          lcmq*nb )
                 CALL sgesd2d( icontxt, np0, m, work, np0,
     $                         mcrow, jccol )
               END IF
*
*           A destination node receives the copied blocks
*
            ELSE IF( myrow.EQ.mcrow .AND. mycol.EQ.jccol ) THEN
              IF( lcmp.EQ.1 .AND. tbeta.EQ.zero ) THEN
                CALL sgerv2d( icontxt, np0, m, c, ldc, iarow, mccol )
              ELSE
                CALL sgerv2d( icontxt, np0, m, work, np0, iarow, mccol )
                CALL pbstr2af( icontxt, 'Col', np-idex, m, nb, work,
     $                         np0, tbeta, c(idex+1,1), ldc, lcmp, lcmq,
     $                         np0 )
              END IF
            END IF
   40     CONTINUE
*
*         Broadcast a column block of WORK in each row of template
*
          IF( iccol.EQ.-1 ) THEN
            IF( mycol.EQ.jccol ) THEN
              CALL sgebs2d( icontxt, 'Row', '1-tree', np, m, c, ldc )
            ELSE
              CALL sgebr2d( icontxt, 'Row', '1-tree', np, m, c, ldc,
     $                       myrow, jccol )
            END IF
          END IF
*
*       When all row procesors have a copy of the row block A,
*
        ELSE
          IF( lcmp.EQ.1 ) np0 = np
*
*         Processors, which have diagonal blocks of A, copy them to
*         WORK array in transposed form
*
          DO 50 i = 0, lcmq-1
            IF( mrrow.EQ.mod(npcol*i+mrcol, nprow) ) THEN
              IF( lcmp.EQ.1.AND.(iccol.EQ.-1.OR.iccol.EQ.mycol) ) THEN
                CALL pbstr2bt( icontxt, 'Row', trans, m, nq-i*nb, nb,
     $                         a(1,i*nb+1), lda, beta, c, ldc,
     $                         lcmq*nb )
              ELSE
                CALL pbstr2bt( icontxt, 'Row', trans, m, nq-i*nb, nb,
     $                         a(1,i*nb+1), lda, zero, work, np0,
     $                         lcmq*nb )
              END IF
            END IF
   50     CONTINUE
*
*         Get diagonal blocks of A for each row of the template
*
          mccol = mod( mod(mrrow, npcol)+iacol, npcol )
          IF( lcmp.GT.1 ) THEN
            mcrow = mod( nprow+myrow-icrow, nprow )
            CALL pbstrget( icontxt, 'Col', np0, m, iceil(n,nb), work,
     $                     np0, mcrow, mccol, igd, myrow, mycol, nprow,
     $                     npcol )
          END IF
*
*         Broadcast a column block of WORK in every column of template
*
          IF( iccol.EQ.-1 ) THEN
            IF( mycol.EQ.mccol ) THEN
              IF( lcmp.GT.1 )
     $          CALL pbstrsrt( icontxt, 'Col', np, m, nb, work, np0,
     $                         beta, c, ldc, lcmp, lcmq, np0 )
              CALL sgebs2d( icontxt, 'Row', '1-tree', np, m, c, ldc )
            ELSE
              CALL sgebr2d( icontxt, 'Row', '1-tree', np, m, c, ldc,
     $                       myrow, mccol )
            END IF
*
*         Send a column block of WORK to the destination column
*
          ELSE
            IF( lcmp.EQ.1 ) THEN
              IF( mycol.EQ.mccol ) THEN
                IF( mycol.NE.iccol )
     $            CALL sgesd2d( icontxt, np, m, work, np, myrow, iccol )
              ELSE IF( mycol.EQ.iccol ) THEN
                IF( beta.EQ.zero ) THEN
                  CALL sgerv2d( icontxt, np, m, c, ldc, myrow, mccol )
                ELSE
                  CALL sgerv2d( icontxt, np, m, work, np, myrow, mccol )
                  CALL pbsmatadd( icontxt, 'G', np, m, one, work, np,
     $                            beta, c, ldc )
                END IF
              END IF
*
            ELSE
              ml = m * min( lcmp, max( 0, iceil(n,nb) - mcrow ) )
              IF( mycol.EQ.mccol ) THEN
                IF( mycol.NE.iccol )
     $            CALL sgesd2d( icontxt, np0, ml, work, np0,
     $                          myrow, iccol )
              ELSE IF( mycol.EQ.iccol ) THEN
                CALL sgerv2d( icontxt, np0, ml, work, np0,
     $                        myrow, mccol )
              END IF
*
              IF( mycol.EQ.iccol )
     $          CALL pbstrsrt( icontxt, 'Col', np, m, nb, work, np0,
     $                         beta, c, ldc, lcmp, lcmq, np0 )
            END IF
          END IF
*
        END IF
      END IF
*
      RETURN
*
*     End of PBSTRAN
*

Here is the call graph for this function:

Here is the caller graph for this function: