Actual source code: sorti.c
2: /*
3: This file contains routines for sorting integers. Values are sorted in place.
4: */
5: #include <petscsys.h> /*I "petscsys.h" I*/
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: /* -----------------------------------------------------------------------*/
24: /*
25: A simple version of quicksort; taken from Kernighan and Ritchie, page 87.
26: Assumes 0 origin for v, number of elements = right+1 (right is index of
27: right-most member).
28: */
29: static void PetscSortInt_Private(PetscInt *v,PetscInt right)
30: {
31: PetscInt i,j,pivot,tmp;
33: if (right <= 1) {
34: if (right == 1) {
35: if (v[0] > v[1]) SWAP(v[0],v[1],tmp);
36: }
37: return;
38: }
39: i = MEDIAN(v,right); /* Choose a pivot */
40: SWAP(v[0],v[i],tmp); /* Move it out of the way */
41: pivot = v[0];
42: for (i=0,j=right+1;;) {
43: while (++i < j && v[i] <= pivot); /* Scan from the left */
44: while (v[--j] > pivot); /* Scan from the right */
45: if (i >= j) break;
46: SWAP(v[i],v[j],tmp);
47: }
48: SWAP(v[0],v[j],tmp); /* Put pivot back in place. */
49: PetscSortInt_Private(v,j-1);
50: PetscSortInt_Private(v+j+1,right-(j+1));
51: }
55: /*@
56: PetscSortInt - Sorts an array of integers in place in increasing order.
58: Not Collective
60: Input Parameters:
61: + n - number of values
62: - i - array of integers
64: Level: intermediate
66: Concepts: sorting^ints
68: .seealso: PetscSortReal(), PetscSortIntWithPermutation()
69: @*/
70: PetscErrorCode PetscSortInt(PetscInt n,PetscInt i[])
71: {
72: PetscInt j,k,tmp,ik;
75: if (n<8) {
76: for (k=0; k<n; k++) {
77: ik = i[k];
78: for (j=k+1; j<n; j++) {
79: if (ik > i[j]) {
80: SWAP(i[k],i[j],tmp);
81: ik = i[k];
82: }
83: }
84: }
85: } else {
86: PetscSortInt_Private(i,n-1);
87: }
88: return(0);
89: }
93: /*@
94: PetscSortRemoveDupsInt - Sorts an array of integers in place in increasing order removes all duplicate entries
96: Not Collective
98: Input Parameters:
99: + n - number of values
100: - ii - array of integers
102: Output Parameter:
103: . n - number of non-redundant values
105: Level: intermediate
107: Concepts: sorting^ints
109: .seealso: PetscSortReal(), PetscSortIntWithPermutation(), PetscSortInt()
110: @*/
111: PetscErrorCode PetscSortRemoveDupsInt(PetscInt *n,PetscInt ii[])
112: {
114: PetscInt i,s = 0,N = *n, b = 0;
117: PetscSortInt(N,ii);
118: for (i=0; i<N-1; i++) {
119: if (ii[b+s+1] != ii[b]) {ii[b+1] = ii[b+s+1]; b++;}
120: else s++;
121: }
122: *n = N - s;
123: return(0);
124: }
127: /* -----------------------------------------------------------------------*/
128: #define SWAP2(a,b,c,d,t) {t=a;a=b;b=t;t=c;c=d;d=t;}
132: /*
133: A simple version of quicksort; taken from Kernighan and Ritchie, page 87.
134: Assumes 0 origin for v, number of elements = right+1 (right is index of
135: right-most member).
136: */
137: static PetscErrorCode PetscSortIntWithArray_Private(PetscInt *v,PetscInt *V,PetscInt right)
138: {
140: PetscInt i,vl,last,tmp;
143: if (right <= 1) {
144: if (right == 1) {
145: if (v[0] > v[1]) SWAP2(v[0],v[1],V[0],V[1],tmp);
146: }
147: return(0);
148: }
149: SWAP2(v[0],v[right/2],V[0],V[right/2],tmp);
150: vl = v[0];
151: last = 0;
152: for (i=1; i<=right; i++) {
153: if (v[i] < vl) {last++; SWAP2(v[last],v[i],V[last],V[i],tmp);}
154: }
155: SWAP2(v[0],v[last],V[0],V[last],tmp);
156: PetscSortIntWithArray_Private(v,V,last-1);
157: PetscSortIntWithArray_Private(v+last+1,V+last+1,right-(last+1));
158: return(0);
159: }
163: /*@
164: PetscSortIntWithArray - Sorts an array of integers in place in increasing order;
165: changes a second array to match the sorted first array.
167: Not Collective
169: Input Parameters:
170: + n - number of values
171: . i - array of integers
172: - I - second array of integers
174: Level: intermediate
176: Concepts: sorting^ints with array
178: .seealso: PetscSortReal(), PetscSortIntPermutation(), PetscSortInt()
179: @*/
180: PetscErrorCode PetscSortIntWithArray(PetscInt n,PetscInt i[],PetscInt Ii[])
181: {
183: PetscInt j,k,tmp,ik;
186: if (n<8) {
187: for (k=0; k<n; k++) {
188: ik = i[k];
189: for (j=k+1; j<n; j++) {
190: if (ik > i[j]) {
191: SWAP2(i[k],i[j],Ii[k],Ii[j],tmp);
192: ik = i[k];
193: }
194: }
195: }
196: } else {
197: PetscSortIntWithArray_Private(i,Ii,n-1);
198: }
199: return(0);
200: }
204: /*
205: A simple version of quicksort; taken from Kernighan and Ritchie, page 87.
206: Assumes 0 origin for v, number of elements = right+1 (right is index of
207: right-most member).
208: */
209: static PetscErrorCode PetscSortMPIIntWithArray_Private(PetscMPIInt *v,PetscMPIInt *V,PetscMPIInt right)
210: {
212: PetscMPIInt i,vl,last,tmp;
215: if (right <= 1) {
216: if (right == 1) {
217: if (v[0] > v[1]) SWAP2(v[0],v[1],V[0],V[1],tmp);
218: }
219: return(0);
220: }
221: SWAP2(v[0],v[right/2],V[0],V[right/2],tmp);
222: vl = v[0];
223: last = 0;
224: for (i=1; i<=right; i++) {
225: if (v[i] < vl) {last++; SWAP2(v[last],v[i],V[last],V[i],tmp);}
226: }
227: SWAP2(v[0],v[last],V[0],V[last],tmp);
228: PetscSortMPIIntWithArray_Private(v,V,last-1);
229: PetscSortMPIIntWithArray_Private(v+last+1,V+last+1,right-(last+1));
230: return(0);
231: }
235: /*@
236: PetscSortMPIIntWithArray - Sorts an array of integers in place in increasing order;
237: changes a second array to match the sorted first array.
239: Not Collective
241: Input Parameters:
242: + n - number of values
243: . i - array of integers
244: - I - second array of integers
246: Level: intermediate
248: Concepts: sorting^ints with array
250: .seealso: PetscSortReal(), PetscSortIntPermutation(), PetscSortInt()
251: @*/
252: PetscErrorCode PetscSortMPIIntWithArray(PetscMPIInt n,PetscMPIInt i[],PetscMPIInt Ii[])
253: {
255: PetscMPIInt j,k,tmp,ik;
258: if (n<8) {
259: for (k=0; k<n; k++) {
260: ik = i[k];
261: for (j=k+1; j<n; j++) {
262: if (ik > i[j]) {
263: SWAP2(i[k],i[j],Ii[k],Ii[j],tmp);
264: ik = i[k];
265: }
266: }
267: }
268: } else {
269: PetscSortMPIIntWithArray_Private(i,Ii,n-1);
270: }
271: return(0);
272: }
274: /* -----------------------------------------------------------------------*/
275: #define SWAP2IntScalar(a,b,c,d,t,ts) {t=a;a=b;b=t;ts=c;c=d;d=ts;}
279: /*
280: Modified from PetscSortIntWithArray_Private().
281: */
282: static PetscErrorCode PetscSortIntWithScalarArray_Private(PetscInt *v,PetscScalar *V,PetscInt right)
283: {
285: PetscInt i,vl,last,tmp;
286: PetscScalar stmp;
289: if (right <= 1) {
290: if (right == 1) {
291: if (v[0] > v[1]) SWAP2IntScalar(v[0],v[1],V[0],V[1],tmp,stmp);
292: }
293: return(0);
294: }
295: SWAP2IntScalar(v[0],v[right/2],V[0],V[right/2],tmp,stmp);
296: vl = v[0];
297: last = 0;
298: for (i=1; i<=right; i++) {
299: if (v[i] < vl) {last++; SWAP2IntScalar(v[last],v[i],V[last],V[i],tmp,stmp);}
300: }
301: SWAP2IntScalar(v[0],v[last],V[0],V[last],tmp,stmp);
302: PetscSortIntWithScalarArray_Private(v,V,last-1);
303: PetscSortIntWithScalarArray_Private(v+last+1,V+last+1,right-(last+1));
304: return(0);
305: }
309: /*@
310: PetscSortIntWithScalarArray - Sorts an array of integers in place in increasing order;
311: changes a second SCALAR array to match the sorted first INTEGER array.
313: Not Collective
315: Input Parameters:
316: + n - number of values
317: . i - array of integers
318: - I - second array of scalars
320: Level: intermediate
322: Concepts: sorting^ints with array
324: .seealso: PetscSortReal(), PetscSortIntPermutation(), PetscSortInt(), PetscSortIntWithArray()
325: @*/
326: PetscErrorCode PetscSortIntWithScalarArray(PetscInt n,PetscInt i[],PetscScalar Ii[])
327: {
329: PetscInt j,k,tmp,ik;
330: PetscScalar stmp;
333: if (n<8) {
334: for (k=0; k<n; k++) {
335: ik = i[k];
336: for (j=k+1; j<n; j++) {
337: if (ik > i[j]) {
338: SWAP2IntScalar(i[k],i[j],Ii[k],Ii[j],tmp,stmp);
339: ik = i[k];
340: }
341: }
342: }
343: } else {
344: PetscSortIntWithScalarArray_Private(i,Ii,n-1);
345: }
346: return(0);
347: }
352: /*@
353: PetscProcessTree - Prepares tree data to be displayed graphically
355: Not Collective
357: Input Parameters:
358: + n - number of values
359: . mask - indicates those entries in the tree, location 0 is always masked
360: - parentid - indicates the parent of each entry
362: Output Parameters:
363: + Nlevels - the number of levels
364: . Level - for each node tells its level
365: . Levelcnts - the number of nodes on each level
366: . Idbylevel - a list of ids on each of the levels, first level followed by second etc
367: - Column - for each id tells its column index
369: Level: intermediate
372: .seealso: PetscSortReal(), PetscSortIntWithPermutation()
373: @*/
374: PetscErrorCode PetscProcessTree(PetscInt n,const PetscBool mask[],const PetscInt parentid[],PetscInt *Nlevels,PetscInt **Level,PetscInt **Levelcnt,PetscInt **Idbylevel,PetscInt **Column)
375: {
376: PetscInt i,j,cnt,nmask = 0,nlevels = 0,*level,*levelcnt,levelmax = 0,*workid,*workparentid,tcnt = 0,*idbylevel,*column;
378: PetscBool done = PETSC_FALSE;
381: if (!mask[0]) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Mask of 0th location must be set");
382: for (i=0; i<n; i++) {
383: if (mask[i]) continue;
384: if (parentid[i] == i) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Node labeled as own parent");
385: if (parentid[i] && mask[parentid[i]]) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Parent is masked");
386: }
389: for (i=0; i<n; i++) {
390: if (!mask[i]) nmask++;
391: }
393: /* determine the level in the tree of each node */
394: PetscMalloc(n*sizeof(PetscInt),&level);
395: PetscMemzero(level,n*sizeof(PetscInt));
396: level[0] = 1;
397: while (!done) {
398: done = PETSC_TRUE;
399: for (i=0; i<n; i++) {
400: if (mask[i]) continue;
401: if (!level[i] && level[parentid[i]]) level[i] = level[parentid[i]] + 1;
402: else if (!level[i]) done = PETSC_FALSE;
403: }
404: }
405: for (i=0; i<n; i++) {
406: level[i]--;
407: nlevels = PetscMax(nlevels,level[i]);
408: }
410: /* count the number of nodes on each level and its max */
411: PetscMalloc(nlevels*sizeof(PetscInt),&levelcnt);
412: PetscMemzero(levelcnt,nlevels*sizeof(PetscInt));
413: for (i=0; i<n; i++) {
414: if (mask[i]) continue;
415: levelcnt[level[i]-1]++;
416: }
417: for (i=0; i<nlevels;i++) {
418: levelmax = PetscMax(levelmax,levelcnt[i]);
419: }
421: /* for each level sort the ids by the parent id */
422: PetscMalloc2(levelmax,PetscInt,&workid,levelmax,PetscInt,&workparentid);
423: PetscMalloc(nmask*sizeof(PetscInt),&idbylevel);
424: for (j=1; j<=nlevels;j++) {
425: cnt = 0;
426: for (i=0; i<n; i++) {
427: if (mask[i]) continue;
428: if (level[i] != j) continue;
429: workid[cnt] = i;
430: workparentid[cnt++] = parentid[i];
431: }
432: /* PetscIntView(cnt,workparentid,0);
433: PetscIntView(cnt,workid,0);
434: PetscSortIntWithArray(cnt,workparentid,workid);
435: PetscIntView(cnt,workparentid,0);
436: PetscIntView(cnt,workid,0);*/
437: PetscMemcpy(idbylevel+tcnt,workid,cnt*sizeof(PetscInt));
438: tcnt += cnt;
439: }
440: if (tcnt != nmask) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"Inconsistent count of unmasked nodes");
441: PetscFree2(workid,workparentid);
443: /* for each node list its column */
444: PetscMalloc(n*sizeof(PetscInt),&column);
445: cnt = 0;
446: for (j=0; j<nlevels; j++) {
447: for (i=0; i<levelcnt[j]; i++) {
448: column[idbylevel[cnt++]] = i;
449: }
450: }
452: *Nlevels = nlevels;
453: *Level = level;
454: *Levelcnt = levelcnt;
455: *Idbylevel = idbylevel;
456: *Column = column;
457: return(0);
458: }