Actual source code: pheap.c
1: #include <petsc/private/petscimpl.h>
2: #include <petscviewer.h>
4: typedef struct {
5: PetscInt id;
6: PetscInt value;
7: } HeapNode;
9: struct _n_PetscHeap {
10: PetscInt end; /* one past the last item */
11: PetscInt alloc; /* length of array */
12: PetscInt stash; /* stash grows down, this points to last item */
13: HeapNode *base;
14: };
16: /*
17: The arity of the heap can be changed via the parameter B below. Consider the B=2 (arity=4 case below)
19: [00 (sentinel); 01 (min node); 10 (unused); 11 (unused); 0100 (first child); 0101; 0110; 0111; ...]
21: Slots 10 and 11 are referred to as the "hole" below in the implementation.
22: */
24: #define B 1 /* log2(ARITY) */
25: #define ARITY (1 << B) /* tree branching factor */
26: static inline PetscInt Parent(PetscInt loc)
27: {
28: PetscInt p = loc >> B;
29: if (p < ARITY) return (PetscInt)(loc != 1); /* Parent(1) is 0, otherwise fix entries ending up in the hole */
30: return p;
31: }
32: #define Value(h, loc) ((h)->base[loc].value)
33: #define Id(h, loc) ((h)->base[loc].id)
35: static inline void Swap(PetscHeap h, PetscInt loc, PetscInt loc2)
36: {
37: PetscInt id, val;
38: id = Id(h, loc);
39: val = Value(h, loc);
40: h->base[loc].id = Id(h, loc2);
41: h->base[loc].value = Value(h, loc2);
42: h->base[loc2].id = id;
43: h->base[loc2].value = val;
44: }
45: static inline PetscInt MinChild(PetscHeap h, PetscInt loc)
46: {
47: PetscInt min, chld, left, right;
48: left = loc << B;
49: right = PetscMin(left + ARITY - 1, h->end - 1);
50: chld = 0;
51: min = PETSC_INT_MAX;
52: for (; left <= right; left++) {
53: PetscInt val = Value(h, left);
54: if (val < min) {
55: min = val;
56: chld = left;
57: }
58: }
59: return chld;
60: }
62: /*@
63: PetscHeapCreate - Creates a `PetscHeap` object, a simple min-heap for `(id, value)` pairs.
65: Not Collective
67: Input Parameter:
68: . maxsize - the maximum number of items the heap can hold at once
70: Output Parameter:
71: . heap - the newly created `PetscHeap` object
73: Level: developer
75: Note:
76: The heap is ordered by `value`; items with equal values may be returned in any order.
78: .seealso: `PetscHeap`, `PetscHeapAdd()`, `PetscHeapPop()`, `PetscHeapPeek()`, `PetscHeapStash()`, `PetscHeapUnstash()`, `PetscHeapView()`, `PetscHeapDestroy()`
79: @*/
80: PetscErrorCode PetscHeapCreate(PetscInt maxsize, PetscHeap *heap)
81: {
82: PetscHeap h;
84: PetscFunctionBegin;
85: *heap = NULL;
86: PetscCall(PetscMalloc1(1, &h));
87: h->end = 1;
88: h->alloc = maxsize + ARITY; /* We waste all but one slot (loc=1) in the first ARITY slots */
89: h->stash = h->alloc;
90: PetscCall(PetscCalloc1(h->alloc, &h->base));
91: h->base[0].id = -1;
92: h->base[0].value = PETSC_INT_MIN;
93: *heap = h;
94: PetscFunctionReturn(PETSC_SUCCESS);
95: }
97: /*@
98: PetscHeapAdd - Insert an item into a `PetscHeap`.
100: Not Collective
102: Input Parameters:
103: + h - the `PetscHeap`
104: . id - the item identifier
105: - val - the value used for heap ordering
107: Level: developer
109: .seealso: `PetscHeap`, `PetscHeapCreate()`, `PetscHeapPop()`, `PetscHeapPeek()`, `PetscHeapStash()`, `PetscHeapUnstash()`, `PetscHeapDestroy()`
110: @*/
111: PetscErrorCode PetscHeapAdd(PetscHeap h, PetscInt id, PetscInt val)
112: {
113: PetscInt loc, par;
115: PetscFunctionBegin;
116: if (1 < h->end && h->end < ARITY) h->end = ARITY;
117: loc = h->end++;
118: PetscCheck(h->end <= h->stash, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Addition would exceed allocation %" PetscInt_FMT " (%" PetscInt_FMT " stashed)", h->alloc, h->alloc - h->stash);
119: h->base[loc].id = id;
120: h->base[loc].value = val;
122: /* move up until heap condition is satisfied */
123: while ((void)(par = Parent(loc)), Value(h, par) > val) {
124: Swap(h, loc, par);
125: loc = par;
126: }
127: PetscFunctionReturn(PETSC_SUCCESS);
128: }
130: /*@
131: PetscHeapPop - Remove and return the minimum item from a `PetscHeap`.
133: Not Collective
135: Input Parameter:
136: . h - the `PetscHeap`
138: Output Parameters:
139: + id - identifier of the popped item, or `-1` if the heap is empty
140: - val - value of the popped item, or `PETSC_INT_MIN` if the heap is empty
142: Level: developer
144: .seealso: `PetscHeap`, `PetscHeapCreate()`, `PetscHeapAdd()`, `PetscHeapPeek()`, `PetscHeapStash()`, `PetscHeapUnstash()`, `PetscHeapDestroy()`
145: @*/
146: PetscErrorCode PetscHeapPop(PetscHeap h, PetscInt *id, PetscInt *val)
147: {
148: PetscInt loc, chld;
150: PetscFunctionBegin;
151: if (h->end == 1) {
152: *id = h->base[0].id;
153: *val = h->base[0].value;
154: PetscFunctionReturn(PETSC_SUCCESS);
155: }
157: *id = h->base[1].id;
158: *val = h->base[1].value;
160: /* rotate last entry into first position */
161: loc = --h->end;
162: if (h->end == ARITY) h->end = 2; /* Skip over hole */
163: h->base[1].id = Id(h, loc);
164: h->base[1].value = Value(h, loc);
166: /* move down until min heap condition is satisfied */
167: loc = 1;
168: while ((chld = MinChild(h, loc)) && Value(h, loc) > Value(h, chld)) {
169: Swap(h, loc, chld);
170: loc = chld;
171: }
172: PetscFunctionReturn(PETSC_SUCCESS);
173: }
175: /*@
176: PetscHeapPeek - Return the minimum item of a `PetscHeap` without removing it.
178: Not Collective
180: Input Parameter:
181: . h - the `PetscHeap`
183: Output Parameters:
184: + id - identifier of the minimum item, or `-1` if the heap is empty
185: - val - value of the minimum item, or `PETSC_INT_MIN` if the heap is empty
187: Level: developer
189: .seealso: `PetscHeap`, `PetscHeapCreate()`, `PetscHeapAdd()`, `PetscHeapPop()`, `PetscHeapDestroy()`
190: @*/
191: PetscErrorCode PetscHeapPeek(PetscHeap h, PetscInt *id, PetscInt *val)
192: {
193: PetscFunctionBegin;
194: if (h->end == 1) {
195: *id = h->base[0].id;
196: *val = h->base[0].value;
197: PetscFunctionReturn(PETSC_SUCCESS);
198: }
200: *id = h->base[1].id;
201: *val = h->base[1].value;
202: PetscFunctionReturn(PETSC_SUCCESS);
203: }
205: /*@
206: PetscHeapStash - Set aside an item in a `PetscHeap` for later insertion via `PetscHeapUnstash()`.
208: Not Collective
210: Input Parameters:
211: + h - the `PetscHeap`
212: . id - the item identifier
213: - val - the value used for heap ordering
215: Level: developer
217: Note:
218: Stashed items are held in the trailing portion of the heap's storage and do not participate in
219: heap ordering until `PetscHeapUnstash()` reinserts them.
221: .seealso: `PetscHeap`, `PetscHeapCreate()`, `PetscHeapAdd()`, `PetscHeapUnstash()`, `PetscHeapDestroy()`
222: @*/
223: PetscErrorCode PetscHeapStash(PetscHeap h, PetscInt id, PetscInt val)
224: {
225: PetscInt loc;
227: PetscFunctionBegin;
228: loc = --h->stash;
229: h->base[loc].id = id;
230: h->base[loc].value = val;
231: PetscFunctionReturn(PETSC_SUCCESS);
232: }
234: /*@
235: PetscHeapUnstash - Reinsert all items previously stashed with `PetscHeapStash()` into the heap.
237: Not Collective
239: Input Parameter:
240: . h - the `PetscHeap`
242: Level: developer
244: .seealso: `PetscHeap`, `PetscHeapCreate()`, `PetscHeapAdd()`, `PetscHeapStash()`, `PetscHeapDestroy()`
245: @*/
246: PetscErrorCode PetscHeapUnstash(PetscHeap h)
247: {
248: PetscFunctionBegin;
249: while (h->stash < h->alloc) {
250: PetscInt id = Id(h, h->stash), value = Value(h, h->stash);
251: h->stash++;
252: PetscCall(PetscHeapAdd(h, id, value));
253: }
254: PetscFunctionReturn(PETSC_SUCCESS);
255: }
257: /*@
258: PetscHeapDestroy - Destroys a `PetscHeap` created with `PetscHeapCreate()`.
260: Not Collective
262: Input Parameter:
263: . heap - the `PetscHeap` to destroy; set to `NULL` on return
265: Level: developer
267: .seealso: `PetscHeap`, `PetscHeapCreate()`
268: @*/
269: PetscErrorCode PetscHeapDestroy(PetscHeap *heap)
270: {
271: PetscFunctionBegin;
272: PetscCall(PetscFree((*heap)->base));
273: PetscCall(PetscFree(*heap));
274: PetscFunctionReturn(PETSC_SUCCESS);
275: }
277: /*@
278: PetscHeapView - View the contents of a `PetscHeap`, including any stashed items.
280: Not Collective
282: Input Parameters:
283: + h - the `PetscHeap`
284: - viewer - a `PetscViewer`, or `NULL` to use `PETSC_VIEWER_STDOUT_SELF`
286: Level: developer
288: .seealso: `PetscHeap`, `PetscHeapCreate()`, `PetscHeapAdd()`, `PetscHeapPop()`
289: @*/
290: PetscErrorCode PetscHeapView(PetscHeap h, PetscViewer viewer)
291: {
292: PetscBool isascii;
294: PetscFunctionBegin;
295: if (!viewer) PetscCall(PetscViewerASCIIGetStdout(PETSC_COMM_SELF, &viewer));
297: PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &isascii));
298: if (isascii) {
299: PetscCall(PetscViewerASCIIPrintf(viewer, "Heap size %" PetscInt_FMT " with %" PetscInt_FMT " stashed\n", h->end - 1, h->alloc - h->stash));
300: PetscCall(PetscViewerASCIIPrintf(viewer, "Heap in (id,value) pairs\n"));
301: PetscCall(PetscIntView(2 * (h->end - 1), (const PetscInt *)(h->base + 1), viewer));
302: PetscCall(PetscViewerASCIIPrintf(viewer, "Stash in (id,value) pairs\n"));
303: PetscCall(PetscIntView(2 * (h->alloc - h->stash), (const PetscInt *)(h->base + h->stash), viewer));
304: }
305: PetscFunctionReturn(PETSC_SUCCESS);
306: }