Actual source code: sorti.c

  1: /*
  2:    This file contains routines for sorting integers. Values are sorted in place.
  3:    One can use src/sys/tests/ex52.c for benchmarking.
  4:  */
  5: #include <petsc/private/petscimpl.h>
  6: #include <petsc/private/hashseti.h>

  8: #define MEDIAN3(v, a, b, c) (v[a] < v[b] ? (v[b] < v[c] ? (b) : (v[a] < v[c] ? (c) : (a))) : (v[c] < v[b] ? (b) : (v[a] < v[c] ? (a) : (c))))

 10: #define MEDIAN(v, right) MEDIAN3(v, right / 4, right / 2, right / 4 * 3)

 12: /* Swap one, two or three pairs. Each pair can have its own type */
 13: #define SWAP1(a, b, t1) \
 14:   do { \
 15:     t1 = a; \
 16:     a  = b; \
 17:     b  = t1; \
 18:   } while (0)
 19: #define SWAP2(a, b, c, d, t1, t2) \
 20:   do { \
 21:     t1 = a; \
 22:     a  = b; \
 23:     b  = t1; \
 24:     t2 = c; \
 25:     c  = d; \
 26:     d  = t2; \
 27:   } while (0)
 28: #define SWAP3(a, b, c, d, e, f, t1, t2, t3) \
 29:   do { \
 30:     t1 = a; \
 31:     a  = b; \
 32:     b  = t1; \
 33:     t2 = c; \
 34:     c  = d; \
 35:     d  = t2; \
 36:     t3 = e; \
 37:     e  = f; \
 38:     f  = t3; \
 39:   } while (0)

 41: /* Swap a & b, *c & *d. c, d, t2 are pointers to a type of size <siz> */
 42: #define SWAP2Data(a, b, c, d, t1, t2, siz) \
 43:   do { \
 44:     t1 = a; \
 45:     a  = b; \
 46:     b  = t1; \
 47:     PetscCall(PetscMemcpy(t2, c, siz)); \
 48:     PetscCall(PetscMemcpy(c, d, siz)); \
 49:     PetscCall(PetscMemcpy(d, t2, siz)); \
 50:   } while (0)

 52: /*
 53:    Partition X[lo,hi] into two parts: X[lo,l) <= pivot; X[r,hi] > pivot

 55:    Input Parameters:
 56:     + X         - array to partition
 57:     . pivot     - a pivot of X[]
 58:     . t1        - temp variable for X
 59:     - lo,hi     - lower and upper bound of the array

 61:    Output Parameters:
 62:     + l,r       - of type PetscInt

 64:    Note:
 65:     The TwoWayPartition2/3 variants also partition other arrays along with X.
 66:     These arrays can have different types, so they provide their own temp t2,t3
 67:  */
 68: #define TwoWayPartition1(X, pivot, t1, lo, hi, l, r) \
 69:   do { \
 70:     l = lo; \
 71:     r = hi; \
 72:     while (1) { \
 73:       while (X[l] < pivot) l++; \
 74:       while (X[r] > pivot) r--; \
 75:       if (l >= r) { \
 76:         r++; \
 77:         break; \
 78:       } \
 79:       SWAP1(X[l], X[r], t1); \
 80:       l++; \
 81:       r--; \
 82:     } \
 83:   } while (0)

 85: /*
 86:    Partition X[lo,hi] into two parts: X[lo,l) >= pivot; X[r,hi] < pivot

 88:    Input Parameters:
 89:     + X         - array to partition
 90:     . pivot     - a pivot of X[]
 91:     . t1        - temp variable for X
 92:     - lo,hi     - lower and upper bound of the array

 94:    Output Parameters:
 95:     + l,r       - of type PetscInt

 97:    Note:
 98:     The TwoWayPartition2/3 variants also partition other arrays along with X.
 99:     These arrays can have different types, so they provide their own temp t2,t3
100:  */
101: #define TwoWayPartitionReverse1(X, pivot, t1, lo, hi, l, r) \
102:   do { \
103:     l = lo; \
104:     r = hi; \
105:     while (1) { \
106:       while (X[l] > pivot) l++; \
107:       while (X[r] < pivot) r--; \
108:       if (l >= r) { \
109:         r++; \
110:         break; \
111:       } \
112:       SWAP1(X[l], X[r], t1); \
113:       l++; \
114:       r--; \
115:     } \
116:   } while (0)

118: #define TwoWayPartition2(X, Y, pivot, t1, t2, lo, hi, l, r) \
119:   do { \
120:     l = lo; \
121:     r = hi; \
122:     while (1) { \
123:       while (X[l] < pivot) l++; \
124:       while (X[r] > pivot) r--; \
125:       if (l >= r) { \
126:         r++; \
127:         break; \
128:       } \
129:       SWAP2(X[l], X[r], Y[l], Y[r], t1, t2); \
130:       l++; \
131:       r--; \
132:     } \
133:   } while (0)

135: #define TwoWayPartition3(X, Y, Z, pivot, t1, t2, t3, lo, hi, l, r) \
136:   do { \
137:     l = lo; \
138:     r = hi; \
139:     while (1) { \
140:       while (X[l] < pivot) l++; \
141:       while (X[r] > pivot) r--; \
142:       if (l >= r) { \
143:         r++; \
144:         break; \
145:       } \
146:       SWAP3(X[l], X[r], Y[l], Y[r], Z[l], Z[r], t1, t2, t3); \
147:       l++; \
148:       r--; \
149:     } \
150:   } while (0)

152: /* Templates for similar functions used below */
153: #define QuickSort1(FuncName, X, n, pivot, t1) \
154:   do { \
155:     PetscCount i, j, p, l, r, hi = n - 1; \
156:     if (n < 8) { \
157:       for (i = 0; i < n; i++) { \
158:         pivot = X[i]; \
159:         for (j = i + 1; j < n; j++) { \
160:           if (pivot > X[j]) { \
161:             SWAP1(X[i], X[j], t1); \
162:             pivot = X[i]; \
163:           } \
164:         } \
165:       } \
166:     } else { \
167:       p     = MEDIAN(X, hi); \
168:       pivot = X[p]; \
169:       TwoWayPartition1(X, pivot, t1, 0, hi, l, r); \
170:       PetscCall(FuncName(l, X)); \
171:       PetscCall(FuncName(hi - r + 1, X + r)); \
172:     } \
173:   } while (0)

