~fluidity-core/fluidity/embedded_models

« back to all changes in this revision

Viewing changes to libmba2d/dsort.f

  • Committer: Timothy Bond
  • Date: 2011-04-14 15:40:24 UTC
  • Revision ID: timothy.bond@imperial.ac.uk-20110414154024-116ci9gq6mwigmaw
Following the move from svn to bzr we change the nature of inclusion of these
four software libraries. Previously, they were included as svn externals and
pulled at latest version for trunk, pinned to specific versions for release
and stable trunk. Since bzr is less elegant at dealing with externals we have
made the decision to include the packages directly into the trunk instead.

At this import the versions are:

libadaptivity: r163
libvtkfortran: r67
libspud: r545
libmba2d: r28

Show diffs side-by-side

added added

removed removed

Lines of Context:
 
1
      SUBROUTINE DSORT (DX, IX, N, KFLAG, IERR)
 
2
C***BEGIN PROLOGUE  DSORT
 
3
C***PURPOSE  Sort an array and optionally make the same interchanges in
 
4
C            an auxiliary array.  The array may be sorted in increasing
 
5
C            or decreasing order.  A slightly modified QUICKSORT
 
6
C            algorithm is used.
 
7
C***LIBRARY   SLATEC
 
8
C***CATEGORY  N6A2B
 
9
C***TYPE      real (SSORT-S, DSORT-D, ISORT-I)
 
10
C***KEYWORDS  SINGLETON QUICKSORT, SORT, SORTING
 
11
C***AUTHOR  Jones, R. E., (SNLA)
 
12
C           Wisniewski, J. A., (SNLA)
 
13
C***DESCRIPTION
 
14
C
 
15
C   DSORT sorts array DX and optionally makes the same interchanges in
 
16
C   array IX.  The array DX may be sorted in increasing order or
 
17
C   decreasing order.  A slightly modified quicksort algorithm is used.
 
18
C
 
19
C   Description of Parameters
 
20
C      DX - array of values to be sorted   (usually abscissas)
 
21
C      IX - array to be (optionally) carried along
 
22
C      N  - number of values in array DX to be sorted
 
23
C      KFLAG - control parameter
 
24
C            =  2  means sort DX in increasing order and carry IX along.
 
25
C            =  1  means sort DX in increasing order (ignoring IX)
 
26
C            = -1  means sort DX in decreasing order (ignoring IX)
 
27
C            = -2  means sort DX in decreasing order and carry IX along.
 
28
C
 
29
C      IERR - the error's indicator
 
30
C           = 0 - Okay
 
31
C           != 0 - errors
 
32
C
 
33
C***REFERENCES  R. C. Singleton, Algorithm 347, An efficient algorithm
 
34
C                 for sorting with minimal storage, Communications of
 
35
C                 the ACM, 12, 3 (1969), pp. 185-187.
 
36
C***ROUTINES CALLED  XERMSG
 
37
C***REVISION HISTORY  (YYMMDD)
 
38
C   761101  DATE WRITTEN
 
39
C   761118  Modified to use the Singleton quicksort algorithm.  (JAW)
 
40
C   890531  Changed all specific intrinsics to generic.  (WRB)
 
41
C   890831  Modified array declarations.  (WRB)
 
42
C   891009  Removed unreferenced statement labels.  (WRB)
 
43
C   891024  Changed category.  (WRB)
 
44
C   891024  REVISION DATE from Version 3.2
 
45
C   891214  Prologue converted to Version 4.0 format.  (BAB)
 
46
C   900315  CALLs to XERROR changed to CALLs to XERMSG.  (THJ)
 
47
C   901012  Declared all variables; changed X,Y to DX,IX; changed
 
48
C           code to parallel SSORT. (M. McClain)
 
49
C   920501  Reformatted the REFERENCES section.  (DWL, WRB)
 
50
C   920519  Clarified error messages.  (DWL)
 
51
C   920801  Declarations section rebuilt and code restructured to use
 
52
C           IF-THEN-ELSE-ENDIF.  (RWC, WRB)
 
53
C***END PROLOGUE  DSORT
 
54
C     .. Scalar Arguments ..
 
55
      INTEGER KFLAG, N
 
56
C     .. Array Arguments ..
 
57
      real DX(*)
 
58
      INTEGER IX(*)
 
59
C     .. Local Scalars ..
 
60
      real R, T, TT
 
61
      INTEGER I, IJ, J, K, KK, L, M, NN
 
62
C     .. Local Arrays ..
 
63
      INTEGER IL(21), IU(21)
 
64
C     .. Intrinsic Functions ..
 
65
      INTRINSIC ABS, INT
 
