Actual source code: sortd.c
petsc-3.8.0 2017-09-26
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>
10: #define SWAP(a,b,t) {t=a;a=b;b=t;}
12: /* A simple version of quicksort; taken from Kernighan and Ritchie, page 87 */
13: static PetscErrorCode PetscSortReal_Private(PetscReal *v,PetscInt right)
14: {
15: PetscInt i,last;
16: PetscReal vl,tmp;
19: if (right <= 1) {
20: if (right == 1) {
21: if (v[0] > v[1]) SWAP(v[0],v[1],tmp);
22: }
23: return(0);
24: }
25: SWAP(v[0],v[right/2],tmp);
26: vl = v[0];
27: last = 0;
28: for (i=1; i<=right; i++) {
29: if (v[i] < vl) {last++; SWAP(v[last],v[i],tmp);}
30: }
31: SWAP(v[0],v[last],tmp);
32: PetscSortReal_Private(v,last-1);
33: PetscSortReal_Private(v+last+1,right-(last+1));
34: return(0);
35: }
37: /*@
38: PetscSortReal - Sorts an array of doubles in place in increasing order.
40: Not Collective
42: Input Parameters:
43: + n - number of values
44: - v - array of doubles
46: Level: intermediate
48: Concepts: sorting^doubles
50: .seealso: PetscSortInt(), PetscSortRealWithPermutation()
51: @*/
52: PetscErrorCode PetscSortReal(PetscInt n,PetscReal v[])
53: {
54: PetscInt j,k;
55: PetscReal tmp,vk;
58: if (n<8) {
59: for (k=0; k<n; k++) {
60: vk = v[k];
61: for (j=k+1; j<n; j++) {
62: if (vk > v[j]) {
63: SWAP(v[k],v[j],tmp);
64: vk = v[k];
65: }
66: }
67: }
68: } else PetscSortReal_Private(v,n-1);
69: return(0);
70: }
73: /*@
74: PetscSortRemoveDupsReal - Sorts an array of doubles in place in increasing order removes all duplicate entries
76: Not Collective
78: Input Parameters:
79: + n - number of values
80: - v - array of doubles
82: Output Parameter:
83: . n - number of non-redundant values
85: Level: intermediate
87: Concepts: sorting^doubles
89: .seealso: PetscSortReal(), PetscSortRemoveDupsInt()
90: @*/
91: PetscErrorCode PetscSortRemoveDupsReal(PetscInt *n,PetscReal v[])
92: {
94: PetscInt i,s = 0,N = *n, b = 0;
97: PetscSortReal(N,v);
98: for (i=0; i<N-1; i++) {
99: if (v[b+s+1] != v[b]) {
100: v[b+1] = v[b+s+1]; b++;
101: } else s++;
102: }
103: *n = N - s;
104: return(0);
105: }
107: /*@
108: PetscSortSplit - Quick-sort split of an array of PetscScalars in place.
110: Not Collective
112: Input Parameters:
113: + ncut - splitig index
114: . n - number of values to sort
115: . a - array of values
116: - idx - index for array a
118: Output Parameters:
119: + a - permuted array of values such that its elements satisfy:
120: abs(a[i]) >= abs(a[ncut-1]) for i < ncut and
121: abs(a[i]) <= abs(a[ncut-1]) for i >= ncut
122: - idx - permuted index of array a
124: Level: intermediate
126: Concepts: sorting^doubles
128: .seealso: PetscSortInt(), PetscSortRealWithPermutation()
129: @*/
130: PetscErrorCode PetscSortSplit(PetscInt ncut,PetscInt n,PetscScalar a[],PetscInt idx[])
131: {
132: PetscInt i,mid,last,itmp,j,first;
133: PetscScalar d,tmp;
134: PetscReal abskey;
137: first = 0;
138: last = n-1;
139: if (ncut < first || ncut > last) return(0);
141: while (1) {
142: mid = first;
143: d = a[mid];
144: abskey = PetscAbsScalar(d);
145: i = last;
146: for (j = first + 1; j <= i; ++j) {
147: d = a[j];
148: if (PetscAbsScalar(d) >= abskey) {
149: ++mid;
150: /* interchange */
151: tmp = a[mid]; itmp = idx[mid];
152: a[mid] = a[j]; idx[mid] = idx[j];
153: a[j] = tmp; idx[j] = itmp;
154: }
155: }
157: /* interchange */
158: tmp = a[mid]; itmp = idx[mid];
159: a[mid] = a[first]; idx[mid] = idx[first];
160: a[first] = tmp; idx[first] = itmp;
162: /* test for while loop */
163: if (mid == ncut) break;
164: else if (mid > ncut) last = mid - 1;
165: else first = mid + 1;
166: }
167: return(0);
168: }
170: /*@
171: PetscSortSplitReal - Quick-sort split of an array of PetscReals in place.
173: Not Collective
175: Input Parameters:
176: + ncut - splitig index
177: . n - number of values to sort
178: . a - array of values in PetscReal
179: - idx - index for array a
181: Output Parameters:
182: + a - permuted array of real values such that its elements satisfy:
183: abs(a[i]) >= abs(a[ncut-1]) for i < ncut and
184: abs(a[i]) <= abs(a[ncut-1]) for i >= ncut
185: - idx - permuted index of array a
187: Level: intermediate
189: Concepts: sorting^doubles
191: .seealso: PetscSortInt(), PetscSortRealWithPermutation()
192: @*/
193: PetscErrorCode PetscSortSplitReal(PetscInt ncut,PetscInt n,PetscReal a[],PetscInt idx[])
194: {
195: PetscInt i,mid,last,itmp,j,first;
196: PetscReal d,tmp;
197: PetscReal abskey;
200: first = 0;
201: last = n-1;
202: if (ncut < first || ncut > last) return(0);
204: while (1) {
205: mid = first;
206: d = a[mid];
207: abskey = PetscAbsReal(d);
208: i = last;
209: for (j = first + 1; j <= i; ++j) {
210: d = a[j];
211: if (PetscAbsReal(d) >= abskey) {
212: ++mid;
213: /* interchange */
214: tmp = a[mid]; itmp = idx[mid];
215: a[mid] = a[j]; idx[mid] = idx[j];
216: a[j] = tmp; idx[j] = itmp;
217: }
218: }
220: /* interchange */
221: tmp = a[mid]; itmp = idx[mid];
222: a[mid] = a[first]; idx[mid] = idx[first];
223: a[first] = tmp; idx[first] = itmp;
225: /* test for while loop */
226: if (mid == ncut) break;
227: else if (mid > ncut) last = mid - 1;
228: else first = mid + 1;
229: }
230: return(0);
231: }