Actual source code: sorti.c

petsc-3.8.3 2017-12-09
Report Typos and Errors

  2: /*
  3:    This file contains routines for sorting integers. Values are sorted in place.
  4:  */
  5:  #include <petsc/private/petscimpl.h>

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

  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: /* -----------------------------------------------------------------------*/

 22: /*
 23:    A simple version of quicksort; taken from Kernighan and Ritchie, page 87.
 24:    Assumes 0 origin for v, number of elements = right+1 (right is index of
 25:    right-most member).
 26: */
 27: static void PetscSortInt_Private(PetscInt *v,PetscInt right)
 28: {
 29:   PetscInt i,j,pivot,tmp;

 31:   if (right <= 1) {
 32:     if (right == 1) {
 33:       if (v[0] > v[1]) SWAP(v[0],v[1],tmp);
 34:     }
 35:     return;
 36:   }
 37:   i = MEDIAN(v,right);          /* Choose a pivot */
 38:   SWAP(v[0],v[i],tmp);          /* Move it out of the way */
 39:   pivot = v[0];
 40:   for (i=0,j=right+1;; ) {
 41:     while (++i < j && v[i] <= pivot) ; /* Scan from the left */
 42:     while (v[--j] > pivot) ;           /* Scan from the right */
 43:     if (i >= j) break;
 44:     SWAP(v[i],v[j],tmp);
 45:   }
 46:   SWAP(v[0],v[j],tmp);          /* Put pivot back in place. */
 47:   PetscSortInt_Private(v,j-1);
 48:   PetscSortInt_Private(v+j+1,right-(j+1));
 49: }

 51: /*@
 52:    PetscSortInt - Sorts an array of integers in place in increasing order.

 54:    Not Collective

 56:    Input Parameters:
 57: +  n  - number of values
 58: -  i  - array of integers

 60:    Level: intermediate

 62:    Concepts: sorting^ints

 64: .seealso: PetscSortReal(), PetscSortIntWithPermutation()
 65: @*/
 66: PetscErrorCode  PetscSortInt(PetscInt n,PetscInt i[])
 67: {
 68:   PetscInt j,k,tmp,ik;

 71:   if (n<8) {
 72:     for (k=0; k<n; k++) {
 73:       ik = i[k];
 74:       for (j=k+1; j<n; j++) {
 75:         if (ik > i[j]) {
 76:           SWAP(i[k],i[j],tmp);
 77:           ik = i[k];
 78:         }
 79:       }
 80:     }
 81:   } else PetscSortInt_Private(i,n-1);
 82:   return(0);
 83: }

 85: /*@
 86:    PetscSortedRemoveDupsInt - Removes all duplicate entries of a sorted input array

 88:    Not Collective

 90:    Input Parameters:
 91: +  n  - number of values
 92: -  ii  - sorted array of integers

 94:    Output Parameter:
 95: .  n - number of non-redundant values

 97:    Level: intermediate

 99:    Concepts: sorting^ints

101: .seealso: PetscSortInt()
102: @*/
103: PetscErrorCode  PetscSortedRemoveDupsInt(PetscInt *n,PetscInt ii[])
104: {
105:   PetscInt i,s = 0,N = *n, b = 0;

108:   for (i=0; i<N-1; i++) {
109:     if (PetscUnlikely(ii[b+s+1] < ii[b])) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"Input array is not sorted");
110:     if (ii[b+s+1] != ii[b]) {
111:       ii[b+1] = ii[b+s+1]; b++;
112:     } else s++;
113:   }
114:   *n = N - s;
115:   return(0);
116: }

118: /*@
119:    PetscSortRemoveDupsInt - Sorts an array of integers in place in increasing order removes all duplicate entries

121:    Not Collective

123:    Input Parameters:
124: +  n  - number of values
125: -  ii  - array of integers

127:    Output Parameter:
128: .  n - number of non-redundant values

130:    Level: intermediate

132:    Concepts: sorting^ints

134: .seealso: PetscSortReal(), PetscSortIntWithPermutation(), PetscSortInt(), PetscSortedRemoveDupsInt()
135: @*/
136: PetscErrorCode  PetscSortRemoveDupsInt(PetscInt *n,PetscInt ii[])
137: {

141:   PetscSortInt(*n,ii);
142:   PetscSortedRemoveDupsInt(n,ii);
143:   return(0);
144: }

