Actual source code: sorti.c
petsc-3.11.4 2019-09-28
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 SWAP2MPIIntInt(a,b,c,d,t1,t2) {t1=a;a=b;b=t1;t2=c;c=d;d=t2;}
531: /*
532: A simple version of quicksort; taken from Kernighan and Ritchie, page 87.
533: Assumes 0 origin for v, number of elements = right+1 (right is index of
534: right-most member).
535: */
536: static PetscErrorCode PetscSortMPIIntWithIntArray_Private(PetscMPIInt *v,PetscInt *V,PetscMPIInt right)
537: {
539: PetscMPIInt i,vl,last,t1;
540: PetscInt t2;
543: if (right <= 1) {
544: if (right == 1) {
545: if (v[0] > v[1]) SWAP2MPIIntInt(v[0],v[1],V[0],V[1],t1,t2);
546: }
547: return(0);
548: }
549: SWAP2MPIIntInt(v[0],v[right/2],V[0],V[right/2],t1,t2);
550: vl = v[0];
551: last = 0;
552: for (i=1; i<=right; i++) {
553: if (v[i] < vl) {last++; SWAP2MPIIntInt(v[last],v[i],V[last],V[i],t1,t2);}
554: }
555: SWAP2MPIIntInt(v[0],v[last],V[0],V[last],t1,t2);
556: PetscSortMPIIntWithIntArray_Private(v,V,last-1);
557: PetscSortMPIIntWithIntArray_Private(v+last+1,V+last+1,right-(last+1));
558: return(0);
559: }
561: /*@
562: PetscSortMPIIntWithIntArray - Sorts an array of MPI integers in place in increasing order;
563: changes a second array of Petsc intergers to match the sorted first array.
565: Not Collective
567: Input Parameters:
568: + n - number of values
569: . i - array of MPI integers
570: - I - second array of Petsc integers
572: Level: intermediate
574: Notes: this routine is useful when one needs to sort MPI ranks with other integer arrays.
576: Concepts: sorting^ints with array
578: .seealso: PetscSortMPIIntWithArray()
579: @*/
580: PetscErrorCode PetscSortMPIIntWithIntArray(PetscMPIInt n,PetscMPIInt i[],PetscInt Ii[])
581: {
583: PetscMPIInt j,k,t1,ik;
584: PetscInt t2;
587: if (n<8) {
588: for (k=0; k<n; k++) {
589: ik = i[k];
590: for (j=k+1; j<n; j++) {
591: if (ik > i[j]) {
592: SWAP2MPIIntInt(i[k],i[j],Ii[k],Ii[j],t1,t2);
593: ik = i[k];
594: }
595: }
596: }
597: } else {
598: PetscSortMPIIntWithIntArray_Private(i,Ii,n-1);
599: }
600: return(0);
601: }
603: /* -----------------------------------------------------------------------*/
604: #define SWAP2IntScalar(a,b,c,d,t,ts) {t=a;a=b;b=t;ts=c;c=d;d=ts;}
606: /*
607: Modified from PetscSortIntWithArray_Private().
608: */
609: static PetscErrorCode PetscSortIntWithScalarArray_Private(PetscInt *v,PetscScalar *V,PetscInt right)
610: {
612: PetscInt i,vl,last,tmp;
613: PetscScalar stmp;
616: if (right <= 1) {
617: if (right == 1) {
618: if (v[0] > v[1]) SWAP2IntScalar(v[0],v[1],V[0],V[1],tmp,stmp);
619: }
620: return(0);
621: }
622: SWAP2IntScalar(v[0],v[right/2],V[0],V[right/2],tmp,stmp);
623: vl = v[0];
624: last = 0;
625: for (i=1; i<=right; i++) {
626: if (v[i] < vl) {last++; SWAP2IntScalar(v[last],v[i],V[last],V[i],tmp,stmp);}
627: }
628: SWAP2IntScalar(v[0],v[last],V[0],V[last],tmp,stmp);
629: PetscSortIntWithScalarArray_Private(v,V,last-1);
630: PetscSortIntWithScalarArray_Private(v+last+1,V+last+1,right-(last+1));
631: return(0);
632: }
634: /*@
635: PetscSortIntWithScalarArray - Sorts an array of integers in place in increasing order;
636: changes a second SCALAR array to match the sorted first INTEGER array.
638: Not Collective
640: Input Parameters:
641: + n - number of values
642: . i - array of integers
643: - I - second array of scalars
645: Level: intermediate
647: Concepts: sorting^ints with array
649: .seealso: PetscSortReal(), PetscSortIntPermutation(), PetscSortInt(), PetscSortIntWithArray()
650: @*/
651: PetscErrorCode PetscSortIntWithScalarArray(PetscInt n,PetscInt i[],PetscScalar Ii[])
652: {
654: PetscInt j,k,tmp,ik;
655: PetscScalar stmp;
658: if (n<8) {
659: for (k=0; k<n; k++) {
660: ik = i[k];
661: for (j=k+1; j<n; j++) {
662: if (ik > i[j]) {
663: SWAP2IntScalar(i[k],i[j],Ii[k],Ii[j],tmp,stmp);
664: ik = i[k];
665: }
666: }
667: }
668: } else {
669: PetscSortIntWithScalarArray_Private(i,Ii,n-1);
670: }
671: return(0);
672: }
674: #define SWAP2IntData(a,b,c,d,t,td,siz) \
675: do { \
676: PetscErrorCode _ierr; \
677: t=a;a=b;b=t; \
678: _PetscMemcpy(td,c,siz);CHKERRQ(_ierr); \
679: _PetscMemcpy(c,d,siz);CHKERRQ(_ierr); \
680: _PetscMemcpy(d,td,siz);CHKERRQ(_ierr); \
681: } while(0)
683: /*
684: Modified from PetscSortIntWithArray_Private().
685: */
686: static PetscErrorCode PetscSortIntWithDataArray_Private(PetscInt *v,char *V,PetscInt right,size_t size,void *work)
687: {
689: PetscInt i,vl,last,tmp;
692: if (right <= 1) {
693: if (right == 1) {
694: if (v[0] > v[1]) SWAP2IntData(v[0],v[1],V,V+size,tmp,work,size);
695: }
696: return(0);
697: }
698: SWAP2IntData(v[0],v[right/2],V,V+size*(right/2),tmp,work,size);
699: vl = v[0];
700: last = 0;
701: for (i=1; i<=right; i++) {
702: if (v[i] < vl) {last++; SWAP2IntData(v[last],v[i],V+size*last,V+size*i,tmp,work,size);}
703: }
704: SWAP2IntData(v[0],v[last],V,V + size*last,tmp,work,size);
705: PetscSortIntWithDataArray_Private(v,V,last-1,size,work);
706: PetscSortIntWithDataArray_Private(v+last+1,V+size*(last+1),right-(last+1),size,work);
707: return(0);
708: }
710: /*@C
711: PetscSortIntWithDataArray - Sorts an array of integers in place in increasing order;
712: changes a second array to match the sorted first INTEGER array. Unlike other sort routines, the user must
713: provide workspace (the size of an element in the data array) to use when sorting.
715: Not Collective
717: Input Parameters:
718: + n - number of values
719: . i - array of integers
720: . Ii - second array of data
721: . size - sizeof elements in the data array in bytes
722: - work - workspace of "size" bytes used when sorting
724: Level: intermediate
726: Concepts: sorting^ints with array
728: .seealso: PetscSortReal(), PetscSortIntPermutation(), PetscSortInt(), PetscSortIntWithArray()
729: @*/
730: PetscErrorCode PetscSortIntWithDataArray(PetscInt n,PetscInt i[],void *Ii,size_t size,void *work)
731: {
732: char *V = (char *) Ii;
734: PetscInt j,k,tmp,ik;
737: if (n<8) {
738: for (k=0; k<n; k++) {
739: ik = i[k];
740: for (j=k+1; j<n; j++) {
741: if (ik > i[j]) {
742: SWAP2IntData(i[k],i[j],V+size*k,V+size*j,tmp,work,size);
743: ik = i[k];
744: }
745: }
746: }
747: } else {
748: PetscSortIntWithDataArray_Private(i,V,n-1,size,work);
749: }
750: return(0);
751: }
754: /*@
755: PetscMergeIntArray - Merges two SORTED integer arrays, removes duplicate elements.
757: Not Collective
759: Input Parameters:
760: + an - number of values in the first array
761: . aI - first sorted array of integers
762: . bn - number of values in the second array
763: - bI - second array of integers
765: Output Parameters:
766: + n - number of values in the merged array
767: - L - merged sorted array, this is allocated if an array is not provided
769: Level: intermediate
771: Concepts: merging^arrays
773: .seealso: PetscSortReal(), PetscSortIntPermutation(), PetscSortInt(), PetscSortIntWithArray()
774: @*/
775: PetscErrorCode PetscMergeIntArray(PetscInt an,const PetscInt aI[], PetscInt bn, const PetscInt bI[], PetscInt *n, PetscInt **L)
776: {
778: PetscInt *L_ = *L, ak, bk, k;
780: if (!L_) {
781: PetscMalloc1(an+bn, L);
782: L_ = *L;
783: }
784: k = ak = bk = 0;
785: while (ak < an && bk < bn) {
786: if (aI[ak] == bI[bk]) {
787: L_[k] = aI[ak];
788: ++ak;
789: ++bk;
790: ++k;
791: } else if (aI[ak] < bI[bk]) {
792: L_[k] = aI[ak];
793: ++ak;
794: ++k;
795: } else {
796: L_[k] = bI[bk];
797: ++bk;
798: ++k;
799: }
800: }
801: if (ak < an) {
802: PetscMemcpy(L_+k,aI+ak,(an-ak)*sizeof(PetscInt));
803: k += (an-ak);
804: }
805: if (bk < bn) {
806: PetscMemcpy(L_+k,bI+bk,(bn-bk)*sizeof(PetscInt));
807: k += (bn-bk);
808: }
809: *n = k;
810: return(0);
811: }
813: /*@
814: PetscMergeIntArrayPair - Merges two SORTED integer arrays that share NO common values along with an additional array of integers.
815: The additional arrays are the same length as sorted arrays and are merged
816: in the order determined by the merging of the sorted pair.
818: Not Collective
820: Input Parameters:
821: + an - number of values in the first array
822: . aI - first sorted array of integers
823: . aJ - first additional array of integers
824: . bn - number of values in the second array
825: . bI - second array of integers
826: - bJ - second additional of integers
828: Output Parameters:
829: + n - number of values in the merged array (== an + bn)
830: . L - merged sorted array
831: - J - merged additional array
833: Notes:
834: 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
835: Level: intermediate
837: Concepts: merging^arrays
839: .seealso: PetscSortReal(), PetscSortIntPermutation(), PetscSortInt(), PetscSortIntWithArray()
840: @*/
841: PetscErrorCode PetscMergeIntArrayPair(PetscInt an,const PetscInt aI[], const PetscInt aJ[], PetscInt bn, const PetscInt bI[], const PetscInt bJ[], PetscInt *n, PetscInt **L, PetscInt **J)
842: {
844: PetscInt n_, *L_, *J_, ak, bk, k;
849: n_ = an + bn;
850: *n = n_;
851: if (!*L) {
852: PetscMalloc1(n_, L);
853: }
854: L_ = *L;
855: if (!*J) {
856: PetscMalloc1(n_, J);
857: }
858: J_ = *J;
859: k = ak = bk = 0;
860: while (ak < an && bk < bn) {
861: if (aI[ak] <= bI[bk]) {
862: L_[k] = aI[ak];
863: J_[k] = aJ[ak];
864: ++ak;
865: ++k;
866: } else {
867: L_[k] = bI[bk];
868: J_[k] = bJ[bk];
869: ++bk;
870: ++k;
871: }
872: }
873: if (ak < an) {
874: PetscMemcpy(L_+k,aI+ak,(an-ak)*sizeof(PetscInt));
875: PetscMemcpy(J_+k,aJ+ak,(an-ak)*sizeof(PetscInt));
876: k += (an-ak);
877: }
878: if (bk < bn) {
879: PetscMemcpy(L_+k,bI+bk,(bn-bk)*sizeof(PetscInt));
880: PetscMemcpy(J_+k,bJ+bk,(bn-bk)*sizeof(PetscInt));
881: }
882: return(0);
883: }
885: /*@
886: PetscMergeMPIIntArray - Merges two SORTED integer arrays.
888: Not Collective
890: Input Parameters:
891: + an - number of values in the first array
892: . aI - first sorted array of integers
893: . bn - number of values in the second array
894: - bI - second array of integers
896: Output Parameters:
897: + n - number of values in the merged array (<= an + bn)
898: - L - merged sorted array, allocated if address of NULL pointer is passed
900: Level: intermediate
902: Concepts: merging^arrays
904: .seealso: PetscSortReal(), PetscSortIntPermutation(), PetscSortInt(), PetscSortIntWithArray()
905: @*/
906: PetscErrorCode PetscMergeMPIIntArray(PetscInt an,const PetscMPIInt aI[],PetscInt bn,const PetscMPIInt bI[],PetscInt *n,PetscMPIInt **L)
907: {
909: PetscInt ai,bi,k;
912: if (!*L) {PetscMalloc1((an+bn),L);}
913: for (ai=0,bi=0,k=0; ai<an || bi<bn; ) {
914: PetscInt t = -1;
915: for ( ; ai<an && (!bn || aI[ai] <= bI[bi]); ai++) (*L)[k++] = t = aI[ai];
916: for ( ; bi<bn && bI[bi] == t; bi++);
917: for ( ; bi<bn && (!an || bI[bi] <= aI[ai]); bi++) (*L)[k++] = t = bI[bi];
918: for ( ; ai<an && aI[ai] == t; ai++);
919: }
920: *n = k;
921: return(0);
922: }
924: /*@C
925: PetscProcessTree - Prepares tree data to be displayed graphically
927: Not Collective
929: Input Parameters:
930: + n - number of values
931: . mask - indicates those entries in the tree, location 0 is always masked
932: - parentid - indicates the parent of each entry
934: Output Parameters:
935: + Nlevels - the number of levels
936: . Level - for each node tells its level
937: . Levelcnts - the number of nodes on each level
938: . Idbylevel - a list of ids on each of the levels, first level followed by second etc
939: - Column - for each id tells its column index
941: Level: developer
943: Notes:
944: This code is not currently used
946: .seealso: PetscSortReal(), PetscSortIntWithPermutation()
947: @*/
948: PetscErrorCode PetscProcessTree(PetscInt n,const PetscBool mask[],const PetscInt parentid[],PetscInt *Nlevels,PetscInt **Level,PetscInt **Levelcnt,PetscInt **Idbylevel,PetscInt **Column)
949: {
950: PetscInt i,j,cnt,nmask = 0,nlevels = 0,*level,*levelcnt,levelmax = 0,*workid,*workparentid,tcnt = 0,*idbylevel,*column;
952: PetscBool done = PETSC_FALSE;
955: if (!mask[0]) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Mask of 0th location must be set");
956: for (i=0; i<n; i++) {
957: if (mask[i]) continue;
958: if (parentid[i] == i) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Node labeled as own parent");
959: if (parentid[i] && mask[parentid[i]]) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Parent is masked");
960: }
962: for (i=0; i<n; i++) {
963: if (!mask[i]) nmask++;
964: }
966: /* determine the level in the tree of each node */
967: PetscCalloc1(n,&level);
969: level[0] = 1;
970: while (!done) {
971: done = PETSC_TRUE;
972: for (i=0; i<n; i++) {
973: if (mask[i]) continue;
974: if (!level[i] && level[parentid[i]]) level[i] = level[parentid[i]] + 1;
975: else if (!level[i]) done = PETSC_FALSE;
976: }
977: }
978: for (i=0; i<n; i++) {
979: level[i]--;
980: nlevels = PetscMax(nlevels,level[i]);
981: }
983: /* count the number of nodes on each level and its max */
984: PetscCalloc1(nlevels,&levelcnt);
985: for (i=0; i<n; i++) {
986: if (mask[i]) continue;
987: levelcnt[level[i]-1]++;
988: }
989: for (i=0; i<nlevels;i++) levelmax = PetscMax(levelmax,levelcnt[i]);
991: /* for each level sort the ids by the parent id */
992: PetscMalloc2(levelmax,&workid,levelmax,&workparentid);
993: PetscMalloc1(nmask,&idbylevel);
994: for (j=1; j<=nlevels;j++) {
995: cnt = 0;
996: for (i=0; i<n; i++) {
997: if (mask[i]) continue;
998: if (level[i] != j) continue;
999: workid[cnt] = i;
1000: workparentid[cnt++] = parentid[i];
1001: }
1002: /* PetscIntView(cnt,workparentid,0);
1003: PetscIntView(cnt,workid,0);
1004: PetscSortIntWithArray(cnt,workparentid,workid);
1005: PetscIntView(cnt,workparentid,0);
1006: PetscIntView(cnt,workid,0);*/
1007: PetscMemcpy(idbylevel+tcnt,workid,cnt*sizeof(PetscInt));
1008: tcnt += cnt;
1009: }
1010: if (tcnt != nmask) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"Inconsistent count of unmasked nodes");
1011: PetscFree2(workid,workparentid);
1013: /* for each node list its column */
1014: PetscMalloc1(n,&column);
1015: cnt = 0;
1016: for (j=0; j<nlevels; j++) {
1017: for (i=0; i<levelcnt[j]; i++) {
1018: column[idbylevel[cnt++]] = i;
1019: }
1020: }
1022: *Nlevels = nlevels;
1023: *Level = level;
1024: *Levelcnt = levelcnt;
1025: *Idbylevel = idbylevel;
1026: *Column = column;
1027: return(0);
1028: }