PTHash: Minimal Perfect Hashing

Imagine we have $n$ distinct keys, known in advance. A perfect hash function (PHF) can map them into a range of size about $n$ with zero collisions. A minimal perfect hash function (MPHF) goes further and maps them onto exactly $0..n-1$.

In this article we'll be building a perfect hash function that uses about 2 bits per key on average, no matter how long the keys are. A billion keys (URLs, for example) can fit into a fraction of a gigabyte in a data structure that can answer the question "which slot is this key's data in?".

The only catch is that the function only promises anything for the $n$ keys it was built on. A key outside the set still returns some slot number, without any error. A perfect hash function doesn't store the actual keys (obviously), so it can't tell whether a key was in the original set. Any use of such a structure should either tolerate that or keep the keys elsewhere.

I implemented PTHash from the paper by Giulio Ermanno Pibiri and Roberto Trani, PTHash: Revisiting FCH Minimal Perfect Hashing, plus one piece of the follow-up PHOBIC paper.

This article walks through my pthash.zig and builds a working toy version along the way. We'll also use the data structures from the Elias-Fano and Rank-Select articles to achieve the required compression levels.

General idea

A whole PTHash query fits in five lines. First we hash the key, then pick a bucket, then read that bucket's pilot (a small integer found at construction time), then hash the pilot, xor and reduce:

const hash = hashKey(key);
const bucket = getBucket(hash);
const pilot = pilots[bucket];
const position = (hash ^ hashPilot(pilot)) % table_size;
Copy

For a regular PHF, the structure stores only the pilots array and a few numbers like the seed and the table size. For MPHF we also need to store a small array of rearranged positions, which gives a slightly worse bits-per-key ratio. The keys themselves are not stored anywhere.

Construction is the problem of finding pilots that make the table collision-free. Keys land in buckets, each bucket has a set of assigned keys defined by the seed, and all keys of a bucket share one pilot. We can change the pilot and re-hash the whole bucket to completely new positions, since the pilot's hash gets xor-ed with the key hash. Our task is to find pilots that produce no collisions in the table. Each attempt places the bucket's $s$ keys at new independent random positions, and we keep trying until all of them land in free slots.

That's the general idea of the algorithm. Everything else in the paper is about making the search fast and the pilots small.

Construction

Buckets

A key's bucket must depend only on its hash, because the query has to recompute it. We'll get to that in detail below, but for now let's keep it simple and define our bucket assignment function like this.

fn getBucket(hash: u64) u64 {
  return hash % num_buckets;
}
Copy

Building the structure starts by hashing every key. We'll use Wyhash for string keys, then sort the (bucket, hash) pairs by bucket size. We need to sort to process the largest buckets first to simplify pilot search: smaller buckets are usually easier to map.

The next step is to find pilots for each bucket so the mapped key positions don't collide.

Pilot Search

Pilot search in our case is relatively simple. We process each bucket in order of decreasing size to leave the easiest bucket for the end, when the table is almost full. For each bucket, we're looking for a pilot that makes all the keys in this bucket map to positions not yet occupied in the table. If the current pilot causes collisions, we just increment it and try the next one.

Here is the search algorithm:

pilot_search: while (true) : (pilot += 1) {
    const hashed_pilot = hashPilot(seed, pilot);

    positions.clearRetainingCapacity();
    for (hashes[bucket.i_start..bucket.i_end]) |it| {
        const p = (it.hash ^ hashed_pilot) % table_size;
        if (table[p]) continue :pilot_search;
        positions.appendAssumeCapacity(p);
    }

    std.mem.sort(usize, positions.items, {}, std.sort.asc(usize));
    for (1..positions.items.len) |j| {
        if (positions.items[j - 1] == positions.items[j]) continue :pilot_search;
    }

    pilots[bucket.bucket] = pilot;
    for (positions.items) |p| table[p] = true;
    break;
}
Copy

An attempt can fail in one of two ways: a key hits a slot taken by an earlier bucket, or two keys of the same bucket hit the same free slot, which is why we sort the positions and scan them for duplicates. On collision the pilot increments and everything is tried again.

Here's the simplified runnable version without Elias-Fano compression.

const std = @import("std");
const print = std.debug.print;

fn hashKey(key: []const u8) u64 {
    return std.hash.Wyhash.hash(0x42B, key);
}

fn hashPilot(pilot: u64) u64 {
    return std.hash.murmur.Murmur2_64.hashWithSeed(std.mem.asBytes(&pilot), 0x42C);
}