146: /*@
147:   PetscFindInt - Finds integer in a sorted array of integers

149:    Not Collective

151:    Input Parameters:
152: +  key - the integer to locate
153: .  n   - number of values in the array
154: -  ii  - array of integers

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

159:    Level: intermediate

161:    Concepts: sorting^ints

163: .seealso: PetscSortInt(), PetscSortIntWithArray(), PetscSortRemoveDupsInt()
164: @*/
165: PetscErrorCode PetscFindInt(PetscInt key, PetscInt n, const PetscInt ii[], PetscInt *loc)
166: {
167:   PetscInt lo = 0,hi = n;

171:   if (!n) {*loc = -1; return(0);}
173:   while (hi - lo > 1) {
174:     PetscInt mid = lo + (hi - lo)/2;
175:     if (key < ii[mid]) hi = mid;
176:     else               lo = mid;
177:   }
178:   *loc = key == ii[lo] ? lo : -(lo + (key > ii[lo]) + 1);
179:   return(0);
180: }

182: /*@
183:   PetscFindMPIInt - Finds MPI integer in a sorted array of integers

185:    Not Collective

187:    Input Parameters:
188: +  key - the integer to locate
189: .  n   - number of values in the array
190: -  ii  - array of integers

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

195:    Level: intermediate

197:    Concepts: sorting^ints

199: .seealso: PetscSortInt(), PetscSortIntWithArray(), PetscSortRemoveDupsInt()
200: @*/
201: PetscErrorCode PetscFindMPIInt(PetscMPIInt key, PetscInt n, const PetscMPIInt ii[], PetscInt *loc)
202: {
203:   PetscInt lo = 0,hi = n;

207:   if (!n) {*loc = -1; return(0);}
209:   while (hi - lo > 1) {
210:     PetscInt mid = lo + (hi - lo)/2;
211:     if (key < ii[mid]) hi = mid;
212:     else               lo = mid;
213:   }
214:   *loc = key == ii[lo] ? lo : -(lo + (key > ii[lo]) + 1);
215:   return(0);
216: }

218: /* -----------------------------------------------------------------------*/
219: #define SWAP2(a,b,c,d,t) {t=a;a=b;b=t;t=c;c=d;d=t;}

221: /*
222:    A simple version of quicksort; taken from Kernighan and Ritchie, page 87.
223:    Assumes 0 origin for v, number of elements = right+1 (right is index of
224:    right-most member).
225: */
226: static PetscErrorCode PetscSortIntWithArray_Private(PetscInt *v,PetscInt *V,PetscInt right)
227: {
229:   PetscInt       i,vl,last,tmp;

232:   if (right <= 1) {
233:     if (right == 1) {
234:       if (v[0] > v[1]) SWAP2(v[0],v[1],V[0],V[1],tmp);
235:     }
236:     return(0);
237:   }
238:   SWAP2(v[0],v[right/2],V[0],V[right/2],tmp);
239:   vl   = v[0];
240:   last = 0;
241:   for (i=1; i<=right; i++) {
242:     if (v[i] < vl) {last++; SWAP2(v[last],v[i],V[last],V[i],tmp);}
243:   }
244:   SWAP2(v[0],v[last],V[0],V[last],tmp);
245:   PetscSortIntWithArray_Private(v,V,last-1);
246:   PetscSortIntWithArray_Private(v+last+1,V+last+1,right-(last+1));
247:   return(0);
248: }

250: /*@
251:    PetscSortIntWithArray - Sorts an array of integers in place in increasing order;
252:        changes a second array to match the sorted first array.

254:    Not Collective

256:    Input Parameters:
257: +  n  - number of values
258: .  i  - array of integers
259: -  I - second array of integers

261:    Level: intermediate

263:    Concepts: sorting^ints with array

265: .seealso: PetscSortReal(), PetscSortIntPermutation(), PetscSortInt()
266: @*/
267: PetscErrorCode  PetscSortIntWithArray(PetscInt n,PetscInt i[],PetscInt Ii[])
268: {
270:   PetscInt       j,k,tmp,ik;

273:   if (n<8) {
274:     for (k=0; k<n; k++) {
275:       ik = i[k];
276:       for (j=k+1; j<n; j++) {
277:         if (ik > i[j]) {
278:           SWAP2(i[k],i[j],Ii[k],Ii[j],tmp);
279:           ik = i[k];
280:         }
281:       }
282:     }
283:   } else {
284:     PetscSortIntWithArray_Private(i,Ii,n-1);
285:   }
286:   return(0);
287: }


