Actual source code: sorti.c
2: /*
3: This file contains routines for sorting integers. Values are sorted in place.
4: One can use src/sys/tests/ex52.c for benchmarking.
5: */
6: #include <petsc/private/petscimpl.h>
7: #include <petsc/private/hashseti.h>
9: #define MEDIAN3(v,a,b,c) \
10: (v[a]<v[b] \
11: ? (v[b]<v[c] \
12: ? (b) \
13: : (v[a]<v[c] ? (c) : (a))) \
14: : (v[c]<v[b] \
15: ? (b) \
16: : (v[a]<v[c] ? (a) : (c))))
18: #define MEDIAN(v,right) MEDIAN3(v,right/4,right/2,right/4*3)
20: /* Swap one, two or three pairs. Each pair can have its own type */
21: #define SWAP1(a,b,t1) do {t1=a;a=b;b=t1;} while (0)
22: #define SWAP2(a,b,c,d,t1,t2) do {t1=a;a=b;b=t1; t2=c;c=d;d=t2;} while (0)
23: #define SWAP3(a,b,c,d,e,f,t1,t2,t3) do {t1=a;a=b;b=t1; t2=c;c=d;d=t2; t3=e;e=f;f=t3;} while (0)
25: /* Swap a & b, *c & *d. c, d, t2 are pointers to a type of size <siz> */
26: #define SWAP2Data(a,b,c,d,t1,t2,siz) \
27: do { \
28: t1=a; a=b; b=t1; \
29: PetscMemcpy(t2,c,siz); \
30: PetscMemcpy(c,d,siz); \
31: PetscMemcpy(d,t2,siz); \
32: } while (0)
34: /*
35: Partition X[lo,hi] into two parts: X[lo,l) <= pivot; X[r,hi] > pivot
37: Input Parameters:
38: + X - array to partition
39: . pivot - a pivot of X[]
40: . t1 - temp variable for X
41: - lo,hi - lower and upper bound of the array
43: Output Parameters:
44: + l,r - of type PetscInt
46: Notes:
47: The TwoWayPartition2/3 variants also partition other arrays along with X.
48: These arrays can have different types, so they provide their own temp t2,t3
49: */
50: #define TwoWayPartition1(X,pivot,t1,lo,hi,l,r) \
51: do { \
52: l = lo; \
53: r = hi; \
54: while (1) { \
55: while (X[l] < pivot) l++; \
56: while (X[r] > pivot) r--; \
57: if (l >= r) {r++; break;} \
58: SWAP1(X[l],X[r],t1); \
59: l++; \
60: r--; \
61: } \
62: } while (0)
64: /*
65: Partition X[lo,hi] into two parts: X[lo,l) >= pivot; X[r,hi] < pivot
67: Input Parameters:
68: + X - array to partition
69: . pivot - a pivot of X[]
70: . t1 - temp variable for X
71: - lo,hi - lower and upper bound of the array
73: Output Parameters:
74: + l,r - of type PetscInt
76: Notes:
77: The TwoWayPartition2/3 variants also partition other arrays along with X.
78: These arrays can have different types, so they provide their own temp t2,t3
79: */
80: #define TwoWayPartitionReverse1(X,pivot,t1,lo,hi,l,r) \
81: do { \
82: l = lo; \
83: r = hi; \
84: while (1) { \
85: while (X[l] > pivot) l++; \
86: while (X[r] < pivot) r--; \
87: if (l >= r) {r++; break;} \
88: SWAP1(X[l],X[r],t1); \
89: l++; \
90: r--; \
91: } \
92: } while (0)
94: #define TwoWayPartition2(X,Y,pivot,t1,t2,lo,hi,l,r) \
95: do { \
96: l = lo; \
97: r = hi; \
98: while (1) { \
99: while (X[l] < pivot) l++; \
100: while (X[r] > pivot) r--; \
101: if (l >= r) {r++; break;} \
102: SWAP2(X[l],X[r],Y[l],Y[r],t1,t2); \
103: l++; \
104: r--; \
105: } \
106: } while (0)
108: #define TwoWayPartition3(X,Y,Z,pivot,t1,t2,t3,lo,hi,l,r) \
109: do { \
110: l = lo; \
111: r = hi; \
112: while (1) { \
113: while (X[l] < pivot) l++; \
114: while (X[r] > pivot) r--; \
115: if (l >= r) {r++; break;} \
116: SWAP3(X[l],X[r],Y[l],Y[r],Z[l],Z[r],t1,t2,t3); \
117: l++; \
118: r--; \
119: } \
120: } while (0)
122: /* Templates for similar functions used below */
123: #define QuickSort1(FuncName,X,n,pivot,t1) \
124: do { \
125: PetscCount i,j,p,l,r,hi=n-1; \
126: if (n < 8) { \
127: for (i=0; i<n; i++) { \
128: pivot = X[i]; \
129: for (j=i+1; j<n; j++) { \
130: if (pivot > X[j]) { \
131: SWAP1(X[i],X[j],t1); \
132: pivot = X[i]; \
133: } \
134: } \
135: } \
136: } else { \
137: p = MEDIAN(X,hi); \
138: pivot = X[p]; \
139: TwoWayPartition1(X,pivot,t1,0,hi,l,r); \
140: FuncName(l,X); \
141: FuncName(hi-r+1,X+r); \
142: } \
143: } while (0)
145: /* Templates for similar functions used below */
146: #define QuickSortReverse1(FuncName,X,n,pivot,t1) \
147: do { \
148: PetscCount i,j,p,l,r,hi=n-1; \
149: if (n < 8) { \
150: for (i=0; i<n; i++) { \
151: pivot = X[i]; \
152: for (j=i+1; j<n; j++) { \
153: if (pivot < X[j]) { \
154: SWAP1(X[i],X[j],t1); \
155: pivot = X[i]; \
156: } \
157: } \
158: } \
159: } else { \
160: p = MEDIAN(X,hi); \
161: pivot = X[p]; \
162: TwoWayPartitionReverse1(X,pivot,t1,0,hi,l,r); \
163: FuncName(l,X); \
164: FuncName(hi-r+1,X+r); \
165: } \
166: } while (0)
168: #define QuickSort2(FuncName,X,Y,n,pivot,t1,t2) \
169: do { \
170: PetscCount i,j,p,l,r,hi=n-1; \
171: if (n < 8) { \
172: for (i=0; i<n; i++) { \
173: pivot = X[i]; \
174: for (j=i+1; j<n; j++) { \
175: if (pivot > X[j]) { \
176: SWAP2(X[i],X[j],Y[i],Y[j],t1,t2); \
177: pivot = X[i]; \
178: } \
179: } \
180: } \
181: } else { \
182: p = MEDIAN(X,hi); \
183: pivot = X[p]; \
184: TwoWayPartition2(X,Y,pivot,t1,t2,0,hi,l,r); \
185: FuncName(l,X,Y); \
186: FuncName(hi-r+1,X+r,Y+r); \
187: } \
188: } while (0)
190: #define QuickSort3(FuncName,X,Y,Z,n,pivot,t1,t2,t3) \
191: do { \
192: PetscCount i,j,p,l,r,hi=n-1; \
193: if (n < 8) { \
194: for (i=0; i<n; i++) { \
195: pivot = X[i]; \
196: for (j=i+1; j<n; j++) { \
197: if (pivot > X[j]) { \
198: SWAP3(X[i],X[j],Y[i],Y[j],Z[i],Z[j],t1,t2,t3); \
199: pivot = X[i]; \
200: } \
201: } \
202: } \
203: } else { \
204: p = MEDIAN(X,hi); \
205: pivot = X[p]; \
206: TwoWayPartition3(X,Y,Z,pivot,t1,t2,t3,0,hi,l,r); \
207: FuncName(l,X,Y,Z); \
208: FuncName(hi-r+1,X+r,Y+r,Z+r); \
209: } \
210: } while (0)
212: /*@
213: PetscSortedInt - Determines whether the array is sorted.
215: Not Collective
217: Input Parameters:
218: + n - number of values
219: - X - array of integers
221: Output Parameters:
222: . sorted - flag whether the array is sorted
224: Level: intermediate
226: .seealso: PetscSortInt(), PetscSortedMPIInt(), PetscSortedReal()
227: @*/
228: PetscErrorCode PetscSortedInt(PetscInt n,const PetscInt X[],PetscBool *sorted)
229: {
232: PetscSorted(n,X,*sorted);
233: return 0;
234: }
236: /*@
237: PetscSortInt - Sorts an array of integers in place in increasing order.
239: Not Collective
241: Input Parameters:
242: + n - number of values
243: - X - array of integers
245: Notes:
246: This function serves as an alternative to PetscIntSortSemiOrdered(), and may perform faster especially if the array
247: is completely random. There are exceptions to this and so it is __highly__ recommended that the user benchmark their
248: code to see which routine is fastest.
250: Level: intermediate
252: .seealso: PetscIntSortSemiOrdered(), PetscSortReal(), PetscSortIntWithPermutation()
253: @*/
254: PetscErrorCode PetscSortInt(PetscInt n,PetscInt X[])
255: {
256: PetscInt pivot,t1;
259: QuickSort1(PetscSortInt,X,n,pivot,t1);
260: return 0;
261: }
263: /*@
264: PetscSortReverseInt - Sorts an array of integers in place in decreasing order.
266: Not Collective
268: Input Parameters:
269: + n - number of values
270: - X - array of integers
272: Level: intermediate
274: .seealso: PetscIntSortSemiOrdered(), PetscSortInt(), PetscSortIntWithPermutation()
275: @*/
276: PetscErrorCode PetscSortReverseInt(PetscInt n,PetscInt X[])
277: {
278: PetscInt pivot,t1;
281: QuickSortReverse1(PetscSortReverseInt,X,n,pivot,t1);
282: return 0;
283: }
285: /*@
286: PetscSortedRemoveDupsInt - Removes all duplicate entries of a sorted input array
288: Not Collective
290: Input Parameters:
291: + n - number of values
292: - X - sorted array of integers
294: Output Parameter:
295: . n - number of non-redundant values
297: Level: intermediate
299: .seealso: PetscSortInt()
300: @*/
301: PetscErrorCode PetscSortedRemoveDupsInt(PetscInt *n,PetscInt X[])
302: {
303: PetscInt i,s = 0,N = *n, b = 0;
307: for (i=0; i<N-1; i++) {
308: if (X[b+s+1] != X[b]) {
309: X[b+1] = X[b+s+1]; b++;
310: } else s++;
311: }
312: *n = N - s;
313: return 0;
314: }
316: /*@
317: PetscSortedCheckDupsInt - Checks if a sorted integer array has duplicates
319: Not Collective
321: Input Parameters:
322: + n - number of values
323: - X - sorted array of integers
325: Output Parameter:
326: . dups - True if the array has dups, otherwise false
328: Level: intermediate
331: @*/
332: PetscErrorCode PetscSortedCheckDupsInt(PetscInt n,const PetscInt X[],PetscBool *flg)
333: {
334: PetscInt i;
337: *flg = PETSC_FALSE;
338: for (i=0; i<n-1; i++) {
339: if (X[i+1] == X[i]) {
340: *flg = PETSC_TRUE;
341: break;
342: }
343: }
344: return 0;
345: }
347: /*@
348: PetscSortRemoveDupsInt - Sorts an array of integers in place in increasing order removes all duplicate entries
350: Not Collective
352: Input Parameters:
353: + n - number of values
354: - X - array of integers
356: Output Parameter:
357: . n - number of non-redundant values
359: Level: intermediate
361: .seealso: PetscIntSortSemiOrdered(), PetscSortReal(), PetscSortIntWithPermutation(), PetscSortInt(), PetscSortedRemoveDupsInt()
362: @*/
363: PetscErrorCode PetscSortRemoveDupsInt(PetscInt *n,PetscInt X[])
364: {
366: PetscSortInt(*n,X);
367: PetscSortedRemoveDupsInt(n,X);
368: return 0;
369: }
371: /*@
372: PetscFindInt - Finds integer in a sorted array of integers
374: Not Collective
376: Input Parameters:
377: + key - the integer to locate
378: . n - number of values in the array
379: - X - array of integers
381: Output Parameter:
382: . loc - the location if found, otherwise -(slot+1) where slot is the place the value would go
384: Level: intermediate
386: .seealso: PetscIntSortSemiOrdered(), PetscSortInt(), PetscSortIntWithArray(), PetscSortRemoveDupsInt()
387: @*/
388: PetscErrorCode PetscFindInt(PetscInt key, PetscInt n, const PetscInt X[], PetscInt *loc)
389: {
390: PetscInt lo = 0,hi = n;
393: if (!n) {*loc = -1; return 0;}
396: while (hi - lo > 1) {
397: PetscInt mid = lo + (hi - lo)/2;
398: if (key < X[mid]) hi = mid;
399: else lo = mid;
400: }
401: *loc = key == X[lo] ? lo : -(lo + (key > X[lo]) + 1);
402: return 0;
403: }
405: /*@
408: Not Collective
410: Input Parameters:
411: + n - number of values in the array
412: - X - array of integers
414: Output Parameter:
415: . dups - True if the array has dups, otherwise false
417: Level: intermediate
419: .seealso: PetscSortRemoveDupsInt(), PetscSortedCheckDupsInt()
420: @*/
422: {
423: PetscInt i;
424: PetscHSetI ht;
425: PetscBool missing;
429: *dups = PETSC_FALSE;
430: if (n > 1) {
431: PetscHSetICreate(&ht);
432: PetscHSetIResize(ht,n);
433: for (i=0; i<n; i++) {
434: PetscHSetIQueryAdd(ht,X[i],&missing);
435: if (!missing) {*dups = PETSC_TRUE; break;}
436: }
437: PetscHSetIDestroy(&ht);
438: }
439: return 0;
440: }
442: /*@
443: PetscFindMPIInt - Finds MPI integer in a sorted array of integers
445: Not Collective
447: Input Parameters:
448: + key - the integer to locate
449: . n - number of values in the array
450: - X - array of integers
452: Output Parameter:
453: . loc - the location if found, otherwise -(slot+1) where slot is the place the value would go
455: Level: intermediate
457: .seealso: PetscMPIIntSortSemiOrdered(), PetscSortInt(), PetscSortIntWithArray(), PetscSortRemoveDupsInt()
458: @*/
459: PetscErrorCode PetscFindMPIInt(PetscMPIInt key, PetscInt n, const PetscMPIInt X[], PetscInt *loc)
460: {
461: PetscInt lo = 0,hi = n;
464: if (!n) {*loc = -1; return 0;}
467: while (hi - lo > 1) {
468: PetscInt mid = lo + (hi - lo)/2;
469: if (key < X[mid]) hi = mid;
470: else lo = mid;
471: }
472: *loc = key == X[lo] ? lo : -(lo + (key > X[lo]) + 1);
473: return 0;
474: }
476: /*@
477: PetscSortIntWithArray - Sorts an array of integers in place in increasing order;
478: changes a second array to match the sorted first array.
480: Not Collective
482: Input Parameters:
483: + n - number of values
484: . X - array of integers
485: - Y - second array of integers
487: Level: intermediate
489: .seealso: PetscIntSortSemiOrderedWithArray(), PetscSortReal(), PetscSortIntWithPermutation(), PetscSortInt(), PetscSortIntWithCountArray()
490: @*/
491: PetscErrorCode PetscSortIntWithArray(PetscInt n,PetscInt X[],PetscInt Y[])
492: {
493: PetscInt pivot,t1,t2;
495: QuickSort2(PetscSortIntWithArray,X,Y,n,pivot,t1,t2);
496: return 0;
497: }
499: /*@
500: PetscSortIntWithArrayPair - Sorts an array of integers in place in increasing order;
501: changes a pair of integer arrays to match the sorted first array.
503: Not Collective
505: Input Parameters:
506: + n - number of values
507: . X - array of integers
508: . Y - second array of integers (first array of the pair)
509: - Z - third array of integers (second array of the pair)
511: Level: intermediate
513: .seealso: PetscSortReal(), PetscSortIntWithPermutation(), PetscSortIntWithArray(), PetscIntSortSemiOrdered(), PetscSortIntWithIntCountArrayPair()
514: @*/
515: PetscErrorCode PetscSortIntWithArrayPair(PetscInt n,PetscInt X[],PetscInt Y[],PetscInt Z[])
516: {
517: PetscInt pivot,t1,t2,t3;
519: QuickSort3(PetscSortIntWithArrayPair,X,Y,Z,n,pivot,t1,t2,t3);
520: return 0;
521: }
523: /*@
524: PetscSortIntWithCountArray - Sorts an array of integers in place in increasing order;
525: changes a second array to match the sorted first array.
527: Not Collective
529: Input Parameters:
530: + n - number of values
531: . X - array of integers
532: - Y - second array of PetscCounts (signed integers)
534: Level: intermediate
536: .seealso: PetscIntSortSemiOrderedWithArray(), PetscSortReal(), PetscSortIntPermutation(), PetscSortInt(), PetscSortIntWithArray()
537: @*/
538: PetscErrorCode PetscSortIntWithCountArray(PetscCount n,PetscInt X[],PetscCount Y[])
539: {
540: PetscInt pivot,t1;
541: PetscCount t2;
543: QuickSort2(PetscSortIntWithCountArray,X,Y,n,pivot,t1,t2);
544: return 0;
545: }
547: /*@
548: PetscSortIntWithIntCountArrayPair - Sorts an array of integers in place in increasing order;
549: changes an integer array and a PetscCount array to match the sorted first array.
551: Not Collective
553: Input Parameters:
554: + n - number of values
555: . X - array of integers
556: . Y - second array of integers (first array of the pair)
557: - Z - third array of PetscCounts (second array of the pair)
559: Level: intermediate
561: Notes:
562: Usually X, Y are matrix row/column indices, and Z is a permutation array and therefore Z's type is PetscCount to allow 2B+ nonzeros even with 32-bit PetscInt.
564: .seealso: PetscSortReal(), PetscSortIntPermutation(), PetscSortIntWithArray(), PetscIntSortSemiOrdered(), PetscSortIntWithArrayPair()
565: @*/
566: PetscErrorCode PetscSortIntWithIntCountArrayPair(PetscCount n,PetscInt X[],PetscInt Y[],PetscCount Z[])
567: {
568: PetscInt pivot,t1,t2; /* pivot is take from X[], so its type is still PetscInt */
569: PetscCount t3; /* temp for Z[] */
571: QuickSort3(PetscSortIntWithIntCountArrayPair,X,Y,Z,n,pivot,t1,t2,t3);
572: return 0;
573: }
575: /*@
576: PetscSortedMPIInt - Determines whether the array is sorted.
578: Not Collective
580: Input Parameters:
581: + n - number of values
582: - X - array of integers
584: Output Parameters:
585: . sorted - flag whether the array is sorted
587: Level: intermediate
589: .seealso: PetscMPIIntSortSemiOrdered(), PetscSortMPIInt(), PetscSortedInt(), PetscSortedReal()
590: @*/
591: PetscErrorCode PetscSortedMPIInt(PetscInt n,const PetscMPIInt X[],PetscBool *sorted)
592: {
593: PetscSorted(n,X,*sorted);
594: return 0;
595: }
597: /*@
598: PetscSortMPIInt - Sorts an array of MPI integers in place in increasing order.
600: Not Collective
602: Input Parameters:
603: + n - number of values
604: - X - array of integers
606: Level: intermediate
608: Notes:
609: This function serves as an alternative to PetscMPIIntSortSemiOrdered(), and may perform faster especially if the array
610: is completely random. There are exceptions to this and so it is __highly__ recommended that the user benchmark their
611: code to see which routine is fastest.
613: .seealso: PetscMPIIntSortSemiOrdered(), PetscSortReal(), PetscSortIntWithPermutation()
614: @*/
615: PetscErrorCode PetscSortMPIInt(PetscInt n,PetscMPIInt X[])
616: {
617: PetscMPIInt pivot,t1;
619: QuickSort1(PetscSortMPIInt,X,n,pivot,t1);
620: return 0;
621: }
623: /*@
624: PetscSortRemoveDupsMPIInt - Sorts an array of MPI integers in place in increasing order removes all duplicate entries
626: Not Collective
628: Input Parameters:
629: + n - number of values
630: - X - array of integers
632: Output Parameter:
633: . n - number of non-redundant values
635: Level: intermediate
637: .seealso: PetscSortReal(), PetscSortIntWithPermutation(), PetscSortInt()
638: @*/
639: PetscErrorCode PetscSortRemoveDupsMPIInt(PetscInt *n,PetscMPIInt X[])
640: {
641: PetscInt s = 0,N = *n,b = 0;
643: PetscSortMPIInt(N,X);
644: for (PetscInt i=0; i<N-1; i++) {
645: if (X[b+s+1] != X[b]) {
646: X[b+1] = X[b+s+1]; b++;
647: } else s++;
648: }
649: *n = N - s;
650: return 0;
651: }
653: /*@
654: PetscSortMPIIntWithArray - Sorts an array of integers in place in increasing order;
655: changes a second array to match the sorted first array.
657: Not Collective
659: Input Parameters:
660: + n - number of values
661: . X - array of integers
662: - Y - second array of integers
664: Level: intermediate
666: .seealso: PetscMPIIntSortSemiOrderedWithArray(), PetscSortReal(), PetscSortIntWithPermutation(), PetscSortInt()
667: @*/
668: PetscErrorCode PetscSortMPIIntWithArray(PetscMPIInt n,PetscMPIInt X[],PetscMPIInt Y[])
669: {
670: PetscMPIInt pivot,t1,t2;
672: QuickSort2(PetscSortMPIIntWithArray,X,Y,n,pivot,t1,t2);
673: return 0;
674: }
676: /*@
677: PetscSortMPIIntWithIntArray - Sorts an array of MPI integers in place in increasing order;
678: changes a second array of Petsc intergers to match the sorted first array.
680: Not Collective
682: Input Parameters:
683: + n - number of values
684: . X - array of MPI integers
685: - Y - second array of Petsc integers
687: Level: intermediate
689: Notes: this routine is useful when one needs to sort MPI ranks with other integer arrays.
691: .seealso: PetscSortMPIIntWithArray(), PetscIntSortSemiOrderedWithArray(), PetscTimSortWithArray()
692: @*/
693: PetscErrorCode PetscSortMPIIntWithIntArray(PetscMPIInt n,PetscMPIInt X[],PetscInt Y[])
694: {
695: PetscMPIInt pivot,t1;
696: PetscInt t2;
698: QuickSort2(PetscSortMPIIntWithIntArray,X,Y,n,pivot,t1,t2);
699: return 0;
700: }
702: /*@
703: PetscSortIntWithScalarArray - Sorts an array of integers in place in increasing order;
704: changes a second SCALAR array to match the sorted first INTEGER array.
706: Not Collective
708: Input Parameters:
709: + n - number of values
710: . X - array of integers
711: - Y - second array of scalars
713: Level: intermediate
715: .seealso: PetscTimSortWithArray(), PetscSortReal(), PetscSortIntWithPermutation(), PetscSortInt(), PetscSortIntWithArray()
716: @*/
717: PetscErrorCode PetscSortIntWithScalarArray(PetscInt n,PetscInt X[],PetscScalar Y[])
718: {
719: PetscInt pivot,t1;
720: PetscScalar t2;
722: QuickSort2(PetscSortIntWithScalarArray,X,Y,n,pivot,t1,t2);
723: return 0;
724: }
726: /*@C
727: PetscSortIntWithDataArray - Sorts an array of integers in place in increasing order;
728: changes a second array to match the sorted first INTEGER array. Unlike other sort routines, the user must
729: provide workspace (the size of an element in the data array) to use when sorting.
731: Not Collective
733: Input Parameters:
734: + n - number of values
735: . X - array of integers
736: . Y - second array of data
737: . size - sizeof elements in the data array in bytes
738: - t2 - workspace of "size" bytes used when sorting
740: Level: intermediate
742: .seealso: PetscTimSortWithArray(), PetscSortReal(), PetscSortIntWithPermutation(), PetscSortInt(), PetscSortIntWithArray()
743: @*/
744: PetscErrorCode PetscSortIntWithDataArray(PetscInt n,PetscInt X[],void *Y,size_t size,void *t2)
745: {
746: char *YY = (char*)Y;
747: PetscInt t1,pivot,hi = n-1;
749: if (n<8) {
750: for (PetscInt i=0; i<n; i++) {
751: pivot = X[i];
752: for (PetscInt j=i+1; j<n; j++) {
753: if (pivot > X[j]) {
754: SWAP2Data(X[i],X[j],YY+size*i,YY+size*j,t1,t2,size);
755: pivot = X[i];
756: }
757: }
758: }
759: } else {
760: /* Two way partition */
761: PetscInt l = 0,r = hi;
763: pivot = X[MEDIAN(X,hi)];
764: while (1) {
765: while (X[l] < pivot) l++;
766: while (X[r] > pivot) r--;
767: if (l >= r) {r++; break;}
768: SWAP2Data(X[l],X[r],YY+size*l,YY+size*r,t1,t2,size);
769: l++;
770: r--;
771: }
772: PetscSortIntWithDataArray(l,X,Y,size,t2);
773: PetscSortIntWithDataArray(hi-r+1,X+r,YY+size*r,size,t2);
774: }
775: return 0;
776: }
778: /*@
779: PetscMergeIntArray - Merges two SORTED integer arrays, removes duplicate elements.
781: Not Collective
783: Input Parameters:
784: + an - number of values in the first array
785: . aI - first sorted array of integers
786: . bn - number of values in the second array
787: - bI - second array of integers
789: Output Parameters:
790: + n - number of values in the merged array
791: - L - merged sorted array, this is allocated if an array is not provided
793: Level: intermediate
795: .seealso: PetscSortReal(), PetscSortIntWithPermutation(), PetscSortInt(), PetscSortIntWithArray()
796: @*/
797: PetscErrorCode PetscMergeIntArray(PetscInt an,const PetscInt aI[], PetscInt bn, const PetscInt bI[], PetscInt *n, PetscInt **L)
798: {
799: PetscInt *L_ = *L, ak, bk, k;
801: if (!L_) {
802: PetscMalloc1(an+bn, L);
803: L_ = *L;
804: }
805: k = ak = bk = 0;
806: while (ak < an && bk < bn) {
807: if (aI[ak] == bI[bk]) {
808: L_[k] = aI[ak];
809: ++ak;
810: ++bk;
811: ++k;
812: } else if (aI[ak] < bI[bk]) {
813: L_[k] = aI[ak];
814: ++ak;
815: ++k;
816: } else {
817: L_[k] = bI[bk];
818: ++bk;
819: ++k;
820: }
821: }
822: if (ak < an) {
823: PetscArraycpy(L_+k,aI+ak,an-ak);
824: k += (an-ak);
825: }
826: if (bk < bn) {
827: PetscArraycpy(L_+k,bI+bk,bn-bk);
828: k += (bn-bk);
829: }
830: *n = k;
831: return 0;
832: }
834: /*@
835: PetscMergeIntArrayPair - Merges two SORTED integer arrays that share NO common values along with an additional array of integers.
836: The additional arrays are the same length as sorted arrays and are merged
837: in the order determined by the merging of the sorted pair.
839: Not Collective
841: Input Parameters:
842: + an - number of values in the first array
843: . aI - first sorted array of integers
844: . aJ - first additional array of integers
845: . bn - number of values in the second array
846: . bI - second array of integers
847: - bJ - second additional of integers
849: Output Parameters:
850: + n - number of values in the merged array (== an + bn)
851: . L - merged sorted array
852: - J - merged additional array
854: Notes:
855: if L or J point to non-null arrays then this routine will assume they are of the approproate size and use them, otherwise this routine will allocate space for them
856: Level: intermediate
858: .seealso: PetscIntSortSemiOrdered(), PetscSortReal(), PetscSortIntWithPermutation(), PetscSortInt(), PetscSortIntWithArray()
859: @*/
860: PetscErrorCode PetscMergeIntArrayPair(PetscInt an,const PetscInt aI[], const PetscInt aJ[], PetscInt bn, const PetscInt bI[], const PetscInt bJ[], PetscInt *n, PetscInt **L, PetscInt **J)
861: {
862: PetscInt n_, *L_, *J_, ak, bk, k;
866: n_ = an + bn;
867: *n = n_;
868: if (!*L) {
869: PetscMalloc1(n_, L);
870: }
871: L_ = *L;
872: if (!*J) {
873: PetscMalloc1(n_, J);
874: }
875: J_ = *J;
876: k = ak = bk = 0;
877: while (ak < an && bk < bn) {
878: if (aI[ak] <= bI[bk]) {
879: L_[k] = aI[ak];
880: J_[k] = aJ[ak];
881: ++ak;
882: ++k;
883: } else {
884: L_[k] = bI[bk];
885: J_[k] = bJ[bk];
886: ++bk;
887: ++k;
888: }
889: }
890: if (ak < an) {
891: PetscArraycpy(L_+k,aI+ak,an-ak);
892: PetscArraycpy(J_+k,aJ+ak,an-ak);
893: k += (an-ak);
894: }
895: if (bk < bn) {
896: PetscArraycpy(L_+k,bI+bk,bn-bk);
897: PetscArraycpy(J_+k,bJ+bk,bn-bk);
898: }
899: return 0;
900: }
902: /*@
903: PetscMergeMPIIntArray - Merges two SORTED integer arrays.
905: Not Collective
907: Input Parameters:
908: + an - number of values in the first array
909: . aI - first sorted array of integers
910: . bn - number of values in the second array
911: - bI - second array of integers
913: Output Parameters:
914: + n - number of values in the merged array (<= an + bn)
915: - L - merged sorted array, allocated if address of NULL pointer is passed
917: Level: intermediate
919: .seealso: PetscIntSortSemiOrdered(), PetscSortReal(), PetscSortIntWithPermutation(), PetscSortInt(), PetscSortIntWithArray()
920: @*/
921: PetscErrorCode PetscMergeMPIIntArray(PetscInt an,const PetscMPIInt aI[],PetscInt bn,const PetscMPIInt bI[],PetscInt *n,PetscMPIInt **L)
922: {
923: PetscInt ai,bi,k;
925: if (!*L) PetscMalloc1((an+bn),L);
926: for (ai=0,bi=0,k=0; ai<an || bi<bn;) {
927: PetscInt t = -1;
928: for (; ai<an && (!bn || aI[ai] <= bI[bi]); ai++) (*L)[k++] = t = aI[ai];
929: for (; bi<bn && bI[bi] == t; bi++);
930: for (; bi<bn && (!an || bI[bi] <= aI[ai]); bi++) (*L)[k++] = t = bI[bi];
931: for (; ai<an && aI[ai] == t; ai++);
932: }
933: *n = k;
934: return 0;
935: }
937: /*@C
938: PetscProcessTree - Prepares tree data to be displayed graphically
940: Not Collective
942: Input Parameters:
943: + n - number of values
944: . mask - indicates those entries in the tree, location 0 is always masked
945: - parentid - indicates the parent of each entry
947: Output Parameters:
948: + Nlevels - the number of levels
949: . Level - for each node tells its level
950: . Levelcnts - the number of nodes on each level
951: . Idbylevel - a list of ids on each of the levels, first level followed by second etc
952: - Column - for each id tells its column index
954: Level: developer
956: Notes:
957: This code is not currently used
959: .seealso: PetscSortReal(), PetscSortIntWithPermutation()
960: @*/
961: PetscErrorCode PetscProcessTree(PetscInt n,const PetscBool mask[],const PetscInt parentid[],PetscInt *Nlevels,PetscInt **Level,PetscInt **Levelcnt,PetscInt **Idbylevel,PetscInt **Column)
962: {
963: PetscInt i,j,cnt,nmask = 0,nlevels = 0,*level,*levelcnt,levelmax = 0,*workid,*workparentid,tcnt = 0,*idbylevel,*column;
964: PetscBool done = PETSC_FALSE;
967: for (i=0; i<n; i++) {
968: if (mask[i]) continue;
971: }
973: for (i=0; i<n; i++) {
974: if (!mask[i]) nmask++;
975: }
977: /* determine the level in the tree of each node */
978: PetscCalloc1(n,&level);
980: level[0] = 1;
981: while (!done) {
982: done = PETSC_TRUE;
983: for (i=0; i<n; i++) {
984: if (mask[i]) continue;
985: if (!level[i] && level[parentid[i]]) level[i] = level[parentid[i]] + 1;
986: else if (!level[i]) done = PETSC_FALSE;
987: }
988: }
989: for (i=0; i<n; i++) {
990: level[i]--;
991: nlevels = PetscMax(nlevels,level[i]);
992: }
994: /* count the number of nodes on each level and its max */
995: PetscCalloc1(nlevels,&levelcnt);
996: for (i=0; i<n; i++) {
997: if (mask[i]) continue;
998: levelcnt[level[i]-1]++;
999: }
1000: for (i=0; i<nlevels;i++) levelmax = PetscMax(levelmax,levelcnt[i]);
1002: /* for each level sort the ids by the parent id */
1003: PetscMalloc2(levelmax,&workid,levelmax,&workparentid);
1004: PetscMalloc1(nmask,&idbylevel);
1005: for (j=1; j<=nlevels;j++) {
1006: cnt = 0;
1007: for (i=0; i<n; i++) {
1008: if (mask[i]) continue;
1009: if (level[i] != j) continue;
1010: workid[cnt] = i;
1011: workparentid[cnt++] = parentid[i];
1012: }
1013: /* PetscIntView(cnt,workparentid,0);
1014: PetscIntView(cnt,workid,0);
1015: PetscSortIntWithArray(cnt,workparentid,workid);
1016: PetscIntView(cnt,workparentid,0);
1017: PetscIntView(cnt,workid,0);*/
1018: PetscArraycpy(idbylevel+tcnt,workid,cnt);
1019: tcnt += cnt;
1020: }
1022: PetscFree2(workid,workparentid);
1024: /* for each node list its column */
1025: PetscMalloc1(n,&column);
1026: cnt = 0;
1027: for (j=0; j<nlevels; j++) {
1028: for (i=0; i<levelcnt[j]; i++) {
1029: column[idbylevel[cnt++]] = i;
1030: }
1031: }
1033: *Nlevels = nlevels;
1034: *Level = level;
1035: *Levelcnt = levelcnt;
1036: *Idbylevel = idbylevel;
1037: *Column = column;
1038: return 0;
1039: }
1041: /*@
1042: PetscParallelSortedInt - Check whether an integer array, distributed over a communicator, is globally sorted.
1044: Collective
1046: Input Parameters:
1047: + comm - the MPI communicator
1048: . n - the local number of integers
1049: - keys - the local array of integers
1051: Output Parameters:
1052: . is_sorted - whether the array is globally sorted
1054: Level: developer
1056: .seealso: PetscParallelSortInt()
1057: @*/
1058: PetscErrorCode PetscParallelSortedInt(MPI_Comm comm, PetscInt n, const PetscInt keys[], PetscBool *is_sorted)
1059: {
1060: PetscBool sorted;
1061: PetscInt i, min, max, prevmax;
1062: PetscMPIInt rank;
1064: sorted = PETSC_TRUE;
1065: min = PETSC_MAX_INT;
1066: max = PETSC_MIN_INT;
1067: if (n) {
1068: min = keys[0];
1069: max = keys[0];
1070: }
1071: for (i = 1; i < n; i++) {
1072: if (keys[i] < keys[i - 1]) break;
1073: min = PetscMin(min,keys[i]);
1074: max = PetscMax(max,keys[i]);
1075: }
1076: if (i < n) sorted = PETSC_FALSE;
1077: prevmax = PETSC_MIN_INT;
1078: MPI_Exscan(&max, &prevmax, 1, MPIU_INT, MPI_MAX, comm);
1079: MPI_Comm_rank(comm, &rank);
1080: if (rank == 0) prevmax = PETSC_MIN_INT;
1081: if (prevmax > min) sorted = PETSC_FALSE;
1082: MPI_Allreduce(&sorted, is_sorted, 1, MPIU_BOOL, MPI_LAND, comm);
1083: return 0;
1084: }