175: /* Templates for similar functions used below */
176: #define QuickSortReverse1(FuncName, X, n, pivot, t1) \
177:   do { \
178:     PetscCount i, j, p, l, r, hi = n - 1; \
179:     if (n < 8) { \
180:       for (i = 0; i < n; i++) { \
181:         pivot = X[i]; \
182:         for (j = i + 1; j < n; j++) { \
183:           if (pivot < X[j]) { \
184:             SWAP1(X[i], X[j], t1); \
185:             pivot = X[i]; \
186:           } \
187:         } \
188:       } \
189:     } else { \
190:       p     = MEDIAN(X, hi); \
191:       pivot = X[p]; \
192:       TwoWayPartitionReverse1(X, pivot, t1, 0, hi, l, r); \
193:       PetscCall(FuncName(l, X)); \
194:       PetscCall(FuncName(hi - r + 1, X + r)); \
195:     } \
196:   } while (0)

198: #define QuickSort2(FuncName, X, Y, n, pivot, t1, t2) \
199:   do { \
200:     PetscCount i, j, p, l, r, hi = n - 1; \
201:     if (n < 8) { \
202:       for (i = 0; i < n; i++) { \
203:         pivot = X[i]; \
204:         for (j = i + 1; j < n; j++) { \
205:           if (pivot > X[j]) { \
206:             SWAP2(X[i], X[j], Y[i], Y[j], t1, t2); \
207:             pivot = X[i]; \
208:           } \
209:         } \
210:       } \
211:     } else { \
212:       p     = MEDIAN(X, hi); \
213:       pivot = X[p]; \
214:       TwoWayPartition2(X, Y, pivot, t1, t2, 0, hi, l, r); \
215:       PetscCall(FuncName(l, X, Y)); \
216:       PetscCall(FuncName(hi - r + 1, X + r, Y + r)); \
217:     } \
218:   } while (0)

220: #define QuickSort3(FuncName, X, Y, Z, n, pivot, t1, t2, t3) \
221:   do { \
222:     PetscCount i, j, p, l, r, hi = n - 1; \
223:     if (n < 8) { \
224:       for (i = 0; i < n; i++) { \
225:         pivot = X[i]; \
226:         for (j = i + 1; j < n; j++) { \
227:           if (pivot > X[j]) { \
228:             SWAP3(X[i], X[j], Y[i], Y[j], Z[i], Z[j], t1, t2, t3); \
229:             pivot = X[i]; \
230:           } \
231:         } \
232:       } \
233:     } else { \
234:       p     = MEDIAN(X, hi); \
235:       pivot = X[p]; \
236:       TwoWayPartition3(X, Y, Z, pivot, t1, t2, t3, 0, hi, l, r); \
237:       PetscCall(FuncName(l, X, Y, Z)); \
238:       PetscCall(FuncName(hi - r + 1, X + r, Y + r, Z + r)); \
239:     } \
240:   } while (0)

242: /*@
243:   PetscSortedInt - Determines whether the `PetscInt` array is sorted.

245:   Not Collective

247:   Input Parameters:
248: + n - number of values
249: - X - array of integers

251:   Output Parameter:
252: . sorted - flag whether the array is sorted

254:   Level: intermediate

256: .seealso: `PetscSortInt()`, `PetscSortedMPIInt()`, `PetscSortedReal()`
257: @*/
258: PetscErrorCode PetscSortedInt(PetscInt n, const PetscInt X[], PetscBool *sorted)
259: {
260:   PetscFunctionBegin;
261:   if (n) PetscAssertPointer(X, 2);
262:   PetscAssertPointer(sorted, 3);
263:   PetscSorted(n, X, *sorted);
264:   PetscFunctionReturn(PETSC_SUCCESS);
265: }

267: /*@
268:   PetscSortedInt64 - Determines whether the `PetscInt64` array is sorted.

270:   Not Collective

272:   Input Parameters:
273: + n - number of values
274: - X - array of integers

276:   Output Parameter:
277: . sorted - flag whether the array is sorted

279:   Level: intermediate

281: .seealso: `PetscSortInt64()`, `PetscSortInt()`, `PetscSortedMPIInt()`, `PetscSortedReal()`
282: @*/
283: PetscErrorCode PetscSortedInt64(PetscInt n, const PetscInt64 X[], PetscBool *sorted)
284: {
285:   PetscFunctionBegin;
286:   if (n) PetscAssertPointer(X, 2);
287:   PetscAssertPointer(sorted, 3);
288:   PetscSorted(n, X, *sorted);
289:   PetscFunctionReturn(PETSC_SUCCESS);
290: }

292: /*@
293:   PetscSortInt - Sorts an array of `PetscInt` in place in increasing order.

295:   Not Collective

297:   Input Parameters:
298: + n - number of values
299: - X - array of integers

301:   Note:
302:   This function serves as an alternative to `PetscIntSortSemiOrdered()`, and may perform faster especially if the array
303:   is completely random. There are exceptions to this and so it is __highly__ recommended that the user benchmark their
304:   code to see which routine is fastest.

306:   Level: intermediate

308: .seealso: `PetscIntSortSemiOrdered()`, `PetscSortReal()`, `PetscSortIntWithPermutation()`
309: @*/
310: PetscErrorCode PetscSortInt(PetscInt n, PetscInt X[])
311: {
312:   PetscInt pivot, t1;

314:   PetscFunctionBegin;
315:   if (n) PetscAssertPointer(X, 2);
316:   QuickSort1(PetscSortInt, X, n, pivot, t1);
317:   PetscFunctionReturn(PETSC_SUCCESS);
318: }