const keys = [_][]const u8{
    "mercury", "venus", "earth", "mars", "jupiter", "saturn",
    "uranus",  "neptune", "pluto", "iam", "out", "of", "planets"
};

fn getBucket(hash: u64, num_buckets: u64) u64 {
  return hash % num_buckets;
}

pub fn main() !void {
    const avg_bucket_size = 3;
    const table_size = (keys.len * 10) / 9 + 1;
    const num_buckets = table_size / avg_bucket_size;
    
    print("table_size={} num_keys={} num_buckets={}\n\n", .{ table_size, keys.len, num_buckets });
    
    var hashes: [keys.len]u64 = undefined;
    var bucket_of: [keys.len]u64 = undefined;
    for (keys, 0..) |key, i| {
        hashes[i] = hashKey(key);
        bucket_of[i] = getBucket(hashes[i], num_buckets);
        print("{s: <9} hash=0x{X:0>16} bucket={}\n", .{ key, hashes[i], bucket_of[i] });
    }

    var sizes = [_]u64{0} ** num_buckets;
    for (bucket_of) |b| sizes[b] += 1;

    var order: [num_buckets]usize = undefined;
    for (0..num_buckets) |i| order[i] = i;
    std.mem.sort(usize, &order, &sizes, struct {
        fn func(s: *const [num_buckets]u64, a: usize, b: usize) bool {
            return s[a] > s[b];
        }
    }.func);

    print("\n", .{});

    var taken = [_]u64{0} ** (table_size / 64 + 1); // bitmap to save some space
    var pilots = [_]u64{0} ** num_buckets;
    for (order) |b| {
        if (sizes[b] == 0) continue;
        var pilot: u64 = 0;
        search: while (true) : (pilot += 1) {
            const hp = hashPilot(pilot);
            var positions: [table_size]u64 = undefined;
            var n: usize = 0;
            
            for (hashes, bucket_of) |h, hb| {
                if (hb != b) continue;
                const p = (h ^ hp) % table_size;
                if (taken[p / 64] & (@as(u64, 1) << @intCast(p % 64)) > 0) continue :search;
                for (positions[0..n]) |q| if (q == p) continue :search;
                positions[n] = p;
                n += 1;
            }
            
            for (positions[0..n]) |p| taken[p / 64] |= @as(u64, 1) << @intCast(p % 64);
            
            pilots[b] = pilot;
            print("bucket {} (size {}): pilot = {: <3} ({} attempts)\n", .{ b, sizes[b], pilot, pilot + 1 });
            break;
        }
    }

    print("\n", .{});
    for (keys) |key| {
        const h = hashKey(key);
        const p = (h ^ hashPilot(pilots[getBucket(h, num_buckets)])) % table_size;
        print("{s: <9} -> {: <2}\n", .{ key, p });
    }
}
Copy

Every key gets a distinct slot in 0..15. We allocate slightly more slots than we have keys, otherwise, pilot search can take much longer to complete, and the pilots will be much bigger. That is bad for compression, so by sacrificing a few additional slots, we can still win during compression because the pilots become smaller.

Making It Minimal

With $\alpha < 1$ the function maps $n$ keys into about $n/\alpha$ slots, and some keys land at positions $\ge n$. To get a bijection onto $0..n-1$ those keys must be redirected into the free slots below $n$, and there are exactly as many free slots below $n$ as there are occupied positions above it.

To do that, we walk the overflow positions $n, n+1, \dots$, and for each occupied one we record the next unused free slot below $n$. Lookups then remap in one step:

if (p >= n) p = free_slots[p - n];
Copy

If stored raw that's a u64 per overflow position. The good thing about these array values is that the values are sorted which will allow us to use compress it. Overflow positions that aren't occupied are never queried, so we can store anything there. We repeat the previous value and the whole array stays non-decreasing to allow us to use Elias-Fano encoding.

// ...

const n = keys.len;
var free_slots: [table_size - n]u64 = undefined;
var next_free: u64 = 0;
var prev: u64 = 0;
for (0..free_slots.len) |i| {
    if (taken[n + i]) {
        while (taken[next_free]) next_free += 1;
        free_slots[i] = next_free;
        prev = next_free;
        next_free += 1;
    } else {
        // repeat the previous value to keep the sequence monotone for EliasFano
        free_slots[i] = prev;
    }
}

const fs = try EliasFano.init(allocator, free_slots);
Copy

Storing Pilots

The way how we store pilots affects the bits per key ratio more than anything else. Pilots are small non-negative integers with a skewed distribution, most are small and a few are huge, and they need constant-time random access. The PTHash paper benchmarks a whole family of encodings, from plain fixed-width arrays to dictionary coding to Elias-Fano.

