I have made a generic pairing heap library in C. Pairing heaps are one of the several heap variants with better asymptotic running times than standard binary heaps (others include Fibonacci heaps and binomial heaps). However, pairing heaps are the only ones that really do better than binary heaps according to Wikipedia.
A pairing heap is either empty or consists of an element (the minimum) and a list of non empty pairing heaps. (The name "pairing" is not derived from this definition whatsoever and actually refers to the normal implementation of the pop operation, which I have improved.)
I have a struct, pheap_ft, that holds function pointers and state data for allocating, comparing, and freeing nodes and so forth. The pheap functions then support new, top, push, pop, meld, and delete. It is often also useful to provide an update key operation for a heap, but since my nodes do not have parent pointers this is not practical (and anyway heaps do not support finding an element by value efficiently; an avl tree would probably suit such usage better).
new and top are trivial. meld just adds the input heap with the larger root (minimum) element as the first child of the input heap with the smaller root. push is just new then meld. The only complicated operation is pop (typically, update key would be implemented as popping a non root node and then pushing it back onto the heap, but again I have omitted this because I do not think it is worth adding parent pointers).
pop removes the root (minimum) element and then merges its children. The standard merge procedure pairs off the the child heaps and merges them in twos, creating a queue, the elements of which are then merged into each other iteratively. I have instead used a variation of merge sort that does not need to know the input length. I merge the first 2^0 elements, then 2^0, then 2^1, then 2^2, then 2^3, and so on, as long as there are elements. Calls to merge 2^(n+1) will call merge 2^n twice, and merge 2^0 just splits off one child from the list (a single element is trivially sorted). Once all of the child heaps have been merged, no more recursion will occur even if the 2^n is still not 2^0. This could be further improved using natural merge sort which keeps monotonic runs that happen to occur.
I have made the pairing heap generic by using void pointers for parameter and return types, but the actual nodes use flexible length char arrays, which has a couple of problems I will discuss after the code.
Here is the main header file pheap.h
//pheap.h
//Pairing heap library by Hacatu
#ifndef __PHEAP_H__
#define __PHEAP_H__
#include <stddef.h>
#include "fla.h"
typedef struct pheap_node pheap_node;
struct pheap_node{
pheap_node *first_child, *sibling;
char data[];
};
typedef struct{//function table to hold function pointers and related data
size_t size;
int (*cmp)(const void *data, const void*, const void*);//some cmps are unsequenced so data really must not be modified here
void *(*alloc)(void *data);
void (*free)(void *data, void*);
void *data;
} pheap_ft;
//Allocate a new pairing heap node using ft->alloc(ft->data), copy ft->size bytes from key to
//(return)->data, set (return)->first_child and (return)->sibling, and (return) the new node.
pheap_node *pheap_new(void *key, pheap_node *first_child, pheap_node *sibling, pheap_ft *ft);
//Return a pointer to the minimal data element of a pairing heap, or NULL if the heap is empty
void *pheap_top(pheap_node *r);
//Combine two pairing heaps. The nodes of the existing heaps are reused and re-linked so no new
//allocations are done, but the original heaps are mutated. Returns a pointer to the root element
//of the melded heap
pheap_node *pheap_meld(pheap_node *a, pheap_node *b, pheap_ft *ft);
//Insert an item into a pairing heap, returning the new root element on success or NULL on failure.
//Failure only occurs if allocation fails: items with duplicate keys are permitted with unspecified
//storage order
pheap_node *pheap_push(pheap_node *r, void *key, pheap_ft *ft);
//Remove the root (minimal) element from a pairing heap. The input pointer is modified to indicate the
//new root. The old root is returned, but its links may not be nulled. If pop is called on an empty heap,
//nothing is done and NULL is returned.
pheap_node *pheap_pop(pheap_node **r, pheap_ft *ft);
void pheap_delete(pheap_node *r, pheap_ft *ft);
//macro for obtaining data field (without using a function like pheap_top) despite the strict aliasing
//rule. The pairing heap nodes store data in a flexible length char array because there is no such
//thing as a flexible length void array. This means that the data has a type of char and so to read it
//as another type without optimization messing it up (since pointers to different types "can't" refer
//to the same object) we have to use a union. This macro does that, but it requires the __typeof__ extension.
#define PHEAP_DATA(T, n) FLA_CAST(T, (char*)((pheap_node*)(n))->data)
#endif
Here is the header file fla.h. It contains a macro for working around a particular limitation of flexible length arrays:
//fla.h
#ifndef __FLA_H__
#define __FLA_H__
#include <stddef.h>
#define FLA_CAST(T, p) (((union{__typeof__(p) data; T a;})(p)).a)
#endif
And here is the source code (pheap.c):
//pheap.c
#include <stdlib.h>
#include <string.h>
#include "pheap.h"
void *pheap_top(pheap_node *r){
return r ? r->data : NULL;
}
pheap_node *pheap_new(void *key, pheap_node *first_child, pheap_node *sibling, pheap_ft *ft){
pheap_node *n = ft->alloc(ft->data);
if(n){
*n = (pheap_node){.first_child=first_child, .sibling=sibling};
memcpy(n->data, key, ft->size);
}
return n;
}
pheap_node *pheap_meld(pheap_node *a, pheap_node *b, pheap_ft *ft){
if(!a){
return b;
}else if(!b){
return a;
}else if(ft->cmp(ft->data, a->data, b->data) > 0){
a->sibling = b->first_child;
b->first_child = a;
return b;
}
b->sibling = a->first_child;
a->first_child = b;
return a;
}
pheap_node *pheap_push(pheap_node *r, void *key, pheap_ft *ft){
pheap_node *n = pheap_new(key, NULL, NULL, ft);
return n ? pheap_meld(r, n, ft) : NULL;
}
//helper function for pheap_pop. Uses mergesort to combine a pheap_node's children
//linked list into one pairing heap. Since the length of the list of children is
//not stored, this function accepts a recursion depth and is called by another helper
//function (pheap_merge_exponential) with depth 0, 0, 1, 2, 3, 4, 5, 6, ... until all
//children are merged. If the list is empty, nothing is done. If the depth is 0,
//the first child in the remaining part of the list is removed and returned (1 node is sorted).
static pheap_node *pheap_merge_binary(pheap_node **n, size_t depth, pheap_ft *ft){
pheap_node *a, *b;
if(!*n){
return NULL;
}else if(!depth){
a = *n;
*n = (*n)->sibling;
a->sibling = NULL;
return a;
}
a = pheap_merge_binary(n, depth - 1, ft);
b = pheap_merge_binary(n, depth - 1, ft);
return pheap_meld(a, b, ft);
}
//helper function for pheap_pop. Uses mergesort to combine a pheap_node's children
//linked list into one pairing heap.
static pheap_node *pheap_merge_exponential(pheap_node **n, pheap_ft *ft){
pheap_node *a = pheap_merge_binary(n, 0, ft);
for(size_t depth = 0; *n; ++depth){
a = pheap_meld(a, pheap_merge_binary(n, depth, ft), ft);
}
return a;
}
pheap_node *pheap_pop(pheap_node **r, pheap_ft *ft){
if(!*r){
return NULL;
}
pheap_node *ret = *r;
*r = pheap_merge_exponential(&(*r)->first_child, ft);
return ret;
}
//attempts to be as iterative as possible but is still singly recursive.
void pheap_delete(pheap_node *r, pheap_ft *ft){
for(pheap_node *t; r;){
while(r->sibling){
pheap_delete(r->sibling->first_child, ft);
t = r->sibling;
r->sibling = t->sibling;
ft->free(ft->data, t);
}
t = r;
r = r->first_child;
ft->free(ft->data, t);
}
}
I used a heap based sieve of Eratosthenes and Dijkstra's algorithm (specifically Project Euler 83) to test the library. Both of these tests use a slab allocator and can be found in this gist. The prime number sieve is actually 6% slower than one using a binary heap, but that is probably because it does not use any of the operations that are faster.
The problem with using flexible length char arrays to hold data is that C has the "strict aliasing rule" which says that an object of type a and one of type b cannot have the same address, so basically just casting the char[] to a struct foo* will be optimized to nonsense. The workaround is to use a union since union members are exempt from the strict aliasing rule. I cast through a union using a macro since I can't know the user's type. However, the char[] is not very aligned and because of this there could be alignment issues on some platforms. I am not sure if this is an important issue. I do not want to use an inverted scheme and containerof because it requires extra type definitions.
There are a few aspects of pheap.h, pheap.c, and fla.h that I would like reviewed. Is the approach that I used to make the heap generic good? Are the flexible length char arrays a good idea and if not how should I replace them? Can I make this faster? Thank you.
