Actual source code: sortd.c

petsc-3.5.4 2015-05-23
Report Typos and Errors
```  2: /*
3:    This file contains routines for sorting doubles.  Values are sorted in place.
4:    These are provided because the general sort routines incur a great deal
5:    of overhead in calling the comparision routines.

7:  */
8: #include <petscsys.h>                /*I  "petscsys.h"  I*/

10: #define SWAP(a,b,t) {t=a;a=b;b=t;}

14: /* A simple version of quicksort; taken from Kernighan and Ritchie, page 87 */
15: static PetscErrorCode PetscSortReal_Private(PetscReal *v,PetscInt right)
16: {
17:   PetscInt  i,last;
18:   PetscReal vl,tmp;

21:   if (right <= 1) {
22:     if (right == 1) {
23:       if (v[0] > v[1]) SWAP(v[0],v[1],tmp);
24:     }
25:     return(0);
26:   }
27:   SWAP(v[0],v[right/2],tmp);
28:   vl   = v[0];
29:   last = 0;
30:   for (i=1; i<=right; i++) {
31:     if (v[i] < vl) {last++; SWAP(v[last],v[i],tmp);}
32:   }
33:   SWAP(v[0],v[last],tmp);
34:   PetscSortReal_Private(v,last-1);
35:   PetscSortReal_Private(v+last+1,right-(last+1));
36:   return(0);
37: }

41: /*@
42:    PetscSortReal - Sorts an array of doubles in place in increasing order.

44:    Not Collective

46:    Input Parameters:
47: +  n  - number of values
48: -  v  - array of doubles

50:    Level: intermediate

52:    Concepts: sorting^doubles

54: .seealso: PetscSortInt(), PetscSortRealWithPermutation()
55: @*/
56: PetscErrorCode  PetscSortReal(PetscInt n,PetscReal v[])
57: {
58:   PetscInt  j,k;
59:   PetscReal tmp,vk;

62:   if (n<8) {
63:     for (k=0; k<n; k++) {
64:       vk = v[k];
65:       for (j=k+1; j<n; j++) {
66:         if (vk > v[j]) {
67:           SWAP(v[k],v[j],tmp);
68:           vk = v[k];
69:         }
70:       }
71:     }
72:   } else PetscSortReal_Private(v,n-1);
73:   return(0);
74: }

79: /*@
80:    PetscSortRemoveDupsReal - Sorts an array of doubles in place in increasing order removes all duplicate entries

82:    Not Collective

84:    Input Parameters:
85: +  n  - number of values
86: -  v  - array of doubles

88:    Output Parameter:
89: .  n - number of non-redundant values

91:    Level: intermediate

93:    Concepts: sorting^doubles

95: .seealso: PetscSortReal(), PetscSortRemoveDupsInt()
96: @*/
97: PetscErrorCode  PetscSortRemoveDupsReal(PetscInt *n,PetscReal v[])
98: {
100:   PetscInt       i,s = 0,N = *n, b = 0;

103:   PetscSortReal(N,v);
104:   for (i=0; i<N-1; i++) {
105:     if (v[b+s+1] != v[b]) {
106:       v[b+1] = v[b+s+1]; b++;
107:     } else s++;
108:   }
109:   *n = N - s;
110:   return(0);
111: }

115: /*@
116:    PetscSortSplit - Quick-sort split of an array of PetscScalars in place.

118:    Not Collective

120:    Input Parameters:
121: +  ncut  - splitig index
122: .  n     - number of values to sort
123: .  a     - array of values
124: -  idx   - index for array a

126:    Output Parameters:
127: +  a     - permuted array of values such that its elements satisfy:
128:            abs(a[i]) >= abs(a[ncut-1]) for i < ncut and
129:            abs(a[i]) <= abs(a[ncut-1]) for i >= ncut
130: -  idx   - permuted index of array a

132:    Level: intermediate

134:    Concepts: sorting^doubles

136: .seealso: PetscSortInt(), PetscSortRealWithPermutation()
137: @*/
138: PetscErrorCode  PetscSortSplit(PetscInt ncut,PetscInt n,PetscScalar a[],PetscInt idx[])
139: {
140:   PetscInt    i,mid,last,itmp,j,first;
141:   PetscScalar d,tmp;
142:   PetscReal   abskey;

145:   first = 0;
146:   last  = n-1;
147:   if (ncut < first || ncut > last) return(0);

149:   while (1) {
150:     mid    = first;
151:     abskey = (d = a[mid],PetscAbsScalar(d));
152:     i      = last;
153:     for (j = first + 1; j <= i; ++j) {
154:       if ((d = a[j],PetscAbsScalar(d)) >= abskey) {
155:         ++mid;
156:         /* interchange */
157:         tmp = a[mid];  itmp = idx[mid];
158:         a[mid] = a[j]; idx[mid] = idx[j];
159:         a[j] = tmp;    idx[j] = itmp;
160:       }
161:     }

163:     /* interchange */
164:     tmp = a[mid];      itmp = idx[mid];
165:     a[mid] = a[first]; idx[mid] = idx[first];
166:     a[first] = tmp;    idx[first] = itmp;

168:     /* test for while loop */
169:     if (mid == ncut) break;
170:     else if (mid > ncut) last = mid - 1;
171:     else first = mid + 1;
172:   }
173:   return(0);
174: }

178: /*@
179:    PetscSortSplitReal - Quick-sort split of an array of PetscReals in place.

181:    Not Collective

183:    Input Parameters:
184: +  ncut  - splitig index
185: .  n     - number of values to sort
186: .  a     - array of values in PetscReal
187: -  idx   - index for array a

189:    Output Parameters:
190: +  a     - permuted array of real values such that its elements satisfy:
191:            abs(a[i]) >= abs(a[ncut-1]) for i < ncut and
192:            abs(a[i]) <= abs(a[ncut-1]) for i >= ncut
193: -  idx   - permuted index of array a

195:    Level: intermediate

197:    Concepts: sorting^doubles

199: .seealso: PetscSortInt(), PetscSortRealWithPermutation()
200: @*/
201: PetscErrorCode  PetscSortSplitReal(PetscInt ncut,PetscInt n,PetscReal a[],PetscInt idx[])
202: {
203:   PetscInt  i,mid,last,itmp,j,first;
204:   PetscReal d,tmp;
205:   PetscReal abskey;

208:   first = 0;
209:   last  = n-1;
210:   if (ncut < first || ncut > last) return(0);

212:   while (1) {
213:     mid    = first;
214:     abskey = (d = a[mid],PetscAbsReal(d));
215:     i      = last;
216:     for (j = first + 1; j <= i; ++j) {
217:       if ((d = a[j],PetscAbsReal(d)) >= abskey) {
218:         ++mid;
219:         /* interchange */
220:         tmp = a[mid];  itmp = idx[mid];
221:         a[mid] = a[j]; idx[mid] = idx[j];
222:         a[j] = tmp;    idx[j] = itmp;
223:       }
224:     }

226:     /* interchange */
227:     tmp = a[mid];      itmp = idx[mid];
228:     a[mid] = a[first]; idx[mid] = idx[first];
229:     a[first] = tmp;    idx[first] = itmp;

231:     /* test for while loop */
232:     if (mid == ncut) break;
233:     else if (mid > ncut) last = mid - 1;
234:     else first = mid + 1;
235:   }
236:   return(0);
237: }

```