290: #define SWAP3(a,b,c,d,e,f,t) {t=a;a=b;b=t;t=c;c=d;d=t;t=e;e=f;f=t;}

292: /*
293:    A simple version of quicksort; taken from Kernighan and Ritchie, page 87.
294:    Assumes 0 origin for v, number of elements = right+1 (right is index of
295:    right-most member).
296: */
297: static PetscErrorCode PetscSortIntWithArrayPair_Private(PetscInt *L,PetscInt *J, PetscInt *K, PetscInt right)
298: {
300:   PetscInt       i,vl,last,tmp;

303:   if (right <= 1) {
304:     if (right == 1) {
305:       if (L[0] > L[1]) SWAP3(L[0],L[1],J[0],J[1],K[0],K[1],tmp);
306:     }
307:     return(0);
308:   }
309:   SWAP3(L[0],L[right/2],J[0],J[right/2],K[0],K[right/2],tmp);
310:   vl   = L[0];
311:   last = 0;
312:   for (i=1; i<=right; i++) {
313:     if (L[i] < vl) {last++; SWAP3(L[last],L[i],J[last],J[i],K[last],K[i],tmp);}
314:   }
315:   SWAP3(L[0],L[last],J[0],J[last],K[0],K[last],tmp);
316:   PetscSortIntWithArrayPair_Private(L,J,K,last-1);
317:   PetscSortIntWithArrayPair_Private(L+last+1,J+last+1,K+last+1,right-(last+1));
318:   return(0);
319: }

321: /*@
322:    PetscSortIntWithArrayPair - Sorts an array of integers in place in increasing order;
323:        changes a pair of integer arrays to match the sorted first array.

325:    Not Collective

327:    Input Parameters:
328: +  n  - number of values
329: .  I  - array of integers
330: .  J  - second array of integers (first array of the pair)
331: -  K  - third array of integers  (second array of the pair)

333:    Level: intermediate

335:    Concepts: sorting^ints with array pair

337: .seealso: PetscSortReal(), PetscSortIntPermutation(), PetscSortIntWithArray()
338: @*/
339: PetscErrorCode  PetscSortIntWithArrayPair(PetscInt n,PetscInt L[],PetscInt J[], PetscInt K[])
340: {
342:   PetscInt       j,k,tmp,ik;

345:   if (n<8) {
346:     for (k=0; k<n; k++) {
347:       ik = L[k];
348:       for (j=k+1; j<n; j++) {
349:         if (ik > L[j]) {
350:           SWAP3(L[k],L[j],J[k],J[j],K[k],K[j],tmp);
351:           ik = L[k];
352:         }
353:       }
354:     }
355:   } else {
356:     PetscSortIntWithArrayPair_Private(L,J,K,n-1);
357:   }
358:   return(0);
359: }

361: /*
362:    A simple version of quicksort; taken from Kernighan and Ritchie, page 87.
363:    Assumes 0 origin for v, number of elements = right+1 (right is index of
364:    right-most member).
365: */
366: static void PetscSortMPIInt_Private(PetscMPIInt *v,PetscInt right)
367: {
368:   PetscInt          i,j;
369:   PetscMPIInt       pivot,tmp;

371:   if (right <= 1) {
372:     if (right == 1) {
373:       if (v[0] > v[1]) SWAP(v[0],v[1],tmp);
374:     }
375:     return;
376:   }
377:   i = MEDIAN(v,right);          /* Choose a pivot */
378:   SWAP(v[0],v[i],tmp);          /* Move it out of the way */
379:   pivot = v[0];
380:   for (i=0,j=right+1;; ) {
381:     while (++i < j && v[i] <= pivot) ; /* Scan from the left */
382:     while (v[--j] > pivot) ;           /* Scan from the right */
383:     if (i >= j) break;
384:     SWAP(v[i],v[j],tmp);
385:   }
386:   SWAP(v[0],v[j],tmp);          /* Put pivot back in place. */
387:   PetscSortMPIInt_Private(v,j-1);
388:   PetscSortMPIInt_Private(v+j+1,right-(j+1));
389: }

