Actual source code: pheap.c
petsc-3.13.6 2020-09-29
2: #include <petsc/private/petscimpl.h>
3: #include <petscviewer.h>
5: typedef struct {
6: PetscInt id;
7: PetscInt value;
8: } HeapNode;
10: struct _PetscHeap {
11: PetscInt end; /* one past the last item */
12: PetscInt alloc; /* length of array */
13: PetscInt stash; /* stash grows down, this points to last item */
14: HeapNode *base;
15: };
17: /*
18: * The arity of the heap can be changed via the parameter B below. Consider the B=2 (arity=4 case below)
19: *
20: * [00 (sentinel); 01 (min node); 10 (unused); 11 (unused); 0100 (first child); 0101; 0110; 0111; ...]
21: *
22: * Slots 10 and 11 are referred to as the "hole" below in the implementation.
23: */
25: #define B 1 /* log2(ARITY) */
26: #define ARITY (1 << B) /* tree branching factor */
27: PETSC_STATIC_INLINE PetscInt Parent(PetscInt loc)
28: {
29: PetscInt p = loc >> B;
30: if (p < ARITY) return (PetscInt)(loc != 1); /* Parent(1) is 0, otherwise fix entries ending up in the hole */
31: return p;
32: }
33: #define Value(h,loc) ((h)->base[loc].value)
34: #define Id(h,loc) ((h)->base[loc].id)
36: PETSC_STATIC_INLINE void Swap(PetscHeap h,PetscInt loc,PetscInt loc2)
37: {
38: PetscInt id,val;
39: id = Id(h,loc);
40: val = Value(h,loc);
41: h->base[loc].id = Id(h,loc2);
42: h->base[loc].value = Value(h,loc2);
43: h->base[loc2].id = id;
44: h->base[loc2].value = val;
45: }
46: PETSC_STATIC_INLINE PetscInt MinChild(PetscHeap h,PetscInt loc)
47: {
48: PetscInt min,chld,left,right;
49: left = loc << B;
50: right = PetscMin(left+ARITY-1,h->end-1);
51: chld = 0;
52: min = PETSC_MAX_INT;
53: for (; left <= right; left++) {
54: PetscInt val = Value(h,left);
55: if (val < min) {
56: min = val;
57: chld = left;
58: }
59: }
60: return chld;
61: }
63: PetscErrorCode PetscHeapCreate(PetscInt maxsize,PetscHeap *heap)
64: {
66: PetscHeap h;
69: *heap = NULL;
70: PetscMalloc1(1,&h);
71: h->end = 1;
72: h->alloc = maxsize+ARITY; /* We waste all but one slot (loc=1) in the first ARITY slots */
73: h->stash = h->alloc;
74: PetscCalloc1(h->alloc,&h->base);
75: h->base[0].id = -1;
76: h->base[0].value = PETSC_MIN_INT;
77: *heap = h;
78: return(0);
79: }
81: PetscErrorCode PetscHeapAdd(PetscHeap h,PetscInt id,PetscInt val)
82: {
83: PetscInt loc,par;
86: if (1 < h->end && h->end < ARITY) h->end = ARITY;
87: loc = h->end++;
88: if (h->end > h->stash) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Addition would exceed allocation %D (%D stashed)",h->alloc,(h->alloc-h->stash));
89: h->base[loc].id = id;
90: h->base[loc].value = val;
92: /* move up until heap condition is satisfied */
93: while ((void)(par = Parent(loc)), Value(h,par) > val) {
94: Swap(h,loc,par);
95: loc = par;
96: }
97: return(0);
98: }
100: PetscErrorCode PetscHeapPop(PetscHeap h,PetscInt *id,PetscInt *val)
101: {
102: PetscInt loc,chld;
105: if (h->end == 1) {
106: *id = h->base[0].id;
107: *val = h->base[0].value;
108: return(0);
109: }
111: *id = h->base[1].id;
112: *val = h->base[1].value;
114: /* rotate last entry into first position */
115: loc = --h->end;
116: if (h->end == ARITY) h->end = 2; /* Skip over hole */
117: h->base[1].id = Id(h,loc);
118: h->base[1].value = Value(h,loc);
120: /* move down until min heap condition is satisfied */
121: loc = 1;
122: while ((chld = MinChild(h,loc)) && Value(h,loc) > Value(h,chld)) {
123: Swap(h,loc,chld);
124: loc = chld;
125: }
126: return(0);
127: }
129: PetscErrorCode PetscHeapPeek(PetscHeap h,PetscInt *id,PetscInt *val)
130: {
132: if (h->end == 1) {
133: *id = h->base[0].id;
134: *val = h->base[0].value;
135: return(0);
136: }
138: *id = h->base[1].id;
139: *val = h->base[1].value;
140: return(0);
141: }
143: PetscErrorCode PetscHeapStash(PetscHeap h,PetscInt id,PetscInt val)
144: {
145: PetscInt loc;
148: loc = --h->stash;
149: h->base[loc].id = id;
150: h->base[loc].value = val;
151: return(0);
152: }
154: PetscErrorCode PetscHeapUnstash(PetscHeap h)
155: {
159: while (h->stash < h->alloc) {
160: PetscInt id = Id(h,h->stash),value = Value(h,h->stash);
161: h->stash++;
162: PetscHeapAdd(h,id,value);
163: }
164: return(0);
165: }
167: PetscErrorCode PetscHeapDestroy(PetscHeap *heap)
168: {
172: PetscFree((*heap)->base);
173: PetscFree(*heap);
174: return(0);
175: }
177: PetscErrorCode PetscHeapView(PetscHeap h,PetscViewer viewer)
178: {
180: PetscBool iascii;
183: if (!viewer) {
184: PetscViewerASCIIGetStdout(PETSC_COMM_SELF,&viewer);
185: }
187: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);
188: if (iascii) {
189: PetscViewerASCIIPrintf(viewer,"Heap size %D with %D stashed\n",h->end-1,h->alloc-h->stash);
190: PetscViewerASCIIPrintf(viewer,"Heap in (id,value) pairs\n");
191: PetscIntView(2*(h->end-1),(const PetscInt*)(h->base+1),viewer);
192: PetscViewerASCIIPrintf(viewer,"Stash in (id,value) pairs\n");
193: PetscIntView(2*(h->alloc-h->stash),(const PetscInt*)(h->base+h->stash),viewer);
194: }
195: return(0);
196: }