I only implemented the Elias-Fano option since it provides a good balance between performance and space taken. Pilots aren't monotone, so to make it work with Elias-Fano encoding we need to store their prefix sums. In this case lookup takes the difference of two adjacent values.

To get the result from stored prefix sums, we need to do one select (rank-select article) on the EF upper bits, followed by a quick forward scan for the next one-bit.

Parameters

An attempt succeeds when all $m$ keys of the bucket fall into free slots and don't pair up. With a fraction $f$ of the table already taken, each key survives with probability about $1 - t$, so an attempt succeeds with probability about $(1-f)^m$ and the expected number of attempts is $(1-f)^{-m}$, exponential in the bucket size. Everything about PTHash's tuning follows from this one expression.

$\lambda$ (lambda) is the average bucket size, $\lambda = n / \text{num\_buckets}$. Raising it makes every bucket bigger, which increases $(1-f)^{-s}$ and with it the build time, but there are fewer pilots to store, so that's the space win.

$\alpha$ (alpha) is the load factor. The table gets $\left\lfloor \frac{n}{\alpha} + 1 \right\rfloor$ slots, so at $\alpha = 0.97$ about 3% of slots are never filled. Those spare slots cap how bad $f$ gets for the last buckets. With $\alpha = 1$ the final single-key buckets are throwing keys at a table with one free slot, and the expected number of attempts to find the pilot will be around $n$. A few percent of spare slots remove that worst part of the search. In exchange, a perfect hash now maps into $n/\alpha$ slots instead of $n$, which is no longer minimal, but we'll fix that with the free slots remapping.

To get some real numbers I ran 30k random keys through the toy search:

const std = @import("std");
const print = std.debug.print;

const BucketHash = struct { bucket: u64, hash: u64 };

const Result = struct {
    attempts: u64,
    max_pilot: u64,
    pilot_bits_per_key: f64,
    remap_bits_per_key: f64,
};

fn hashPilot(pilot: u64) u64 {
    return std.hash.murmur.Murmur2_64.hashWithSeed(std.mem.asBytes(&pilot), 0xDEAD);
}

fn eliasFanoSizeInBits(n_values: f64, universe: f64) f64 {
    // N * (ceil(log2(U / N)) + 2), the Elias-Fano size formula
    const low = @max(0, @ceil(@log2(universe / n_values)));
    return n_values * (low + 2);
}


fn buildStats(allocator: std.mem.Allocator, hashes: []BucketHash, num_buckets: u64, table_size: u64) !Result {
    const n = hashes.len;

    const p1: u64 = @max(1, num_buckets * 3 / 10);
    const p1_top: u64 = @intFromFloat(0.6 * @as(f64, @floatFromInt(std.math.maxInt(u64))));
    for (hashes) |*bh| {
        bh.bucket = if (bh.hash < p1_top) bh.hash % p1 else p1 + bh.hash % (num_buckets - p1);
    }

    std.mem.sort(BucketHash, hashes, {}, struct {
        fn func(_: void, a: BucketHash, b: BucketHash) bool {
            return if (a.bucket == b.bucket) a.hash < b.hash else a.bucket < b.bucket;
        }
    }.func);

    var groups: std.ArrayList(struct { start: usize, len: usize }) = .empty;
    defer groups.deinit(allocator);
    var start: usize = 0;
    for (1..n + 1) |i| {
        if (i == n or hashes[i].bucket != hashes[i - 1].bucket) {
            try groups.append(allocator, .{ .start = start, .len = i - start });
            start = i;
        }
    }
    std.mem.sort(@TypeOf(groups.items[0]), groups.items, {}, struct {
        fn func(_: void, a: @TypeOf(groups.items[0]), b: @TypeOf(groups.items[0])) bool {
            return a.len > b.len;
        }
    }.func);

    var taken = try allocator.alloc(u64, (table_size / 64 + 1));
    defer allocator.free(taken);
    @memset(taken, 0);

    var positions: std.ArrayList(u64) = .empty;
    defer positions.deinit(allocator);

    var attempts: u64 = 0;
    var max_pilot: u64 = 0;
    var pilot_sum: f64 = 0;
    var free_count: f64 = 0;
    for (groups.items) |g| {
        var pilot: u64 = 0;
        search: while (true) : (pilot += 1) {
            attempts += 1;
            const hp = hashPilot(pilot);
            
            positions.clearRetainingCapacity();
            for (hashes[g.start..g.start + g.len]) |bh| {
                const p = (bh.hash ^ hp) % table_size;
                if (taken[p / 64] & (@as(u64, 1) << @intCast(p % 64)) > 0) continue :search;
                for (positions.items) |q| if (q == p) continue :search;
                try positions.append(allocator, p);
            }
            
            for (positions.items) |p| {
              taken[p / 64] |= @as(u64, 1) << @intCast(p % 64);
              if (p >= n) free_count += 1;
            }
            
            max_pilot = @max(max_pilot, pilot);
            pilot_sum += @floatFromInt(pilot);
            break;
        }
    }

    const nf: f64 = @floatFromInt(n);
    const nb: f64 = @floatFromInt(num_buckets);
    return .{
        .attempts = attempts,
        .max_pilot = max_pilot,
        .pilot_bits_per_key = eliasFanoSizeInBits(nb + 1, pilot_sum) / nf,
        .remap_bits_per_key = eliasFanoSizeInBits(free_count, nf) / nf,
    };
}


