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: }