320: /*@
321:   PetscSortInt64 - Sorts an array of `PetscInt64` in place in increasing order.

323:   Not Collective

325:   Input Parameters:
326: + n - number of values
327: - X - array of integers

329:   Notes:
330:   This function sorts `PetscInt64`s assumed to be in completely random order

332:   Level: intermediate

334: .seealso: `PetscSortInt()`
335: @*/
336: PetscErrorCode PetscSortInt64(PetscInt n, PetscInt64 X[])
337: {
338:   PetscCount pivot, t1;

340:   PetscFunctionBegin;
341:   if (n) PetscAssertPointer(X, 2);
342:   QuickSort1(PetscSortInt64, X, n, pivot, t1);
343:   PetscFunctionReturn(PETSC_SUCCESS);
344: }

346: /*@
347:   PetscSortCount - Sorts an array of integers in place in increasing order.

349:   Not Collective

351:   Input Parameters:
352: + n - number of values
353: - X - array of integers

355:   Notes:
356:   This function sorts `PetscCount`s assumed to be in completely random order

358:   Level: intermediate

360: .seealso: `PetscSortInt()`
361: @*/
362: PetscErrorCode PetscSortCount(PetscInt n, PetscCount X[])
363: {
364:   PetscCount pivot, t1;

366:   PetscFunctionBegin;
367:   if (n) PetscAssertPointer(X, 2);
368:   QuickSort1(PetscSortCount, X, n, pivot, t1);
369:   PetscFunctionReturn(PETSC_SUCCESS);
370: }

372: /*@
373:   PetscSortReverseInt - Sorts an array of integers in place in decreasing order.

375:   Not Collective

377:   Input Parameters:
378: + n - number of values
379: - X - array of integers

381:   Level: intermediate

383: .seealso: `PetscIntSortSemiOrdered()`, `PetscSortInt()`, `PetscSortIntWithPermutation()`
384: @*/
385: PetscErrorCode PetscSortReverseInt(PetscInt n, PetscInt X[])
386: {
387:   PetscInt pivot, t1;

389:   PetscFunctionBegin;
390:   if (n) PetscAssertPointer(X, 2);
391:   QuickSortReverse1(PetscSortReverseInt, X, n, pivot, t1);
392:   PetscFunctionReturn(PETSC_SUCCESS);
393: }

395: /*@
396:   PetscSortedRemoveDupsInt - Removes all duplicate entries of a sorted `PetscInt` array

398:   Not Collective

400:   Input Parameters:
401: + n - number of values
402: - X - sorted array of integers

404:   Output Parameter:
405: . n - number of non-redundant values

407:   Level: intermediate

409: .seealso: `PetscSortInt()`
410: @*/
411: PetscErrorCode PetscSortedRemoveDupsInt(PetscInt *n, PetscInt X[])
412: {
413:   PetscInt i, s = 0, N = *n, b = 0;

415:   PetscFunctionBegin;
416:   PetscAssertPointer(n, 1);
417:   PetscCheckSorted(*n, X);
418:   for (i = 0; i < N - 1; i++) {
419:     if (X[b + s + 1] != X[b]) {
420:       X[b + 1] = X[b + s + 1];
421:       b++;
422:     } else s++;
423:   }
424:   *n = N - s;
425:   PetscFunctionReturn(PETSC_SUCCESS);
426: }

428: /*@
429:   PetscSortedCheckDupsInt - Checks if a sorted `PetscInt` array has duplicates

431:   Not Collective

433:   Input Parameters:
434: + n - number of values
435: - X - sorted array of integers

437:   Output Parameter:
438: . flg - True if the array has dups, otherwise false

440:   Level: intermediate

442: .seealso: `PetscSortInt()`, `PetscCheckDupsInt()`, `PetscSortRemoveDupsInt()`, `PetscSortedRemoveDupsInt()`
443: @*/
444: PetscErrorCode PetscSortedCheckDupsInt(PetscInt n, const PetscInt X[], PetscBool *flg)
445: {
446:   PetscInt i;

448:   PetscFunctionBegin;
449:   PetscCheckSorted(n, X);
450:   *flg = PETSC_FALSE;
451:   for (i = 0; i < n - 1; i++) {
452:     if (X[i + 1] == X[i]) {
453:       *flg = PETSC_TRUE;
454:       break;
455:     }
456:   }
457:   PetscFunctionReturn(PETSC_SUCCESS);
458: }

460: /*@
461:   PetscSortRemoveDupsInt - Sorts an array of `PetscInt` in place in increasing order removes all duplicate entries

463:   Not Collective

465:   Input Parameters:
466: + n - number of values
467: - X - array of integers

469:   Output Parameter:
470: . n - number of non-redundant values

472:   Level: intermediate

474: .seealso: `PetscIntSortSemiOrdered()`, `PetscSortReal()`, `PetscSortIntWithPermutation()`, `PetscSortInt()`, `PetscSortedRemoveDupsInt()`
475: @*/
476: PetscErrorCode PetscSortRemoveDupsInt(PetscInt *n, PetscInt X[])
477: {
478:   PetscFunctionBegin;
479:   PetscAssertPointer(n, 1);
480:   PetscCall(PetscSortInt(*n, X));
481:   PetscCall(PetscSortedRemoveDupsInt(n, X));
482:   PetscFunctionReturn(PETSC_SUCCESS);
483: }