pub fn main() !void {
    var arena = std.heap.ArenaAllocator.init(std.heap.page_allocator);
    defer arena.deinit();
    const allocator = arena.allocator();

    const n = 30_000;
    var prng = std.Random.DefaultPrng.init(0x2BAD2D0);
    const hashes = try allocator.alloc(BucketHash, n);
    for (hashes) |*bh| bh.hash = prng.random().int(u64);

    print("lambda  alpha   attempts/key  max pilot  bits/key PHF  bits/key MPHF\n", .{});
    const configs = [_]struct { lambda: f64, alpha: f64 }{
        .{ .lambda = 2, .alpha = 0.97 },
        .{ .lambda = 4, .alpha = 0.97 },
        .{ .lambda = 6, .alpha = 0.97 },
        .{ .lambda = 8, .alpha = 0.97 },
        .{ .lambda = 6, .alpha = 0.90 },
        .{ .lambda = 6, .alpha = 0.99 },
    };
    
    for (configs) |cfg| {
        const num_buckets: u64 = @intFromFloat(@as(f64, n) / cfg.lambda);
        const table_size: u64 = @intFromFloat(@as(f64, @floatFromInt(n)) / cfg.alpha + 1);
        const r = try buildStats(allocator, hashes, num_buckets, table_size);
        print("{d: <7} {d: <7} {d: <13.1} {d: <10} {d: <13.2} {d: <14.2}\n", .{
            cfg.lambda,
            cfg.alpha,
            @as(f64, @floatFromInt(r.attempts)) / @as(f64, n), 
            r.max_pilot,
            r.pilot_bits_per_key,
            r.pilot_bits_per_key + r.remap_bits_per_key,
        });
    }
}
Copy

This prints:

lambda  alpha   attempts/key  max pilot  bits/key PHF  bits/key MPHF
2       0.97    4.7           217        3.00          3.24          
4       0.97    16.0          1068       2.00          2.24          
6       0.97    74.6          7273       1.83          2.08          
8       0.97    365.6         45885      1.75          1.99          
6       0.9     28.1          2825       1.67          2.27          
6       0.99    109.8         12079      2.00          2.09
Copy

Going from $\lambda = 2$ to $\lambda = 8$ cuts the pilot storage from 3 bits to 1.75 bits per key while the search does 80x more work. Each step of 2 in $\lambda$ multiplies the attempts by 3-5x, as the exponential predicts. The $\alpha$ represents the table load factor. Dropping to 0.9 makes pilots cheaper at 1.67 bits, but the free-slot remapping eats the savings back. Setting it to 0.99 shrinks the remap to almost nothing while the pilots grow causing EliasFano size to grow.

Bucket Assignment Functions

Early bucket face almost empty table and late buckets face a full one. With uniform assignment all buckets have the same expected size, so the searches are cheap early and very slow at the end. Both PTHash an PHOBIC papers exploit the same observation. Big buckets are cheap to place early and small ones are cheap late, so bucket sizes should decrease over the course of the build. Keeping all buckets small would greatly hurt a bits-per-key ratio so bucket sizes should be dynamic.

In the original PTHash paper this is one hardcoded step. 60% of the hash space goes to the first 30% of buckets and the remaining 40% to the other 70%. So the first 30% of buckets come out twice as big on average, and since the build sorts by size, they mostly go first.

PHOBIC authors replace the step with the function derived to be optimal for this cost model:

$x$ is a hash scaled to $(0, 1)$. Since $\beta$ grows slowly near $0$ and steeply near $1$. The first buckets receive wide ranges of hash space and get big, and the bucket sizes smoothly decrease near to the end instead of in one step. In practice $\beta$ gets blended with the uniform mapping, $\beta_\varepsilon(x) = \varepsilon x + (1 - \varepsilon)\beta(x)$, because pure $\beta$ makes the first buckets too large. In the paper authors calculate $\varepsilon$ using the following formula $\varepsilon = \lambda / (5 * sqrt(P))$ (heuristic) where $P$ is a partition size. We'll use $num\_keys$ instead of $P$ since we don't have partitioning.

Both mappers, plus uniform as the baseline, over a million random hashes and 1000 buckets:

const std = @import("std");
const print = std.debug.print;

const num_buckets = 1000;
const n = 1_000_000;

const p1: u64 = num_buckets * 3 / 10;
const p1_top: u64 = @intFromFloat(0.6 * @as(f64, @floatFromInt(~@as(u64, 0))));

fn skewBucket(hash: u64) u64 {
    return if (hash < p1_top) hash % p1 else p1 + hash % (num_buckets - p1);
}

fn optBucket(hash: u64, eps: f64) u64 {
    const x = @as(f64, @floatFromInt(hash)) / @as(f64, @floatFromInt(~@as(u64, 0)));
    const beta = x + (1.0 - x) * @log(1.0 - x);
    const beta_eps = eps * x + (1.0 - eps) * beta;
    return @intFromFloat(beta_eps * num_buckets);
}

pub fn main() !void {
    var prng = std.Random.DefaultPrng.init(0x5eed);
    const random = prng.random();
    const nf = @as(f64, @floatFromInt(num_buckets));
    const lambda = nf / @as(f64, @floatFromInt(num_buckets));
    const eps = @min(1.0, lambda / (5.0 * std.math.sqrt(nf)));

    var uniform = [_]u64{0} ** num_buckets;
    var skew = [_]u64{0} ** num_buckets;
    var opt = [_]u64{0} ** num_buckets;
    for (0..n) |_| {
        const h = random.int(u64);
        uniform[h % num_buckets] += 1;
        skew[skewBucket(h)] += 1;
        opt[optBucket(h, eps)] += 1;
    }

    print("buckets     avg size: uniform     skew  optimal\n", .{});
    var i: usize = 0;
    while (i < num_buckets) : (i += 100) {
        var sums = [_]u64{ 0, 0, 0 };
        for (i..i + 100) |j| {
            sums[0] += uniform[j];
            sums[1] += skew[j];
            sums[2] += opt[j];
        }
        print("{: >4}..{: <4}   {d: >16.1}  {d: >7.1}  {d: >7.1}\n", .{
            i, i + 99,
            @as(f64, @floatFromInt(sums[0])) / 100.0, @as(f64, @floatFromInt(sums[1])) / 100.0,
            @as(f64, @floatFromInt(sums[2])) / 100.0,
        });
    }
}
Copy

This prints:

buckets     avg size: uniform     skew  optimal
   0..99                998.7   1996.4   4085.6
 100..199              1000.1   1992.3   1494.3
 200..299              1000.5   2005.0   1058.9
 300..399               996.3    570.8    815.0
 400..499              1003.4    569.3    660.3
 500..599              1001.9    573.7    547.1
 600..699              1001.2    572.4    458.2
 700..799               997.7    574.6    375.5
 800..899               998.4    571.0    299.9
 900..999              1001.8    574.5    205.1
Copy

Full Implementation

You can take a look at my full implementation here miara - pthash.zig. It has few differences from the example code presented in the article, and contains several additional optimizations.

Performance [OUTDATED]

I benchmarked against the C++ implementation by Giulio Ermanno Pibiri (jermp/pthash, commit cc4c9c9), same M1 machine, same parameters, EF-encoded pilots in both. BPK is bits per key, BAF is the bucket assignment function:

Query speed is about the same, the zig version wins some rows and loses others. Space is within 4% on the skew rows. My build is 1.2x to 2.6x slower depending on the configuration, so about 2x on average.

I only compared single-threaded non-partitioning versions. Partitioning support is on the TODO list.

Further Reading

• Giulio Ermanno Pibiri, Roberto Trani, PTHash: Revisiting FCH Minimal Perfect Hashing - the paper this article follows
• Stefan Hermann, Hans-Peter Lehmann, Giulio Ermanno Pibiri, Peter Sanders, Stefan Walzer, PHOBIC: Perfect Hashing with Optimized Bucket Sizes and Interleaved Coding - the optimal bucket assignment, and the interleaved pilot coding I haven't ported
• jermp/pthash - the reference C++ implementation
• miara - my Zig implementation