1 SUBROUTINE DLAPST( ID, N, D, INDX, INFO ) 2* 3* -- ScaLAPACK auxiliary routine (version 1.7) -- 4* University of Tennessee, Knoxville, Oak Ridge National Laboratory, 5* and University of California, Berkeley. 6* December 31, 1998 7* 8* .. Scalar Arguments .. 9 CHARACTER ID 10 INTEGER INFO, N 11* .. 12* .. Array Arguments .. 13 INTEGER INDX( * ) 14 DOUBLE PRECISION D( * ) 15* .. 16* 17* Purpose 18* ======= 19* DLAPST is a modified version of the LAPACK routine DLASRT. 20* 21* Define a permutation INDX that sorts the numbers in D 22* in increasing order (if ID = 'I') or 23* in decreasing order (if ID = 'D' ). 24* 25* Use Quick Sort, reverting to Insertion sort on arrays of 26* size <= 20. Dimension of STACK limits N to about 2**32. 27* 28* Arguments 29* ========= 30* 31* ID (input) CHARACTER*1 32* = 'I': sort D in increasing order; 33* = 'D': sort D in decreasing order. 34* 35* N (input) INTEGER 36* The length of the array D. 37* 38* D (input) DOUBLE PRECISION array, dimension (N) 39* The array to be sorted. 40* 41* INDX (ouput) INTEGER array, dimension (N). 42* The permutation which sorts the array D. 43* 44* INFO (output) INTEGER 45* = 0: successful exit 46* < 0: if INFO = -i, the i-th argument had an illegal value 47* 48* ===================================================================== 49* 50* .. Parameters .. 51 INTEGER SELECT 52 PARAMETER ( SELECT = 20 ) 53* .. 54* .. Local Scalars .. 55 INTEGER DIR, ENDD, I, ITMP, J, START, STKPNT 56 DOUBLE PRECISION D1, D2, D3, DMNMX 57* .. 58* .. Local Arrays .. 59 INTEGER STACK( 2, 32 ) 60* .. 61* .. External Functions .. 62 LOGICAL LSAME 63 EXTERNAL LSAME 64* .. 65* .. External Subroutines .. 66 EXTERNAL XERBLA 67* .. 68* .. Executable Statements .. 69* 70* Test the input paramters. 71* 72 INFO = 0 73 DIR = -1 74 IF( LSAME( ID, 'D' ) ) THEN 75 DIR = 0 76 ELSE IF( LSAME( ID, 'I' ) ) THEN 77 DIR = 1 78 END IF 79 IF( DIR.EQ.-1 ) THEN 80 INFO = -1 81 ELSE IF( N.LT.0 ) THEN 82 INFO = -2 83 END IF 84 IF( INFO.NE.0 ) THEN 85 CALL XERBLA( 'DLAPST', -INFO ) 86 RETURN 87 END IF 88* 89* Quick return if possible 90* 91 IF( N.LE.1 ) 92 $ RETURN 93* 94 DO 10 I = 1, N 95 INDX( I ) = I 96 10 CONTINUE 97* 98 STKPNT = 1 99 STACK( 1, 1 ) = 1 100 STACK( 2, 1 ) = N 101 20 CONTINUE 102 START = STACK( 1, STKPNT ) 103 ENDD = STACK( 2, STKPNT ) 104 STKPNT = STKPNT - 1 105 IF( ENDD-START.LE.SELECT .AND. ENDD-START.GT.0 ) THEN 106* 107* Do Insertion sort on D( START:ENDD ) 108* 109 IF( DIR.EQ.0 ) THEN 110* 111* Sort into decreasing order 112* 113 DO 40 I = START + 1, ENDD 114 DO 30 J = I, START + 1, -1 115 IF( D( INDX( J ) ).GT.D( INDX( J-1 ) ) ) THEN 116 ITMP = INDX( J ) 117 INDX( J ) = INDX( J-1 ) 118 INDX( J-1 ) = ITMP 119 ELSE 120 GO TO 40 121 END IF 122 30 CONTINUE 123 40 CONTINUE 124* 125 ELSE 126* 127* Sort into increasing order 128* 129 DO 60 I = START + 1, ENDD 130 DO 50 J = I, START + 1, -1 131 IF( D( INDX( J ) ).LT.D( INDX( J-1 ) ) ) THEN 132 ITMP = INDX( J ) 133 INDX( J ) = INDX( J-1 ) 134 INDX( J-1 ) = ITMP 135 ELSE 136 GO TO 60 137 END IF 138 50 CONTINUE 139 60 CONTINUE 140* 141 END IF 142* 143 ELSE IF( ENDD-START.GT.SELECT ) THEN 144* 145* Partition D( START:ENDD ) and stack parts, largest one first 146* 147* Choose partition entry as median of 3 148* 149 D1 = D( INDX( START ) ) 150 D2 = D( INDX( ENDD ) ) 151 I = ( START+ENDD ) / 2 152 D3 = D( INDX( I ) ) 153 IF( D1.LT.D2 ) THEN 154 IF( D3.LT.D1 ) THEN 155 DMNMX = D1 156 ELSE IF( D3.LT.D2 ) THEN 157 DMNMX = D3 158 ELSE 159 DMNMX = D2 160 END IF 161 ELSE 162 IF( D3.LT.D2 ) THEN 163 DMNMX = D2 164 ELSE IF( D3.LT.D1 ) THEN 165 DMNMX = D3 166 ELSE 167 DMNMX = D1 168 END IF 169 END IF 170* 171 IF( DIR.EQ.0 ) THEN 172* 173* Sort into decreasing order 174* 175 I = START - 1 176 J = ENDD + 1 177 70 CONTINUE 178 80 CONTINUE 179 J = J - 1 180 IF( D( INDX( J ) ).LT.DMNMX ) 181 $ GO TO 80 182 90 CONTINUE 183 I = I + 1 184 IF( D( INDX( I ) ).GT.DMNMX ) 185 $ GO TO 90 186 IF( I.LT.J ) THEN 187 ITMP = INDX( I ) 188 INDX( I ) = INDX( J ) 189 INDX( J ) = ITMP 190 GO TO 70 191 END IF 192 IF( J-START.GT.ENDD-J-1 ) THEN 193 STKPNT = STKPNT + 1 194 STACK( 1, STKPNT ) = START 195 STACK( 2, STKPNT ) = J 196 STKPNT = STKPNT + 1 197 STACK( 1, STKPNT ) = J + 1 198 STACK( 2, STKPNT ) = ENDD 199 ELSE 200 STKPNT = STKPNT + 1 201 STACK( 1, STKPNT ) = J + 1 202 STACK( 2, STKPNT ) = ENDD 203 STKPNT = STKPNT + 1 204 STACK( 1, STKPNT ) = START 205 STACK( 2, STKPNT ) = J 206 END IF 207 ELSE 208* 209* Sort into increasing order 210* 211 I = START - 1 212 J = ENDD + 1 213 100 CONTINUE 214 110 CONTINUE 215 J = J - 1 216 IF( D( INDX( J ) ).GT.DMNMX ) 217 $ GO TO 110 218 120 CONTINUE 219 I = I + 1 220 IF( D( INDX( I ) ).LT.DMNMX ) 221 $ GO TO 120 222 IF( I.LT.J ) THEN 223 ITMP = INDX( I ) 224 INDX( I ) = INDX( J ) 225 INDX( J ) = ITMP 226 GO TO 100 227 END IF 228 IF( J-START.GT.ENDD-J-1 ) THEN 229 STKPNT = STKPNT + 1 230 STACK( 1, STKPNT ) = START 231 STACK( 2, STKPNT ) = J 232 STKPNT = STKPNT + 1 233 STACK( 1, STKPNT ) = J + 1 234 STACK( 2, STKPNT ) = ENDD 235 ELSE 236 STKPNT = STKPNT + 1 237 STACK( 1, STKPNT ) = J + 1 238 STACK( 2, STKPNT ) = ENDD 239 STKPNT = STKPNT + 1 240 STACK( 1, STKPNT ) = START 241 STACK( 2, STKPNT ) = J 242 END IF 243 END IF 244 END IF 245 IF( STKPNT.GT.0 ) 246 $ GO TO 20 247 RETURN 248* 249* End of DLAPST 250* 251 END 252