485: /*@
486:   PetscFindInt - Finds `PetscInt` in a sorted array of `PetscInt`

488:   Not Collective

490:   Input Parameters:
491: + key - the integer to locate
492: . n   - number of values in the array
493: - X   - array of integers

495:   Output Parameter:
496: . loc - the location if found, otherwise -(slot+1) where slot is the place the value would go

498:   Level: intermediate

500: .seealso: `PetscIntSortSemiOrdered()`, `PetscSortInt()`, `PetscSortIntWithArray()`, `PetscSortRemoveDupsInt()`
501: @*/
502: PetscErrorCode PetscFindInt(PetscInt key, PetscInt n, const PetscInt X[], PetscInt *loc)
503: {
504:   PetscInt lo = 0, hi = n;

506:   PetscFunctionBegin;
507:   PetscAssertPointer(loc, 4);
508:   if (!n) {
509:     *loc = -1;
510:     PetscFunctionReturn(PETSC_SUCCESS);
511:   }
512:   PetscAssertPointer(X, 3);
513:   PetscCheckSorted(n, X);
514:   while (hi - lo > 1) {
515:     PetscInt mid = lo + (hi - lo) / 2;
516:     if (key < X[mid]) hi = mid;
517:     else lo = mid;
518:   }
519:   *loc = key == X[lo] ? lo : -(lo + (key > X[lo]) + 1);
520:   PetscFunctionReturn(PETSC_SUCCESS);
521: }

523: /*@
524:   PetscCheckDupsInt - Checks if an `PetscInt` array has duplicates

526:   Not Collective

528:   Input Parameters:
529: + n - number of values in the array
530: - X - array of integers

532:   Output Parameter:
533: . dups - True if the array has dups, otherwise false

535:   Level: intermediate

537: .seealso: `PetscSortRemoveDupsInt()`, `PetscSortedCheckDupsInt()`
538: @*/
539: PetscErrorCode PetscCheckDupsInt(PetscInt n, const PetscInt X[], PetscBool *dups)
540: {
541:   PetscInt   i;
542:   PetscHSetI ht;
543:   PetscBool  missing;

545:   PetscFunctionBegin;
546:   if (n) PetscAssertPointer(X, 2);
547:   PetscAssertPointer(dups, 3);
548:   *dups = PETSC_FALSE;
549:   if (n > 1) {
550:     PetscCall(PetscHSetICreate(&ht));
551:     PetscCall(PetscHSetIResize(ht, n));
552:     for (i = 0; i < n; i++) {
553:       PetscCall(PetscHSetIQueryAdd(ht, X[i], &missing));
554:       if (!missing) {
555:         *dups = PETSC_TRUE;
556:         break;
557:       }
558:     }
559:     PetscCall(PetscHSetIDestroy(&ht));
560:   }
561:   PetscFunctionReturn(PETSC_SUCCESS);
562: }

564: /*@
565:   PetscFindMPIInt - Finds `PetscMPIInt` in a sorted array of `PetscMPIInt`

567:   Not Collective

569:   Input Parameters:
570: + key - the integer to locate
571: . n   - number of values in the array
572: - X   - array of integers

574:   Output Parameter:
575: . loc - the location if found, otherwise -(slot+1) where slot is the place the value would go

577:   Level: intermediate

579: .seealso: `PetscMPIIntSortSemiOrdered()`, `PetscSortInt()`, `PetscSortIntWithArray()`, `PetscSortRemoveDupsInt()`
580: @*/
581: PetscErrorCode PetscFindMPIInt(PetscMPIInt key, PetscInt n, const PetscMPIInt X[], PetscInt *loc)
582: {
583:   PetscInt lo = 0, hi = n;

585:   PetscFunctionBegin;
586:   PetscAssertPointer(loc, 4);
587:   if (!n) {
588:     *loc = -1;
589:     PetscFunctionReturn(PETSC_SUCCESS);
590:   }
591:   PetscAssertPointer(X, 3);
592:   PetscCheckSorted(n, X);
593:   while (hi - lo > 1) {
594:     PetscInt mid = lo + (hi - lo) / 2;
595:     if (key < X[mid]) hi = mid;
596:     else lo = mid;
597:   }
598:   *loc = key == X[lo] ? lo : -(lo + (key > X[lo]) + 1);
599:   PetscFunctionReturn(PETSC_SUCCESS);
600: }

602: /*@
603:   PetscSortIntWithArray - Sorts an array of `PetscInt` in place in increasing order;
604:   changes a second array of `PetscInt` to match the sorted first array.

606:   Not Collective

608:   Input Parameters:
609: + n - number of values
610: . X - array of integers
611: - Y - second array of integers

613:   Level: intermediate

615: .seealso: `PetscIntSortSemiOrderedWithArray()`, `PetscSortReal()`, `PetscSortIntWithPermutation()`, `PetscSortInt()`, `PetscSortIntWithCountArray()`
616: @*/
617: PetscErrorCode PetscSortIntWithArray(PetscInt n, PetscInt X[], PetscInt Y[])
618: {
619:   PetscInt pivot, t1, t2;

621:   PetscFunctionBegin;
622:   QuickSort2(PetscSortIntWithArray, X, Y, n, pivot, t1, t2);
623:   PetscFunctionReturn(PETSC_SUCCESS);
624: }

626: /*@
627:   PetscSortIntWithArrayPair - Sorts an array of `PetscInt` in place in increasing order;
628:   changes a pair of `PetscInt` arrays to match the sorted first array.

630:   Not Collective

632:   Input Parameters:
633: + n - number of values
634: . X - array of integers
635: . Y - second array of integers (first array of the pair)
636: - Z - third array of integers  (second array of the pair)

638:   Level: intermediate

640: .seealso: `PetscSortReal()`, `PetscSortIntWithPermutation()`, `PetscSortIntWithArray()`, `PetscIntSortSemiOrdered()`, `PetscSortIntWithIntCountArrayPair()`
641: @*/
642: PetscErrorCode PetscSortIntWithArrayPair(PetscInt n, PetscInt X[], PetscInt Y[], PetscInt Z[])
643: {
644:   PetscInt pivot, t1, t2, t3;

646:   PetscFunctionBegin;
647:   QuickSort3(PetscSortIntWithArrayPair, X, Y, Z, n, pivot, t1, t2, t3);
648:   PetscFunctionReturn(PETSC_SUCCESS);
649: }