391: /*@
392:    PetscSortMPIInt - Sorts an array of MPI integers in place in increasing order.

394:    Not Collective

396:    Input Parameters:
397: +  n  - number of values
398: -  i  - array of integers

400:    Level: intermediate

402:    Concepts: sorting^ints

404: .seealso: PetscSortReal(), PetscSortIntWithPermutation()
405: @*/
406: PetscErrorCode  PetscSortMPIInt(PetscInt n,PetscMPIInt i[])
407: {
408:   PetscInt    j,k;
409:   PetscMPIInt tmp,ik;

412:   if (n<8) {
413:     for (k=0; k<n; k++) {
414:       ik = i[k];
415:       for (j=k+1; j<n; j++) {
416:         if (ik > i[j]) {
417:           SWAP(i[k],i[j],tmp);
418:           ik = i[k];
419:         }
420:       }
421:     }
422:   } else PetscSortMPIInt_Private(i,n-1);
423:   return(0);
424: }

426: /*@
427:    PetscSortRemoveDupsMPIInt - Sorts an array of MPI integers in place in increasing order removes all duplicate entries

429:    Not Collective

431:    Input Parameters:
432: +  n  - number of values
433: -  ii  - array of integers

435:    Output Parameter:
436: .  n - number of non-redundant values

438:    Level: intermediate

440:    Concepts: sorting^ints

442: .seealso: PetscSortReal(), PetscSortIntWithPermutation(), PetscSortInt()
443: @*/
444: PetscErrorCode  PetscSortRemoveDupsMPIInt(PetscInt *n,PetscMPIInt ii[])
445: {
447:   PetscInt       i,s = 0,N = *n, b = 0;

450:   PetscSortMPIInt(N,ii);
451:   for (i=0; i<N-1; i++) {
452:     if (ii[b+s+1] != ii[b]) {
453:       ii[b+1] = ii[b+s+1]; b++;
454:     } else s++;
455:   }
456:   *n = N - s;
457:   return(0);
458: }

460: /*
461:    A simple version of quicksort; taken from Kernighan and Ritchie, page 87.
462:    Assumes 0 origin for v, number of elements = right+1 (right is index of
463:    right-most member).
464: */
465: static PetscErrorCode PetscSortMPIIntWithArray_Private(PetscMPIInt *v,PetscMPIInt *V,PetscMPIInt right)
466: {
468:   PetscMPIInt    i,vl,last,tmp;

471:   if (right <= 1) {
472:     if (right == 1) {
473:       if (v[0] > v[1]) SWAP2(v[0],v[1],V[0],V[1],tmp);
474:     }
475:     return(0);
476:   }
477:   SWAP2(v[0],v[right/2],V[0],V[right/2],tmp);
478:   vl   = v[0];
479:   last = 0;
480:   for (i=1; i<=right; i++) {
481:     if (v[i] < vl) {last++; SWAP2(v[last],v[i],V[last],V[i],tmp);}
482:   }
483:   SWAP2(v[0],v[last],V[0],V[last],tmp);
484:   PetscSortMPIIntWithArray_Private(v,V,last-1);
485:   PetscSortMPIIntWithArray_Private(v+last+1,V+last+1,right-(last+1));
486:   return(0);
487: }