66
C***FIRST EXECUTABLE STATEMENT  DSORT
 
67
      IERR = 0
 
68
 
 
69
      NN = N
 
70
      IF (NN .LT. 1) THEN
 
71
         IERR = 1201
 
72
         RETURN
 
73
      ENDIF
 
74
C
 
75
      KK = ABS(KFLAG)
 
76
      IF (KK.NE.1 .AND. KK.NE.2) THEN
 
77
         IERR = 1202
 
78
         RETURN
 
79
      ENDIF
 
80
C
 
81
C     Initializing IX
 
82
C
 
83
C      DO 5 I=1,NN
 
84
C        IX(I)=I
 
85
C 5    CONTINUE
 
86
 
 
87
C
 
88
C     Alter array DX to get decreasing order if needed
 
89
C
 
90
      IF (KFLAG .LE. -1) THEN
 
91
         DO 10 I=1,NN
 
92
            DX(I) = -DX(I)
 
93
   10    CONTINUE
 
94
      ENDIF
 
95
C
 
96
      IF (KK .EQ. 2) GO TO 100
 
97
C
 
98
C     Sort DX only
 
99
C
 
100
      M = 1
 
101
      I = 1
 
102
      J = NN
 
103
      R = 0.375D0
 
104
C
 
105
   20 IF (I .EQ. J) GO TO 60
 
106
      IF (R .LE. 0.5898437D0) THEN
 
107
         R = R+3.90625D-2
 
108
      ELSE
 
109
         R = R-0.21875D0
 
110
      ENDIF
 
111
C
 
112
   30 K = I
 
113
C
 
114
C     Select a central element of the array and save it in location T
 
115
C
 
116
      IJ = I + INT((J-I)*R)
 
117
      T = DX(IJ)
 
118
C
 
119
C     If first element of array is greater than T, interchange with T
 
120
C
 
121
      IF (DX(I) .GT. T) THEN
 
122
         DX(IJ) = DX(I)
 
123
         DX(I) = T
 
124
         T = DX(IJ)
 
125
      ENDIF
 
126
      L = J
 
127
C
 
128
C     If last element of array is less than than T, interchange with T
 
129
C
 
130
      IF (DX(J) .LT. T) THEN
 
131
         DX(IJ) = DX(J)
 
132
         DX(J) = T
 
133
         T = DX(IJ)
 
134
C
 
135
C        If first element of array is greater than T, interchange with T
 
136
C
 
137
         IF (DX(I) .GT. T) THEN
 
138
            DX(IJ) = DX(I)
 
139
            DX(I) = T
 
140
            T = DX(IJ)
 
141
         ENDIF
 
142
      ENDIF
 
143
C
 
144
C     Find an element in the second half of the array which is smaller
 
145
C     than T
 
146
C
 
147
   40 L = L-1
 
148
      IF (DX(L) .GT. T) GO TO 40
 
149
C
 
150
C     Find an element in the first half of the array which is greater
 
151
C     than T
 
152
C
 
153
   50 K = K+1
 
154
      IF (DX(K) .LT. T) GO TO 50
 
155
C
 
156
C     Interchange these elements
 
157
C
 
158
      IF (K .LE. L) THEN
 
159
         TT = DX(L)
 
160
         DX(L) = DX(K)
 
161
         DX(K) = TT
 
162
         GO TO 40
 
163
      ENDIF
 
164
C
 
165
C     Save upper and lower subscripts of the array yet to be sorted
 
166
C
 
167
      IF (L-I .GT. J-K) THEN
 
168
         IL(M) = I
 
169
         IU(M) = L
 
170
         I = K
 
171
         M = M+1
 
172
      ELSE
 
173
         IL(M) = K
 
174
         IU(M) = J
 
175
         J = L
 
176
         M = M+1
 
177
      ENDIF
 
178
      GO TO 70
 
179
C
 
180
C     Begin again on another portion of the unsorted array
 
181
C
 
182
   60 M = M-1
 
183
      IF (M .EQ. 0) GO TO 190
 
184
      I = IL(M)
 
185
      J = IU(M)
 
186
C
 
187
   70 IF (J-I .GE. 1) GO TO 30
 
188
      IF (I .EQ. 1) GO TO 20
 
189
      I = I-1
 
190
C
 
191
   80 I = I+1
 
192
      IF (I .EQ. J) GO TO 60
 
193
      T = DX(I+1)
 
194
      IF (DX(I) .LE. T) GO TO 80
 
195
      K = I
 
196
C
 
197
   90 DX(K+1) = DX(K)
 
198
      K = K-1
 
199
      IF (T .LT. DX(K)) GO TO 90
 
200
      DX(K+1) = T
 
201
      GO TO 80
 