651: /*@
652:   PetscSortIntWithCountArray - Sorts an array of `PetscInt` in place in increasing order;
653:   changes a second array of `PetscCount` to match the sorted first array.

655:   Not Collective

657:   Input Parameters:
658: + n - number of values
659: . X - array of integers
660: - Y - second array of PetscCounts (signed integers)

662:   Level: intermediate

664: .seealso: `PetscIntSortSemiOrderedWithArray()`, `PetscSortReal()`, `PetscSortIntPermutation()`, `PetscSortInt()`, `PetscSortIntWithArray()`
665: @*/
666: PetscErrorCode PetscSortIntWithCountArray(PetscCount n, PetscInt X[], PetscCount Y[])
667: {
668:   PetscInt   pivot, t1;
669:   PetscCount t2;

671:   PetscFunctionBegin;
672:   QuickSort2(PetscSortIntWithCountArray, X, Y, n, pivot, t1, t2);
673:   PetscFunctionReturn(PETSC_SUCCESS);
674: }

676: /*@
677:   PetscSortIntWithIntCountArrayPair - Sorts an array of `PetscInt` in place in increasing order;
678:   changes a `PetscInt`  array and a `PetscCount` array to match the sorted first array.

680:   Not Collective

682:   Input Parameters:
683: + n - number of values
684: . X - array of integers
685: . Y - second array of integers (first array of the pair)
686: - Z - third array of PetscCounts  (second array of the pair)

688:   Level: intermediate

690:   Note:
691:   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.

693: .seealso: `PetscSortReal()`, `PetscSortIntPermutation()`, `PetscSortIntWithArray()`, `PetscIntSortSemiOrdered()`, `PetscSortIntWithArrayPair()`
694: @*/
695: PetscErrorCode PetscSortIntWithIntCountArrayPair(PetscCount n, PetscInt X[], PetscInt Y[], PetscCount Z[])
696: {
697:   PetscInt   pivot, t1, t2; /* pivot is take from X[], so its type is still PetscInt */
698:   PetscCount t3;            /* temp for Z[] */

700:   PetscFunctionBegin;
701:   QuickSort3(PetscSortIntWithIntCountArrayPair, X, Y, Z, n, pivot, t1, t2, t3);
702:   PetscFunctionReturn(PETSC_SUCCESS);
703: }

705: /*@
706:   PetscSortedMPIInt - Determines whether the `PetscMPIInt` array is sorted.

708:   Not Collective

710:   Input Parameters:
711: + n - number of values
712: - X - array of integers

714:   Output Parameter:
715: . sorted - flag whether the array is sorted

717:   Level: intermediate

719: .seealso: `PetscMPIIntSortSemiOrdered()`, `PetscSortMPIInt()`, `PetscSortedInt()`, `PetscSortedReal()`
720: @*/
721: PetscErrorCode PetscSortedMPIInt(PetscInt n, const PetscMPIInt X[], PetscBool *sorted)
722: {
723:   PetscFunctionBegin;
724:   PetscSorted(n, X, *sorted);
725:   PetscFunctionReturn(PETSC_SUCCESS);
726: }

728: /*@
729:   PetscSortMPIInt - Sorts an array of `PetscMPIInt` in place in increasing order.

731:   Not Collective

733:   Input Parameters:
734: + n - number of values
735: - X - array of integers

737:   Level: intermediate

739:   Note:
740:   This function serves as an alternative to PetscMPIIntSortSemiOrdered(), and may perform faster especially if the array
741:   is completely random. There are exceptions to this and so it is __highly__ recommended that the user benchmark their
742:   code to see which routine is fastest.

744: .seealso: `PetscMPIIntSortSemiOrdered()`, `PetscSortReal()`, `PetscSortIntWithPermutation()`
745: @*/
746: PetscErrorCode PetscSortMPIInt(PetscInt n, PetscMPIInt X[])
747: {
748:   PetscMPIInt pivot, t1;

750:   PetscFunctionBegin;
751:   QuickSort1(PetscSortMPIInt, X, n, pivot, t1);
752:   PetscFunctionReturn(PETSC_SUCCESS);
753: }

755: /*@
756:   PetscSortRemoveDupsMPIInt - Sorts an array of `PetscMPIInt` in place in increasing order removes all duplicate entries

758:   Not Collective

760:   Input Parameters:
761: + n - number of values
762: - X - array of integers

764:   Output Parameter:
765: . n - number of non-redundant values

767:   Level: intermediate

769: .seealso: `PetscSortReal()`, `PetscSortIntWithPermutation()`, `PetscSortInt()`
770: @*/
771: PetscErrorCode PetscSortRemoveDupsMPIInt(PetscInt *n, PetscMPIInt X[])
772: {
773:   PetscInt s = 0, N = *n, b = 0;

775:   PetscFunctionBegin;
776:   PetscCall(PetscSortMPIInt(N, X));
777:   for (PetscInt i = 0; i < N - 1; i++) {
778:     if (X[b + s + 1] != X[b]) {
779:       X[b + 1] = X[b + s + 1];
780:       b++;
781:     } else s++;
782:   }
783:   *n = N - s;
784:   PetscFunctionReturn(PETSC_SUCCESS);
785: }

787: /*@
788:   PetscSortMPIIntWithArray - Sorts an array of `PetscMPIInt` in place in increasing order;
789:   changes a second `PetscMPIInt` array to match the sorted first array.

791:   Not Collective

793:   Input Parameters:
794: + n - number of values
795: . X - array of integers
796: - Y - second array of integers

798:   Level: intermediate

800: .seealso: `PetscMPIIntSortSemiOrderedWithArray()`, `PetscSortReal()`, `PetscSortIntWithPermutation()`, `PetscSortInt()`
801: @*/
802: PetscErrorCode PetscSortMPIIntWithArray(PetscMPIInt n, PetscMPIInt X[], PetscMPIInt Y[])
803: {
804:   PetscMPIInt pivot, t1, t2;

806:   PetscFunctionBegin;
807:   QuickSort2(PetscSortMPIIntWithArray, X, Y, n, pivot, t1, t2);
808:   PetscFunctionReturn(PETSC_SUCCESS);
809: }