489: /*@
490:    PetscSortMPIIntWithArray - Sorts an array of integers in place in increasing order;
491:        changes a second array to match the sorted first array.

493:    Not Collective

495:    Input Parameters:
496: +  n  - number of values
497: .  i  - array of integers
498: -  I - second array of integers

500:    Level: intermediate

502:    Concepts: sorting^ints with array

504: .seealso: PetscSortReal(), PetscSortIntPermutation(), PetscSortInt()
505: @*/
506: PetscErrorCode  PetscSortMPIIntWithArray(PetscMPIInt n,PetscMPIInt i[],PetscMPIInt Ii[])
507: {
509:   PetscMPIInt    j,k,tmp,ik;

512:   if (n<8) {
513:     for (k=0; k<n; k++) {
514:       ik = i[k];
515:       for (j=k+1; j<n; j++) {
516:         if (ik > i[j]) {
517:           SWAP2(i[k],i[j],Ii[k],Ii[j],tmp);
518:           ik = i[k];
519:         }
520:       }
521:     }
522:   } else {
523:     PetscSortMPIIntWithArray_Private(i,Ii,n-1);
524:   }
525:   return(0);
526: }

528: /* -----------------------------------------------------------------------*/
529: #define SWAP2IntScalar(a,b,c,d,t,ts) {t=a;a=b;b=t;ts=c;c=d;d=ts;}

531: /*
532:    Modified from PetscSortIntWithArray_Private().
533: */
534: static PetscErrorCode PetscSortIntWithScalarArray_Private(PetscInt *v,PetscScalar *V,PetscInt right)
535: {
537:   PetscInt       i,vl,last,tmp;
538:   PetscScalar    stmp;

541:   if (right <= 1) {
542:     if (right == 1) {
543:       if (v[0] > v[1]) SWAP2IntScalar(v[0],v[1],V[0],V[1],tmp,stmp);
544:     }
545:     return(0);
546:   }
547:   SWAP2IntScalar(v[0],v[right/2],V[0],V[right/2],tmp,stmp);
548:   vl   = v[0];
549:   last = 0;
550:   for (i=1; i<=right; i++) {
551:     if (v[i] < vl) {last++; SWAP2IntScalar(v[last],v[i],V[last],V[i],tmp,stmp);}
552:   }
553:   SWAP2IntScalar(v[0],v[last],V[0],V[last],tmp,stmp);
554:   PetscSortIntWithScalarArray_Private(v,V,last-1);
555:   PetscSortIntWithScalarArray_Private(v+last+1,V+last+1,right-(last+1));
556:   return(0);
557: }

559: /*@
560:    PetscSortIntWithScalarArray - Sorts an array of integers in place in increasing order;
561:        changes a second SCALAR array to match the sorted first INTEGER array.

563:    Not Collective

565:    Input Parameters:
566: +  n  - number of values
567: .  i  - array of integers
568: -  I - second array of scalars

570:    Level: intermediate

572:    Concepts: sorting^ints with array

574: .seealso: PetscSortReal(), PetscSortIntPermutation(), PetscSortInt(), PetscSortIntWithArray()
575: @*/
576: PetscErrorCode  PetscSortIntWithScalarArray(PetscInt n,PetscInt i[],PetscScalar Ii[])
577: {
579:   PetscInt       j,k,tmp,ik;
580:   PetscScalar    stmp;

583:   if (n<8) {
584:     for (k=0; k<n; k++) {
585:       ik = i[k];
586:       for (j=k+1; j<n; j++) {
587:         if (ik > i[j]) {
588:           SWAP2IntScalar(i[k],i[j],Ii[k],Ii[j],tmp,stmp);
589:           ik = i[k];
590:         }
591:       }
592:     }
593:   } else {
594:     PetscSortIntWithScalarArray_Private(i,Ii,n-1);
595:   }
596:   return(0);
597: }

599: #define SWAP2IntData(a,b,c,d,t,td,siz)          \
600:   do {                                          \
601:   PetscErrorCode _ierr;                         \
602:   t=a;a=b;b=t;                                  \
603:   _PetscMemcpy(td,c,siz);CHKERRQ(_ierr); \
604:   _PetscMemcpy(c,d,siz);CHKERRQ(_ierr);  \
605:   _PetscMemcpy(d,td,siz);CHKERRQ(_ierr); \
606: } while(0)

