|
SRC\pdlasrt.f |
|
| #lines: 254 size: 8 Kb creation: 18/01/2006 23:36:04 last modification: 08/05/2008 18:37:55 attribute: ARCH Find Reload | |
1: 2: 3: 4: 5: 6: 7: 8: 9: 10: 11: 12: 13: 14: 15: 16: 17: 18: 19: 20: 21: 22: 23: 24: 25: 26: 27: 28: 29: 30: 31: 32: 33: 34: 35: 36: 37: 38: 39: 40: 41: 42: 43: 44: 45: 46: 47: 48: 49: 50: 51: 52: 53: 54: 55: 56: 57: 58: 59: 60: 61: 62: 63: 64: 65: 66: 67: 68: 69: 70: 71: 72: 73: 74: 75: 76: 77: 78: 79: 80: 81: 82: 83: 84: 85: 86: 87: 88: 89: 90: 91: 92: 93: 94: 95: 96: 97: 98: 99: 100: 101: 102: 103: 104: 105: 106: 107: 108: 109: 110: 111: 112: 113: 114: 115: 116: 117: 118: 119: 120: 121: 122: 123: 124: 125: 126: 127: 128: 129: 130: 131: 132: 133: 134: 135: 136: 137: 138: 139: 140: 141: 142: 143: 144: 145: 146: 147: 148: 149: 150: 151: 152: 153: 154: 155: 156: 157: 158: 159: 160: 161: 162: 163: 164: 165: 166: 167: 168: 169: 170: 171: 172: 173: 174: 175: 176: 177: 178: 179: 180: 181: 182: 183: 184: 185: 186: 187: 188: 189: 190: 191: 192: 193: 194: 195: 196: 197: 198: 199: 200: 201: 202: 203: 204: 205: 206: 207: 208: 209: 210: 211: 212: 213: 214: 215: 216: 217: 218: 219: 220: 221: 222: 223: 224: 225: 226: 227: 228: 229: 230: 231: 232: 233: 234: 235: 236: 237: 238: 239: 240: 241: 242: 243: 244: 245: 246: 247: 248: 249: 250: 251: 252: 253: 254: |
SUBROUTINE PDLASRT( ID, N, D, Q, IQ, JQ, DESCQ, WORK, LWORK,
$ IWORK, LIWORK, INFO )
*
* -- ScaLAPACK auxiliary routine (version 1.7) --
* University of Tennessee, Knoxville, Oak Ridge National Laboratory,
* and University of California, Berkeley.
* February 22, 2000
*
* .. Scalar Arguments ..
CHARACTER ID
INTEGER INFO, IQ, JQ, LIWORK, LWORK, N
* ..
* .. Array Arguments ..
INTEGER DESCQ( * ), IWORK( * )
DOUBLE PRECISION D( * ), Q( * ), WORK( * )
* ..
*
* Purpose
* =======
*
* PDLASRT Sort the numbers in D in increasing order and the
* corresponding vectors in Q.
*
* Arguments
* =========
*
* ID (global input) CHARACTER*1
* = 'I': sort D in increasing order;
* = 'D': sort D in decreasing order. (NOT IMPLEMENTED YET)
*
* N (global input) INTEGER
* The number of columns to be operated on i.e the number of
* columns of the distributed submatrix sub( Q ). N >= 0.
*
* D (global input/output) DOUBLE PRECISION array, dimmension (N)
* On exit, the number in D are sorted in increasing order.
*
* Q (local input) DOUBLE PRECISION pointer into the local memory
* to an array of dimension (LLD_Q, LOCc(JQ+N-1) ). This array
* contains the local pieces of the distributed matrix sub( A )
* to be copied from.
*
* IQ (global input) INTEGER
* The row index in the global array A indicating the first
* row of sub( Q ).
*
* JQ (global input) INTEGER
* The column index in the global array A indicating the
* first column of sub( Q ).
*
* DESCQ (global and local input) INTEGER array of dimension DLEN_.
* The array descriptor for the distributed matrix A.
*
* WORK (local workspace/local output) DOUBLE PRECISION array,
* dimension (LWORK)
* LWORK (local or global input) INTEGER
* The dimension of the array WORK.
* LWORK = MAX( N, NP * ( NB + NQ ))
* where
* NP = NUMROC( N, NB, MYROW, IAROW, NPROW ),
* NQ = NUMROC( N, NB, MYCOL, DESCQ( CSRC_ ), NPCOL )
*
* IWORK (local workspace/local output) INTEGER array,
* dimension (LIWORK)
*
* LIWORK (local or global input) INTEGER
* The dimension of the array IWORK.
* LIWORK = N + 2*NB + 2*NPCOL
*
* INFO (global output) INTEGER
* = 0: successful exit
* < 0: If the i-th argument is an array and the j-entry had
* an illegal value, then INFO = -(i*100+j), if the i-th
* argument is a scalar and had an illegal value, then
* INFO = -i.
*
* =====================================================================
*
* .. Parameters ..
INTEGER BLOCK_CYCLIC_2D, DLEN_, DTYPE_, CTXT_, M_, N_,
$ MB_, NB_, RSRC_, CSRC_, LLD_
PARAMETER ( BLOCK_CYCLIC_2D = 1, DLEN_ = 9, DTYPE_ = 1,
$ CTXT_ = 2, M_ = 3, N_ = 4, MB_ = 5, NB_ = 6,
$ RSRC_ = 7, CSRC_ = 8, LLD_ = 9 )
* ..
* .. Local Scalars ..
INTEGER CL, COL, DUMMY, I, ICTXT, IID, IIQ, INDCOL,
$ INDX, INDXC, INDXG, IPQ, IPQ2, IPW, IPWORK, J,
$ JJQ, K, L, LDQ, LEND, LIWMIN, LWMIN, MYCOL,
$ MYROW, NB, ND, NP, NPCOL, NPROW, NQ, PSQ, QCOL,
$ QTOT, SBUF
* ..
* .. External Functions ..
LOGICAL LSAME
INTEGER INDXG2L, INDXG2P, NUMROC
EXTERNAL INDXG2L, INDXG2P, LSAME, NUMROC
* ..
* .. External Subroutines ..
EXTERNAL BLACS_GRIDINFO, CHK1MAT, PXERBLA, DCOPY,
$ DGERV2D, DGESD2D, DLACPY, DLAPST
* ..
* .. Intrinsic Functions ..
INTRINSIC MAX, MIN, MOD
* ..
* .. Executable Statements ..
*
* This is just to keep ftnchek and toolpack/1 happy
IF( BLOCK_CYCLIC_2D*CSRC_*CTXT_*DLEN_*DTYPE_*LLD_*MB_*M_*NB_*N_*
$ RSRC_.LT.0 )RETURN
*
IF( N.EQ.0 )
$ RETURN
*
CALL BLACS_GRIDINFO( DESCQ( CTXT_ ), NPROW, NPCOL, MYROW, MYCOL )
*
* Test the input parameters
*
INFO = 0
IF( NPROW.EQ.-1 ) THEN
INFO = -( 600+CTXT_ )
ELSE
CALL CHK1MAT( N, 1, N, 1, IQ, JQ, DESCQ, 6, INFO )
IF( INFO.EQ.0 ) THEN
NB = DESCQ( NB_ )
LDQ = DESCQ( LLD_ )
NP = NUMROC( N, NB, MYROW, DESCQ( RSRC_ ), NPROW )
NQ = NUMROC( N, NB, MYCOL, DESCQ( CSRC_ ), NPCOL )
LWMIN = MAX( N, NP*( NB+NQ ) )
LIWMIN = N + 2*( NB+NPCOL )
IF( .NOT.LSAME( ID, 'I' ) ) THEN
INFO = -1
ELSE IF( N.LT.0 ) THEN
INFO = -2
ELSE IF( LWORK.LT.LWMIN ) THEN
INFO = -9
ELSE IF( LIWORK.LT.LIWMIN ) THEN
INFO = -11
END IF
END IF
END IF
*
IF( INFO.NE.0 ) THEN
CALL PXERBLA( ICTXT, 'PDLASRT', -INFO )
RETURN
END IF
*
* Set Pointers
*
INDXC = 1
INDX = INDXC + N
INDXG = INDX
INDCOL = INDXG + NB
QTOT = INDCOL + NB
PSQ = QTOT + NPCOL
*
IID = 1
IPQ2 = 1
IPW = IPQ2 + NP*NQ
*
DUMMY = 0
IIQ = INDXG2L( IQ, NB, DUMMY, DUMMY, NPROW )
*
* Sort the eigenvalues in D
*
CALL DLAPST( 'I', N, D, IWORK( INDX ), INFO )
*
DO 10 L = 0, N - 1
WORK( IID+L ) = D( IWORK( INDX+L ) )
IWORK( INDXC-1+IWORK( INDX+L ) ) = IID + L
10 CONTINUE
CALL DCOPY( N, WORK, 1, D, 1 )
*
ND = 0
20 CONTINUE
IF( ND.LT.N ) THEN
LEND = MIN( NB, N-ND )
J = JQ + ND
QCOL = INDXG2P( J, NB, DUMMY, DESCQ( CSRC_ ), NPCOL )
K = 0
DO 30 L = 0, LEND - 1
I = JQ - 1 + IWORK( INDXC+ND+L )
CL = INDXG2P( I, NB, DUMMY, DESCQ( CSRC_ ), NPCOL )
IWORK( INDCOL+L ) = CL
IF( MYCOL.EQ.CL ) THEN
IWORK( INDXG+K ) = IWORK( INDXC+ND+L )
K = K + 1
END IF
30 CONTINUE
*
IF( MYCOL.EQ.QCOL ) THEN
DO 40 CL = 0, NPCOL - 1
IWORK( QTOT+CL ) = 0
40 CONTINUE
DO 50 L = 0, LEND - 1
IWORK( QTOT+IWORK( INDCOL+L ) ) = IWORK( QTOT+
$ IWORK( INDCOL+L ) ) + 1
50 CONTINUE
IWORK( PSQ ) = 1
DO 60 CL = 1, NPCOL - 1
IWORK( PSQ+CL ) = IWORK( PSQ+CL-1 ) + IWORK( QTOT+CL-1 )
60 CONTINUE
DO 70 L = 0, LEND - 1
CL = IWORK( INDCOL+L )
I = JQ + ND + L
JJQ = INDXG2L( I, NB, DUMMY, DUMMY, NPCOL )
IPQ = IIQ + ( JJQ-1 )*LDQ
IPWORK = IPW + ( IWORK( PSQ+CL )-1 )*NP
CALL DCOPY( NP, Q( IPQ ), 1, WORK( IPWORK ), 1 )
IWORK( PSQ+CL ) = IWORK( PSQ+CL ) + 1
70 CONTINUE
IWORK( PSQ ) = 1
DO 80 CL = 1, NPCOL - 1
IWORK( PSQ+CL ) = IWORK( PSQ+CL-1 ) + IWORK( QTOT+CL-1 )
80 CONTINUE
DO 90 L = 0, K - 1
I = IWORK( INDXG+L )
JJQ = INDXG2L( I, NB, DUMMY, DUMMY, NPCOL )
IPQ = IPQ2 + ( JJQ-1 )*NP
IPWORK = IPW + ( IWORK( PSQ+MYCOL )-1 )*NP
CALL DCOPY( NP, WORK( IPWORK ), 1, WORK( IPQ ), 1 )
IWORK( PSQ+MYCOL ) = IWORK( PSQ+MYCOL ) + 1
90 CONTINUE
DO 100 CL = 1, NPCOL - 1
COL = MOD( MYCOL+CL, NPCOL )
SBUF = IWORK( QTOT+COL )
IF( SBUF.NE.0 ) THEN
IPWORK = IPW + ( IWORK( PSQ+COL )-1 )*NP
CALL DGESD2D( DESCQ( CTXT_ ), NP, SBUF,
$ WORK( IPWORK ), NP, MYROW, COL )
END IF
100 CONTINUE
*
ELSE
*
IF( K.NE.0 ) THEN
CALL DGERV2D( DESCQ( CTXT_ ), NP, K, WORK( IPW ), NP,
$ MYROW, QCOL )
DO 110 L = 0, K - 1
I = JQ - 1 + IWORK( INDXG+L )
JJQ = INDXG2L( I, NB, DUMMY, DUMMY, NPCOL )
IPQ = 1 + ( JJQ-1 )*NP
IPWORK = IPW + L*NP
CALL DCOPY( NP, WORK( IPWORK ), 1, WORK( IPQ ), 1 )
110 CONTINUE
END IF
END IF
ND = ND + NB
GO TO 20
END IF
CALL DLACPY( 'Full', NP, NQ, WORK, NP, Q( IIQ ), LDQ )
*
* End of PDLASRT
*
END
|