811: /*@
812:   PetscSortMPIIntWithIntArray - Sorts an array of `PetscMPIInt` in place in increasing order;
813:   changes a second array of `PetscInt` to match the sorted first array.

815:   Not Collective

817:   Input Parameters:
818: + n - number of values
819: . X - array of MPI integers
820: - Y - second array of Petsc integers

822:   Level: intermediate

824:   Note:
825:   This routine is useful when one needs to sort MPI ranks with other integer arrays.

827: .seealso: `PetscSortMPIIntWithArray()`, `PetscIntSortSemiOrderedWithArray()`, `PetscTimSortWithArray()`
828: @*/
829: PetscErrorCode PetscSortMPIIntWithIntArray(PetscMPIInt n, PetscMPIInt X[], PetscInt Y[])
830: {
831:   PetscMPIInt pivot, t1;
832:   PetscInt    t2;

834:   PetscFunctionBegin;
835:   QuickSort2(PetscSortMPIIntWithIntArray, X, Y, n, pivot, t1, t2);
836:   PetscFunctionReturn(PETSC_SUCCESS);
837: }

839: /*@
840:   PetscSortIntWithScalarArray - Sorts an array of `PetscInt` in place in increasing order;
841:   changes a second `PetscScalar` array to match the sorted first array.

843:   Not Collective

845:   Input Parameters:
846: + n - number of values
847: . X - array of integers
848: - Y - second array of scalars

850:   Level: intermediate

852: .seealso: `PetscTimSortWithArray()`, `PetscSortReal()`, `PetscSortIntWithPermutation()`, `PetscSortInt()`, `PetscSortIntWithArray()`
853: @*/
854: PetscErrorCode PetscSortIntWithScalarArray(PetscInt n, PetscInt X[], PetscScalar Y[])
855: {
856:   PetscInt    pivot, t1;
857:   PetscScalar t2;

859:   PetscFunctionBegin;
860:   QuickSort2(PetscSortIntWithScalarArray, X, Y, n, pivot, t1, t2);
861:   PetscFunctionReturn(PETSC_SUCCESS);
862: }

864: /*@C
865:   PetscSortIntWithDataArray - Sorts an array of `PetscInt` in place in increasing order;
866:   changes a second array to match the sorted first INTEGER array.  Unlike other sort routines, the user must
867:   provide workspace (the size of an element in the data array) to use when sorting.

869:   Not Collective

871:   Input Parameters:
872: + n    - number of values
873: . X    - array of integers
874: . Y    - second array of data
875: . size - sizeof elements in the data array in bytes
876: - t2   - workspace of "size" bytes used when sorting

878:   Level: intermediate

880: .seealso: `PetscTimSortWithArray()`, `PetscSortReal()`, `PetscSortIntWithPermutation()`, `PetscSortInt()`, `PetscSortIntWithArray()`
881: @*/
882: PetscErrorCode PetscSortIntWithDataArray(PetscInt n, PetscInt X[], void *Y, size_t size, void *t2)
883: {
884:   char    *YY = (char *)Y;
885:   PetscInt t1, pivot, hi = n - 1;

887:   PetscFunctionBegin;
888:   if (n < 8) {
889:     for (PetscInt i = 0; i < n; i++) {
890:       pivot = X[i];
891:       for (PetscInt j = i + 1; j < n; j++) {
892:         if (pivot > X[j]) {
893:           SWAP2Data(X[i], X[j], YY + size * i, YY + size * j, t1, t2, size);
894:           pivot = X[i];
895:         }
896:       }
897:     }
898:   } else {
899:     /* Two way partition */
900:     PetscInt l = 0, r = hi;

902:     pivot = X[MEDIAN(X, hi)];
903:     while (1) {
904:       while (X[l] < pivot) l++;
905:       while (X[r] > pivot) r--;
906:       if (l >= r) {
907:         r++;
908:         break;
909:       }
910:       SWAP2Data(X[l], X[r], YY + size * l, YY + size * r, t1, t2, size);
911:       l++;
912:       r--;
913:     }
914:     PetscCall(PetscSortIntWithDataArray(l, X, Y, size, t2));
915:     PetscCall(PetscSortIntWithDataArray(hi - r + 1, X + r, YY + size * r, size, t2));
916:   }
917:   PetscFunctionReturn(PETSC_SUCCESS);
918: }

920: /*@
921:   PetscMergeIntArray -     Merges two SORTED `PetscInt` arrays, removes duplicate elements.

923:   Not Collective

925:   Input Parameters:
926: + an - number of values in the first array
927: . aI - first sorted array of integers
928: . bn - number of values in the second array
929: - bI - second array of integers

931:   Output Parameters:
932: + n - number of values in the merged array
933: - L - merged sorted array, this is allocated if an array is not provided

935:   Level: intermediate

937: .seealso: `PetscSortReal()`, `PetscSortIntWithPermutation()`, `PetscSortInt()`, `PetscSortIntWithArray()`
938: @*/
939: PetscErrorCode PetscMergeIntArray(PetscInt an, const PetscInt aI[], PetscInt bn, const PetscInt bI[], PetscInt *n, PetscInt **L)
940: {
941:   PetscInt *L_ = *L, ak, bk, k;

943:   PetscFunctionBegin;
944:   if (!L_) {
945:     PetscCall(PetscMalloc1(an + bn, L));
946:     L_ = *L;
947:   }
948:   k = ak = bk = 0;
949:   while (ak < an && bk < bn) {
950:     if (aI[ak] == bI[bk]) {
951:       L_[k] = aI[ak];
952:       ++ak;
953:       ++bk;
954:       ++k;
955:     } else if (aI[ak] < bI[bk]) {
956:       L_[k] = aI[ak];
957:       ++ak;
958:       ++k;
959:     } else {
960:       L_[k] = bI[bk];
961:       ++bk;
962:       ++k;
963:     }
964:   }
965:   if (ak < an) {
966:     PetscCall(PetscArraycpy(L_ + k, aI + ak, an - ak));
967:     k += (an - ak);
968:   }
969:   if (bk < bn) {
970:     PetscCall(PetscArraycpy(L_ + k, bI + bk, bn - bk));
971:     k += (bn - bk);
972:   }
973:   *n = k;
974:   PetscFunctionReturn(PETSC_SUCCESS);
975: }