608: /*
609:    Modified from PetscSortIntWithArray_Private().
610: */
611: static PetscErrorCode PetscSortIntWithDataArray_Private(PetscInt *v,char *V,PetscInt right,size_t size,void *work)
612: {
614:   PetscInt       i,vl,last,tmp;

617:   if (right <= 1) {
618:     if (right == 1) {
619:       if (v[0] > v[1]) SWAP2IntData(v[0],v[1],V,V+size,tmp,work,size);
620:     }
621:     return(0);
622:   }
623:   SWAP2IntData(v[0],v[right/2],V,V+size*(right/2),tmp,work,size);
624:   vl   = v[0];
625:   last = 0;
626:   for (i=1; i<=right; i++) {
627:     if (v[i] < vl) {last++; SWAP2IntData(v[last],v[i],V+size*last,V+size*i,tmp,work,size);}
628:   }
629:   SWAP2IntData(v[0],v[last],V,V + size*last,tmp,work,size);
630:   PetscSortIntWithDataArray_Private(v,V,last-1,size,work);
631:   PetscSortIntWithDataArray_Private(v+last+1,V+size*(last+1),right-(last+1),size,work);
632:   return(0);
633: }

635: /*@C
636:    PetscSortIntWithDataArray - Sorts an array of integers in place in increasing order;
637:        changes a second array to match the sorted first INTEGER array.  Unlike other sort routines, the user must
638:        provide workspace (the size of an element in the data array) to use when sorting.

640:    Not Collective

642:    Input Parameters:
643: +  n  - number of values
644: .  i  - array of integers
645: .  Ii - second array of data
646: .  size - sizeof elements in the data array in bytes
647: -  work - workspace of "size" bytes used when sorting

649:    Level: intermediate

651:    Concepts: sorting^ints with array

653: .seealso: PetscSortReal(), PetscSortIntPermutation(), PetscSortInt(), PetscSortIntWithArray()
654: @*/
655: PetscErrorCode  PetscSortIntWithDataArray(PetscInt n,PetscInt i[],void *Ii,size_t size,void *work)
656: {
657:   char           *V = (char *) Ii;
659:   PetscInt       j,k,tmp,ik;

662:   if (n<8) {
663:     for (k=0; k<n; k++) {
664:       ik = i[k];
665:       for (j=k+1; j<n; j++) {
666:         if (ik > i[j]) {
667:           SWAP2IntData(i[k],i[j],V+size*k,V+size*j,tmp,work,size);
668:           ik = i[k];
669:         }
670:       }
671:     }
672:   } else {
673:     PetscSortIntWithDataArray_Private(i,V,n-1,size,work);
674:   }
675:   return(0);
676: }


679: /*@
680:    PetscMergeIntArray -     Merges two SORTED integer arrays, removes duplicate elements.

682:    Not Collective

684:    Input Parameters:
685: +  an  - number of values in the first array
686: .  aI  - first sorted array of integers
687: .  bn  - number of values in the second array
688: -  bI  - second array of integers

690:    Output Parameters:
691: +  n   - number of values in the merged array
692: -  L   - merged sorted array, this is allocated if an array is not provided

694:    Level: intermediate

696:    Concepts: merging^arrays

698: .seealso: PetscSortReal(), PetscSortIntPermutation(), PetscSortInt(), PetscSortIntWithArray()
699: @*/
700: PetscErrorCode  PetscMergeIntArray(PetscInt an,const PetscInt aI[], PetscInt bn, const PetscInt bI[],  PetscInt *n, PetscInt **L)
701: {
703:   PetscInt       *L_ = *L, ak, bk, k;

705:   if (!L_) {
706:     PetscMalloc1(an+bn, L);
707:     L_   = *L;
708:   }
709:   k = ak = bk = 0;
710:   while (ak < an && bk < bn) {
711:     if (aI[ak] == bI[bk]) {
712:       L_[k] = aI[ak];
713:       ++ak;
714:       ++bk;
715:       ++k;
716:     } else if (aI[ak] < bI[bk]) {
717:       L_[k] = aI[ak];
718:       ++ak;
719:       ++k;
720:     } else {
721:       L_[k] = bI[bk];
722:       ++bk;
723:       ++k;
724:     }
725:   }
726:   if (ak < an) {
727:     PetscMemcpy(L_+k,aI+ak,(an-ak)*sizeof(PetscInt));
728:     k   += (an-ak);
729:   }
730:   if (bk < bn) {
731:     PetscMemcpy(L_+k,bI+bk,(bn-bk)*sizeof(PetscInt));
732:     k   += (bn-bk);
733:   }
734:   *n = k;
735:   return(0);
736: }

