Actual source code: sorti.c
petsc-3.9.4 2018-09-11
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: }