977: /*@
978:   PetscMergeIntArrayPair -     Merges two SORTED `PetscInt` arrays that share NO common values along with an additional array of `PetscInt`.
979:   The additional arrays are the same length as sorted arrays and are merged
980:   in the order determined by the merging of the sorted pair.

982:   Not Collective

984:   Input Parameters:
985: + an - number of values in the first array
986: . aI - first sorted array of integers
987: . aJ - first additional array of integers
988: . bn - number of values in the second array
989: . bI - second array of integers
990: - bJ - second additional of integers

992:   Output Parameters:
993: + n - number of values in the merged array (== an + bn)
994: . L - merged sorted array
995: - J - merged additional array

997:   Note:
998:   if L or J point to non-null arrays then this routine will assume they are of the appropriate size and use them, otherwise this routine will allocate space for them

1000:   Level: intermediate

1002: .seealso: `PetscIntSortSemiOrdered()`, `PetscSortReal()`, `PetscSortIntWithPermutation()`, `PetscSortInt()`, `PetscSortIntWithArray()`
1003: @*/
1004: PetscErrorCode PetscMergeIntArrayPair(PetscInt an, const PetscInt aI[], const PetscInt aJ[], PetscInt bn, const PetscInt bI[], const PetscInt bJ[], PetscInt *n, PetscInt **L, PetscInt **J)
1005: {
1006:   PetscInt n_, *L_, *J_, ak, bk, k;

1008:   PetscFunctionBegin;
1009:   PetscAssertPointer(L, 8);
1010:   PetscAssertPointer(J, 9);
1011:   n_ = an + bn;
1012:   *n = n_;
1013:   if (!*L) PetscCall(PetscMalloc1(n_, L));
1014:   L_ = *L;
1015:   if (!*J) PetscCall(PetscMalloc1(n_, J));
1016:   J_ = *J;
1017:   k = ak = bk = 0;
1018:   while (ak < an && bk < bn) {
1019:     if (aI[ak] <= bI[bk]) {
1020:       L_[k] = aI[ak];
1021:       J_[k] = aJ[ak];
1022:       ++ak;
1023:       ++k;
1024:     } else {
1025:       L_[k] = bI[bk];
1026:       J_[k] = bJ[bk];
1027:       ++bk;
1028:       ++k;
1029:     }
1030:   }
1031:   if (ak < an) {
1032:     PetscCall(PetscArraycpy(L_ + k, aI + ak, an - ak));
1033:     PetscCall(PetscArraycpy(J_ + k, aJ + ak, an - ak));
1034:     k += (an - ak);
1035:   }
1036:   if (bk < bn) {
1037:     PetscCall(PetscArraycpy(L_ + k, bI + bk, bn - bk));
1038:     PetscCall(PetscArraycpy(J_ + k, bJ + bk, bn - bk));
1039:   }
1040:   PetscFunctionReturn(PETSC_SUCCESS);
1041: }

1043: /*@
1044:   PetscMergeMPIIntArray -     Merges two SORTED `PetscMPIInt` arrays.

1046:   Not Collective

1048:   Input Parameters:
1049: + an - number of values in the first array
1050: . aI - first sorted array of integers
1051: . bn - number of values in the second array
1052: - bI - second array of integers

1054:   Output Parameters:
1055: + n - number of values in the merged array (<= an + bn)
1056: - L - merged sorted array, allocated if address of NULL pointer is passed

1058:   Level: intermediate

1060: .seealso: `PetscIntSortSemiOrdered()`, `PetscSortReal()`, `PetscSortIntWithPermutation()`, `PetscSortInt()`, `PetscSortIntWithArray()`
1061: @*/
1062: PetscErrorCode PetscMergeMPIIntArray(PetscInt an, const PetscMPIInt aI[], PetscInt bn, const PetscMPIInt bI[], PetscInt *n, PetscMPIInt **L)
1063: {
1064:   PetscInt ai, bi, k;

1066:   PetscFunctionBegin;
1067:   if (!*L) PetscCall(PetscMalloc1((an + bn), L));
1068:   for (ai = 0, bi = 0, k = 0; ai < an || bi < bn;) {
1069:     PetscInt t = -1;
1070:     for (; ai < an && (!bn || aI[ai] <= bI[bi]); ai++) (*L)[k++] = t = aI[ai];
1071:     for (; bi < bn && bI[bi] == t; bi++)
1072:       ;
1073:     for (; bi < bn && (!an || bI[bi] <= aI[ai]); bi++) (*L)[k++] = t = bI[bi];
1074:     for (; ai < an && aI[ai] == t; ai++)
1075:       ;
1076:   }
1077:   *n = k;
1078:   PetscFunctionReturn(PETSC_SUCCESS);
1079: }