738: /*@
739:    PetscMergeIntArrayPair -     Merges two SORTED integer arrays that share NO common values along with an additional array of integers.
740:                                 The additional arrays are the same length as sorted arrays and are merged
741:                                 in the order determined by the merging of the sorted pair.

743:    Not Collective

745:    Input Parameters:
746: +  an  - number of values in the first array
747: .  aI  - first sorted array of integers
748: .  aJ  - first additional array of integers
749: .  bn  - number of values in the second array
750: .  bI  - second array of integers
751: -  bJ  - second additional of integers

753:    Output Parameters:
754: +  n   - number of values in the merged array (== an + bn)
755: .  L   - merged sorted array 
756: -  J   - merged additional array

758:    Notes: 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 
759:    Level: intermediate

761:    Concepts: merging^arrays

763: .seealso: PetscSortReal(), PetscSortIntPermutation(), PetscSortInt(), PetscSortIntWithArray()
764: @*/
765: PetscErrorCode  PetscMergeIntArrayPair(PetscInt an,const PetscInt aI[], const PetscInt aJ[], PetscInt bn, const PetscInt bI[], const PetscInt bJ[], PetscInt *n, PetscInt **L, PetscInt **J)
766: {
768:   PetscInt       n_, *L_, *J_, ak, bk, k;

773:   n_ = an + bn;
774:   *n = n_;
775:   if (!*L) {
776:     PetscMalloc1(n_, L);
777:   }
778:   L_ = *L;
779:   if (!*J) {
780:     PetscMalloc1(n_, J);
781:   }
782:   J_   = *J;
783:   k = ak = bk = 0;
784:   while (ak < an && bk < bn) {
785:     if (aI[ak] <= bI[bk]) {
786:       L_[k] = aI[ak];
787:       J_[k] = aJ[ak];
788:       ++ak;
789:       ++k;
790:     } else {
791:       L_[k] = bI[bk];
792:       J_[k] = bJ[bk];
793:       ++bk;
794:       ++k;
795:     }
796:   }
797:   if (ak < an) {
798:     PetscMemcpy(L_+k,aI+ak,(an-ak)*sizeof(PetscInt));
799:     PetscMemcpy(J_+k,aJ+ak,(an-ak)*sizeof(PetscInt));
800:     k   += (an-ak);
801:   }
802:   if (bk < bn) {
803:     PetscMemcpy(L_+k,bI+bk,(bn-bk)*sizeof(PetscInt));
804:     PetscMemcpy(J_+k,bJ+bk,(bn-bk)*sizeof(PetscInt));
805:   }
806:   return(0);
807: }

809: /*@
810:    PetscMergeMPIIntArray -     Merges two SORTED integer arrays.

812:    Not Collective

814:    Input Parameters:
815: +  an  - number of values in the first array
816: .  aI  - first sorted array of integers
817: .  bn  - number of values in the second array
818: -  bI  - second array of integers

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

824:    Level: intermediate

826:    Concepts: merging^arrays

828: .seealso: PetscSortReal(), PetscSortIntPermutation(), PetscSortInt(), PetscSortIntWithArray()
829: @*/
830: PetscErrorCode PetscMergeMPIIntArray(PetscInt an,const PetscMPIInt aI[],PetscInt bn,const PetscMPIInt bI[],PetscInt *n,PetscMPIInt **L)
831: {
833:   PetscInt       ai,bi,k;

836:   if (!*L) {PetscMalloc1((an+bn),L);}
837:   for (ai=0,bi=0,k=0; ai<an || bi<bn; ) {
838:     PetscInt t = -1;
839:     for ( ; ai<an && (!bn || aI[ai] <= bI[bi]); ai++) (*L)[k++] = t = aI[ai];
840:     for ( ; bi<bn && bI[bi] == t; bi++);
841:     for ( ; bi<bn && (!an || bI[bi] <= aI[ai]); bi++) (*L)[k++] = t = bI[bi];
842:     for ( ; ai<an && aI[ai] == t; ai++);
843:   }
844:   *n = k;
845:   return(0);
846: }