202
C
 
203
C     Sort DX and carry IX along
 
204
C
 
205
  100 M = 1
 
206
      I = 1
 
207
      J = NN
 
208
      R = 0.375D0
 
209
C
 
210
  110 IF (I .EQ. J) GO TO 150
 
211
      IF (R .LE. 0.5898437D0) THEN
 
212
         R = R+3.90625D-2
 
213
      ELSE
 
214
         R = R-0.21875D0
 
215
      ENDIF
 
216
C
 
217
  120 K = I
 
218
C
 
219
C     Select a central element of the array and save it in location T
 
220
C
 
221
      IJ = I + INT((J-I)*R)
 
222
      T = DX(IJ)
 
223
      ITX = IX(IJ)
 
224
C
 
225
C     If first element of array is greater than T, interchange with T
 
226
C
 
227
      IF (DX(I) .GT. T) THEN
 
228
         DX(IJ) = DX(I)
 
229
         DX(I) = T
 
230
         T = DX(IJ)
 
231
         IX(IJ) = IX(I)
 
232
         IX(I) = ITX
 
233
         ITX = IX(IJ)
 
234
      ENDIF
 
235
      L = J
 
236
C
 
237
C     If last element of array is less than T, interchange with T
 
238
C
 
239
      IF (DX(J) .LT. T) THEN
 
240
         DX(IJ) = DX(J)
 
241
         DX(J) = T
 
242
         T = DX(IJ)
 
243
         IX(IJ) = IX(J)
 
244
         IX(J) = ITX
 
245
         ITX = IX(IJ)
 
246
C
 
247
C        If first element of array is greater than T, interchange with T
 
248
C
 
249
         IF (DX(I) .GT. T) THEN
 
250
            DX(IJ) = DX(I)
 
251
            DX(I) = T
 
252
            T = DX(IJ)
 
253
            IX(IJ) = IX(I)
 
254
            IX(I) = ITX
 
255
            ITX = IX(IJ)
 
256
         ENDIF
 
257
      ENDIF
 
258
C
 
259
C     Find an element in the second half of the array which is smaller
 
260
C     than T
 
261
C
 
262
  130 L = L-1
 
263
      IF (DX(L) .GT. T) GO TO 130
 
264
C
 
265
C     Find an element in the first half of the array which is greater
 
266
C     than T
 
267
C
 
268
  140 K = K+1
 
269
      IF (DX(K) .LT. T) GO TO 140
 
270
C
 
271
C     Interchange these elements
 
272
C
 
273
      IF (K .LE. L) THEN
 
274
         TT = DX(L)
 
275
         DX(L) = DX(K)
 
276
         DX(K) = TT
 
277
         ITTX = IX(L)
 
278
         IX(L) = IX(K)
 
279
         IX(K) = ITTX
 
280
         GO TO 130
 
281
      ENDIF
 
282
C
 
283
C     Save upper and lower subscripts of the array yet to be sorted
 
284
C
 
285
      IF (L-I .GT. J-K) THEN
 
286
         IL(M) = I
 
287
         IU(M) = L
 
288
         I = K
 
289
         M = M+1
 
290
      ELSE
 
291
         IL(M) = K
 
292
         IU(M) = J
 
293
         J = L
 
294
         M = M+1
 
295
      ENDIF
 
296
      GO TO 160
 
297
C
 
298
C     Begin again on another portion of the unsorted array
 
299
C
 
300
  150 M = M-1
 
301
      IF (M .EQ. 0) GO TO 190
 
302
      I = IL(M)
 
303
      J = IU(M)
 
304
C
 
305
  160 IF (J-I .GE. 1) GO TO 120
 
306
      IF (I .EQ. 1) GO TO 110
 
307
      I = I-1
 
308
C
 
309
  170 I = I+1
 
310
      IF (I .EQ. J) GO TO 150
 
311
      T = DX(I+1)
 
312
      ITX = IX(I+1)
 
313
      IF (DX(I) .LE. T) GO TO 170
 
314
      K = I
 
315
C
 
316
  180 DX(K+1) = DX(K)
 
317
      IX(K+1) = IX(K)
 
318
      K = K-1
 
319
      IF (T .LT. DX(K)) GO TO 180
 
320
      DX(K+1) = T
 
321
      IX(K+1) = ITX
 
322
      GO TO 170
 
323
C
 
324
C     Clean up
 
325
C
 
326
  190 IF (KFLAG .LE. -1) THEN
 
327
         DO 200 I=1,NN
 
328
            DX(I) = -DX(I)
 
329
  200    CONTINUE
 
330
      ENDIF
 
331
      RETURN
 
332
      END
 
333
 
 
334
 
 
335