1081: /*@C
1082:   PetscProcessTree - Prepares tree data to be displayed graphically

1084:   Not Collective

1086:   Input Parameters:
1087: + n        - number of values
1088: . mask     - indicates those entries in the tree, location 0 is always masked
1089: - parentid - indicates the parent of each entry

1091:   Output Parameters:
1092: + Nlevels   - the number of levels
1093: . Level     - for each node tells its level
1094: . Levelcnt  - the number of nodes on each level
1095: . Idbylevel - a list of ids on each of the levels, first level followed by second etc
1096: - Column    - for each id tells its column index

1098:   Level: developer

1100:   Note:
1101:   This code is not currently used

1103: .seealso: `PetscSortReal()`, `PetscSortIntWithPermutation()`
1104: @*/
1105: PetscErrorCode PetscProcessTree(PetscInt n, const PetscBool mask[], const PetscInt parentid[], PetscInt *Nlevels, PetscInt **Level, PetscInt **Levelcnt, PetscInt **Idbylevel, PetscInt **Column)
1106: {
1107:   PetscInt  i, j, cnt, nmask = 0, nlevels = 0, *level, *levelcnt, levelmax = 0, *workid, *workparentid, tcnt = 0, *idbylevel, *column;
1108:   PetscBool done = PETSC_FALSE;

1110:   PetscFunctionBegin;
1111:   PetscCheck(mask[0], PETSC_COMM_SELF, PETSC_ERR_SUP, "Mask of 0th location must be set");
1112:   for (i = 0; i < n; i++) {
1113:     if (mask[i]) continue;
1114:     PetscCheck(parentid[i] != i, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Node labeled as own parent");
1115:     PetscCheck(!parentid[i] || !mask[parentid[i]], PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Parent is masked");
1116:   }

1118:   for (i = 0; i < n; i++) {
1119:     if (!mask[i]) nmask++;
1120:   }

1122:   /* determine the level in the tree of each node */
1123:   PetscCall(PetscCalloc1(n, &level));

1125:   level[0] = 1;
1126:   while (!done) {
1127:     done = PETSC_TRUE;
1128:     for (i = 0; i < n; i++) {
1129:       if (mask[i]) continue;
1130:       if (!level[i] && level[parentid[i]]) level[i] = level[parentid[i]] + 1;
1131:       else if (!level[i]) done = PETSC_FALSE;
1132:     }
1133:   }
1134:   for (i = 0; i < n; i++) {
1135:     level[i]--;
1136:     nlevels = PetscMax(nlevels, level[i]);
1137:   }

1139:   /* count the number of nodes on each level and its max */
1140:   PetscCall(PetscCalloc1(nlevels, &levelcnt));
1141:   for (i = 0; i < n; i++) {
1142:     if (mask[i]) continue;
1143:     levelcnt[level[i] - 1]++;
1144:   }
1145:   for (i = 0; i < nlevels; i++) levelmax = PetscMax(levelmax, levelcnt[i]);

1147:   /* for each level sort the ids by the parent id */
1148:   PetscCall(PetscMalloc2(levelmax, &workid, levelmax, &workparentid));
1149:   PetscCall(PetscMalloc1(nmask, &idbylevel));
1150:   for (j = 1; j <= nlevels; j++) {
1151:     cnt = 0;
1152:     for (i = 0; i < n; i++) {
1153:       if (mask[i]) continue;
1154:       if (level[i] != j) continue;
1155:       workid[cnt]         = i;
1156:       workparentid[cnt++] = parentid[i];
1157:     }
1158:     /*  PetscIntView(cnt,workparentid,0);
1159:     PetscIntView(cnt,workid,0);
1160:     PetscCall(PetscSortIntWithArray(cnt,workparentid,workid));
1161:     PetscIntView(cnt,workparentid,0);
1162:     PetscIntView(cnt,workid,0);*/
1163:     PetscCall(PetscArraycpy(idbylevel + tcnt, workid, cnt));
1164:     tcnt += cnt;
1165:   }
1166:   PetscCheck(tcnt == nmask, PETSC_COMM_SELF, PETSC_ERR_PLIB, "Inconsistent count of unmasked nodes");
1167:   PetscCall(PetscFree2(workid, workparentid));

1169:   /* for each node list its column */
1170:   PetscCall(PetscMalloc1(n, &column));
1171:   cnt = 0;
1172:   for (j = 0; j < nlevels; j++) {
1173:     for (i = 0; i < levelcnt[j]; i++) column[idbylevel[cnt++]] = i;
1174:   }

1176:   *Nlevels   = nlevels;
1177:   *Level     = level;
1178:   *Levelcnt  = levelcnt;
1179:   *Idbylevel = idbylevel;
1180:   *Column    = column;
1181:   PetscFunctionReturn(PETSC_SUCCESS);
1182: }

1184: /*@
1185:   PetscParallelSortedInt - Check whether a `PetscInt` array, distributed over a communicator, is globally sorted.

1187:   Collective

1189:   Input Parameters:
1190: + comm - the MPI communicator
1191: . n    - the local number of integers
1192: - keys - the local array of integers

1194:   Output Parameters:
1195: . is_sorted - whether the array is globally sorted

1197:   Level: developer

1199: .seealso: `PetscParallelSortInt()`
1200: @*/
1201: PetscErrorCode PetscParallelSortedInt(MPI_Comm comm, PetscInt n, const PetscInt keys[], PetscBool *is_sorted)
1202: {
1203:   PetscBool   sorted;
1204:   PetscInt    i, min, max, prevmax;
1205:   PetscMPIInt rank;

1207:   PetscFunctionBegin;
1208:   sorted = PETSC_TRUE;
1209:   min    = PETSC_MAX_INT;
1210:   max    = PETSC_MIN_INT;
1211:   if (n) {
1212:     min = keys[0];
1213:     max = keys[0];
1214:   }
1215:   for (i = 1; i < n; i++) {
1216:     if (keys[i] < keys[i - 1]) break;
1217:     min = PetscMin(min, keys[i]);
1218:     max = PetscMax(max, keys[i]);
1219:   }
1220:   if (i < n) sorted = PETSC_FALSE;
1221:   prevmax = PETSC_MIN_INT;
1222:   PetscCallMPI(MPI_Exscan(&max, &prevmax, 1, MPIU_INT, MPI_MAX, comm));
1223:   PetscCallMPI(MPI_Comm_rank(comm, &rank));
1224:   if (rank == 0) prevmax = PETSC_MIN_INT;
1225:   if (prevmax > min) sorted = PETSC_FALSE;
1226:   PetscCall(MPIU_Allreduce(&sorted, is_sorted, 1, MPIU_BOOL, MPI_LAND, comm));
1227:   PetscFunctionReturn(PETSC_SUCCESS);
1228: }