848: /*@C
849:    PetscProcessTree - Prepares tree data to be displayed graphically

851:    Not Collective

853:    Input Parameters:
854: +  n  - number of values
855: .  mask - indicates those entries in the tree, location 0 is always masked
856: -  parentid - indicates the parent of each entry

858:    Output Parameters:
859: +  Nlevels - the number of levels
860: .  Level - for each node tells its level
861: .  Levelcnts - the number of nodes on each level
862: .  Idbylevel - a list of ids on each of the levels, first level followed by second etc
863: -  Column - for each id tells its column index

865:    Level: developer

867:    Notes: This code is not currently used

869: .seealso: PetscSortReal(), PetscSortIntWithPermutation()
870: @*/
871: PetscErrorCode  PetscProcessTree(PetscInt n,const PetscBool mask[],const PetscInt parentid[],PetscInt *Nlevels,PetscInt **Level,PetscInt **Levelcnt,PetscInt **Idbylevel,PetscInt **Column)
872: {
873:   PetscInt       i,j,cnt,nmask = 0,nlevels = 0,*level,*levelcnt,levelmax = 0,*workid,*workparentid,tcnt = 0,*idbylevel,*column;
875:   PetscBool      done = PETSC_FALSE;

878:   if (!mask[0]) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Mask of 0th location must be set");
879:   for (i=0; i<n; i++) {
880:     if (mask[i]) continue;
881:     if (parentid[i]  == i) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Node labeled as own parent");
882:     if (parentid[i] && mask[parentid[i]]) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Parent is masked");
883:   }

885:   for (i=0; i<n; i++) {
886:     if (!mask[i]) nmask++;
887:   }

889:   /* determine the level in the tree of each node */
890:   PetscCalloc1(n,&level);

892:   level[0] = 1;
893:   while (!done) {
894:     done = PETSC_TRUE;
895:     for (i=0; i<n; i++) {
896:       if (mask[i]) continue;
897:       if (!level[i] && level[parentid[i]]) level[i] = level[parentid[i]] + 1;
898:       else if (!level[i]) done = PETSC_FALSE;
899:     }
900:   }
901:   for (i=0; i<n; i++) {
902:     level[i]--;
903:     nlevels = PetscMax(nlevels,level[i]);
904:   }

906:   /* count the number of nodes on each level and its max */
907:   PetscCalloc1(nlevels,&levelcnt);
908:   for (i=0; i<n; i++) {
909:     if (mask[i]) continue;
910:     levelcnt[level[i]-1]++;
911:   }
912:   for (i=0; i<nlevels;i++) levelmax = PetscMax(levelmax,levelcnt[i]);

914:   /* for each level sort the ids by the parent id */
915:   PetscMalloc2(levelmax,&workid,levelmax,&workparentid);
916:   PetscMalloc1(nmask,&idbylevel);
917:   for (j=1; j<=nlevels;j++) {
918:     cnt = 0;
919:     for (i=0; i<n; i++) {
920:       if (mask[i]) continue;
921:       if (level[i] != j) continue;
922:       workid[cnt]         = i;
923:       workparentid[cnt++] = parentid[i];
924:     }
925:     /*  PetscIntView(cnt,workparentid,0);
926:     PetscIntView(cnt,workid,0);
927:     PetscSortIntWithArray(cnt,workparentid,workid);
928:     PetscIntView(cnt,workparentid,0);
929:     PetscIntView(cnt,workid,0);*/
930:     PetscMemcpy(idbylevel+tcnt,workid,cnt*sizeof(PetscInt));
931:     tcnt += cnt;
932:   }
933:   if (tcnt != nmask) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"Inconsistent count of unmasked nodes");
934:   PetscFree2(workid,workparentid);

936:   /* for each node list its column */
937:   PetscMalloc1(n,&column);
938:   cnt = 0;
939:   for (j=0; j<nlevels; j++) {
940:     for (i=0; i<levelcnt[j]; i++) {
941:       column[idbylevel[cnt++]] = i;
942:     }
943:   }

945:   *Nlevels   = nlevels;
946:   *Level     = level;
947:   *Levelcnt  = levelcnt;
948:   *Idbylevel = idbylevel;
949:   *Column    = column;
950:   return(0);
951: }