diff --git a/BENCHMARKS.md b/BENCHMARKS.md index f6e6375..645c6d6 100644 --- a/BENCHMARKS.md +++ b/BENCHMARKS.md @@ -482,6 +482,85 @@ The rows and columns of the uncompressed layouts are within 20% (chunked column 0.45 -> 0.49 ms, contiguous column 2.55 -> 2.61 ms). This run does not explain the slower windows. +## Concurrent reads + +**Not yet measured.** The harness exists; no numbers are published until it +has been run on an idle machine. The smoke runs used while building it (tiny +files, other jobs compiling on the box) are not results. + +The question: libhdf5's threadsafe build serialises every API call under one +global mutex, and h5py holds a global lock around every call too, so threads +reading through h5py cannot decode in parallel; h5py users scale with +processes. A clawhdf5 `File` is `Send + Sync`, and nothing on the read paths +this harness uses (`read_f32`, `read_f32_selection`) takes a library-wide +lock: the one mutex is the `File`'s chunk cache (keyed per dataset), taken by +full reads of chunked datasets for each chunk's O(1) lookup and insert, never +across a decode; hyperslab reads do not use the cache. How does +decoded throughput scale with threads on one open file, against h5py threads +and h5py processes on the same files? + +Workload (`crates/clawhdf5-bench/src/bin/concurrent_read.rs`; the h5py script +mirrors it): `/deflate.h5` and `/contiguous.h5`, each with 64 `f32` +datasets of 64 MiB decoded (`[16384, 1024]`; the deflate file chunked +`256 x 256`, level 4), written by clawhdf5 on first use and reused while +`manifest.json` matches. The data is a slowly varying ramp plus 8 bits of +noise per element, every value exact in `f32`, so both harnesses check what +they read; it deflates about 3.1x (128 MiB -> 40.7 MiB for two 64 MiB +datasets). For each layout and thread count +(1, 2, 4, 8, 16; fixed total work per repetition, split among the threads): + +- `distinct`: every dataset read in full once, thread `t` taking datasets + `t, t + T, ...`; +- `same`: 1024 random `256 x 256` hyperslabs of `d00` in total, from a seeded + splitmix64 stream that both harnesses generate identically. + +Reported per row: MB/s of decoded (selected) data from the median of the +repetitions, and scaling efficiency `MB/s(T) / (T x MB/s(1))`. Each worker +times itself from a start barrier; a repetition spans the earliest start to +the latest finish. Page cache: warm by default (each file is read once before +timing); `--cold` evicts the files with `posix_fadvise(POSIX_FADV_DONTNEED)` +before every repetition (no root needed; best effort). clawhdf5 opens one +`File` per repetition, shared by all threads; h5py threads share one +`h5py.File`; h5py processes (spawned before timing) each open the file inside +the timed region. + +Decode inside a single clawhdf5 read is itself parallel in this binary +(clawhdf5-format's `parallel` feature, enabled here through clawhdf5-agent; +it is off in the facade's default features), so a 1-thread clawhdf5 full read +of the deflate file already uses the whole rayon pool. Run both +`--decode-threads 1` (each read decodes on its calling thread, like h5py — +this isolates the API's own scaling) and the default pool. + +> **Run** (from the repository root). The default files take about 5.4 GiB +> of disk (4 GiB contiguous + about 1.3 GiB deflate). Generating them is +> memory-hungry because `FileBuilder` holds a whole file in memory: peak RSS +> was 676 MB for `--datasets 2 --mib 64` (2026-09-25, tank, +> `/usr/bin/time -f %M`), about 5x one file's decoded size, so expect about +> 21 GB at the defaults (once; later runs reuse the files). Put `--dir` on a +> real disk, not tmpfs, if `--cold` is to mean anything. +> +> ```bash +> DIR=/path/on/disk/concurrent-read +> BENCH=crates/clawhdf5-bench/scripts +> PY=.venv/bin/python # h5py 3.16 / HDF5 2.0 in this repo +> cargo build --release -p clawhdf5-bench --bin concurrent_read +> B=target/release/concurrent_read +> $B --dir $DIR --json claw-pool.json # generates on first run +> $B --dir $DIR --decode-threads 1 --json claw-1.json +> $PY $BENCH/concurrent_read_h5py.py --dir $DIR --executor threads --json h5py-threads.json +> $PY $BENCH/concurrent_read_h5py.py --dir $DIR --executor processes --json h5py-procs.json +> $PY $BENCH/compare_concurrent_read.py claw-1.json h5py-threads.json h5py-procs.json +> $PY $BENCH/compare_concurrent_read.py claw-pool.json h5py-threads.json h5py-procs.json +> ``` +> +> Cold page cache: add `--cold` to every harness command. Smoke test (seconds): +> `$B --dir /tmp/cr --datasets 4 --mib 1 --threads 1,2,4 --slabs 16 --reps 1` +> and the same `--threads/--slabs/--reps` to the h5py script. + +Other flags (both harnesses): `--threads`, `--reps`, `--slab`, `--slabs`, +`--seed`, `--modes distinct,same`, `--layouts deflate,contiguous`; sizes +(`--datasets`, `--mib`) only on the Rust harness, which writes the files. + ## Search harness baseline (v2.3.0) Produced by `cargo run --release -p clawhdf5-bench --bin search_harness -- --full` diff --git a/CHANGELOG.md b/CHANGELOG.md index 1df18a3..f1934bc 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -305,6 +305,13 @@ README claimed but nothing measured. - `footprint_bench` reports whether it built `float16` or `f32` stores and takes `--f32`; it had kept printing "f32" after the default changed. +- New `concurrent_read` harness, with an h5py counterpart + (`crates/clawhdf5-bench/scripts/concurrent_read_h5py.py`, threads or + processes) and `compare_concurrent_read.py`: decoded read throughput and + scaling efficiency at 1-16 threads on one open file, full reads of distinct + datasets and random hyperslabs of one dataset, deflate and contiguous, warm + or `--cold` page cache, JSON output. Not yet measured — `BENCHMARKS.md` + ("Concurrent reads") has the commands and no numbers. ### Interop - **h5py could not open chunked datasets we wrote with a chunk dimension diff --git a/crates/clawhdf5-bench/Cargo.toml b/crates/clawhdf5-bench/Cargo.toml index 5ac6350..b48afbf 100644 --- a/crates/clawhdf5-bench/Cargo.toml +++ b/crates/clawhdf5-bench/Cargo.toml @@ -34,6 +34,10 @@ path = "src/bin/consolidation_efficiency.rs" name = "ephemeral_perf" path = "src/bin/ephemeral_perf.rs" +[[bin]] +name = "concurrent_read" +path = "src/bin/concurrent_read.rs" + [[bin]] name = "mpi_io_bench" path = "src/bin/mpi_io_bench.rs" @@ -64,6 +68,10 @@ clawhdf5-io = { path = "../clawhdf5-io" } mpi = { version = "0.8", optional = true } serde = { workspace = true } serde_json = "1" +# concurrent_read: size the decode pool (--decode-threads) and evict files +# from the page cache (--cold, posix_fadvise). Both pure Rust / bindings only. +rayon = "1" +libc = "0.2" tempfile = { workspace = true } # Optional: libhdf5 C wrapper for side-by-side comparison (requires system libhdf5). # Enable with: cargo bench -p clawhdf5-bench --features libhdf5-compare diff --git a/crates/clawhdf5-bench/scripts/__pycache__/compare_concurrent_read.cpython-314.pyc b/crates/clawhdf5-bench/scripts/__pycache__/compare_concurrent_read.cpython-314.pyc new file mode 100644 index 0000000..85432af Binary files /dev/null and b/crates/clawhdf5-bench/scripts/__pycache__/compare_concurrent_read.cpython-314.pyc differ diff --git a/crates/clawhdf5-bench/scripts/__pycache__/concurrent_read_h5py.cpython-314.pyc b/crates/clawhdf5-bench/scripts/__pycache__/concurrent_read_h5py.cpython-314.pyc new file mode 100644 index 0000000..4f759f2 Binary files /dev/null and b/crates/clawhdf5-bench/scripts/__pycache__/concurrent_read_h5py.cpython-314.pyc differ diff --git a/crates/clawhdf5-bench/scripts/compare_concurrent_read.py b/crates/clawhdf5-bench/scripts/compare_concurrent_read.py new file mode 100644 index 0000000..b7c75c2 --- /dev/null +++ b/crates/clawhdf5-bench/scripts/compare_concurrent_read.py @@ -0,0 +1,70 @@ +#!/usr/bin/env python3 +"""Tabulate concurrent_read JSON results (clawhdf5, h5py threads/processes). + + python compare_concurrent_read.py clawhdf5.json h5py-threads.json h5py-procs.json + +Prints one Markdown table: for each layout, mode and thread count, every +tool's MB/s and scaling efficiency, and the first file's MB/s relative to each +of the others. Refuses to compare runs whose workload parameters differ. +""" + +import json +import sys + +COMPARED = ("datasets", "rows", "cols", "chunk", "deflate_level", "slab", "slabs", "seed") + + +def main(paths): + if len(paths) < 2: + sys.exit(__doc__) + docs = [] + for p in paths: + with open(p) as fh: + docs.append(json.load(fh)) + ref = docs[0] + for d, p in zip(docs[1:], paths[1:]): + diff = [k for k in COMPARED if d["params"].get(k) != ref["params"].get(k)] + if diff: + sys.exit(f"{p}: workload differs from {paths[0]} in {', '.join(diff)}") + if d["cache"] != ref["cache"]: + print(f"warning: {p} ran {d['cache']!r}, {paths[0]} ran {ref['cache']!r}", + file=sys.stderr) + if d.get("host") != ref.get("host"): + print(f"warning: {p} ran on {d.get('host')}, {paths[0]} on {ref.get('host')}", + file=sys.stderr) + + names = [d["tool"] for d in docs] + for d in docs: + extra = f", HDF5 {d['hdf5_version']}" if "hdf5_version" in d else "" + print(f"- {d['tool']} {d['version']}{extra}: host {d.get('host')}, " + f"{d.get('cpus')} CPUs, cache {d['cache']}, decode threads per read " + f"{d.get('decode_threads')}") + p = ref["params"] + print(f"\n{p['datasets']} datasets of {p['rows']} x {p['cols']} f32, chunks " + f"{p['chunk'][0]} x {p['chunk'][1]} (deflate {p['deflate_level']}); " + f"`same`: {p['slabs']} slabs of {p['slab']} x {p['slab']}\n") + + index = [{(r["layout"], r["mode"], r["threads"]): r for r in d["results"]} for d in docs] + keys = [(r["layout"], r["mode"], r["threads"]) for r in ref["results"]] + + head = ["layout", "mode", "threads"] + head += [f"{n} MB/s (eff)" for n in names] + head += [f"{names[0]} / {n}" for n in names[1:]] + print("| " + " | ".join(head) + " |") + print("|---|---|" + "---:|" * (len(head) - 2)) + for key in keys: + cells = [key[0], key[1], str(key[2])] + rs = [ix.get(key) for ix in index] + for r in rs: + if r is None: + cells.append("-") + else: + eff = "-" if r["efficiency"] is None else f"{r['efficiency']:.2f}" + cells.append(f"{r['mb_s']:.0f} ({eff})") + for r in rs[1:]: + cells.append("-" if r is None else f"{rs[0]['mb_s'] / r['mb_s']:.2f}x") + print("| " + " | ".join(cells) + " |") + + +if __name__ == "__main__": + main(sys.argv[1:]) diff --git a/crates/clawhdf5-bench/scripts/concurrent_read_h5py.py b/crates/clawhdf5-bench/scripts/concurrent_read_h5py.py new file mode 100644 index 0000000..883a850 --- /dev/null +++ b/crates/clawhdf5-bench/scripts/concurrent_read_h5py.py @@ -0,0 +1,265 @@ +#!/usr/bin/env python3 +"""The concurrent_read workload with h5py, on the files concurrent_read wrote. + +libhdf5 serialises every API call under one global lock, and h5py holds its +own global lock around every call as well, so h5py *threads* cannot decode in +parallel. h5py users scale with *processes* instead; ``--executor processes`` +measures that (each worker opens the file itself). + +The workload mirrors ``crates/clawhdf5-bench/src/bin/concurrent_read.rs``: + +* ``distinct``: every dataset read in full once per repetition; worker ``t`` + of ``T`` reads datasets ``t, t + T, ...``. +* ``same``: ``--slabs`` random ``--slab`` x ``--slab`` hyperslabs of ``d00`` + (slab ``j`` to worker ``j % T``), offsets from the same splitmix64 stream. + +Each worker times itself from a start barrier; a repetition spans the earliest +start to the latest finish (CLOCK_MONOTONIC, comparable across processes). +Threads share one ``h5py.File`` per repetition; process workers open the file +inside the timed region (a few ms against reads of many MiB). + +Generate the files first with the Rust harness (it writes ``manifest.json``), +then, for example:: + + python concurrent_read_h5py.py --dir DIR --executor threads --json h5py-threads.json + python concurrent_read_h5py.py --dir DIR --executor processes --json h5py-procs.json +""" + +import argparse +import json +import multiprocessing as mp +import os +import platform +import socket +import sys +import threading +import time + +import h5py +import numpy as np + +M64 = (1 << 64) - 1 + + +def splitmix64(state): + """Return (new_state, value); the same stream as the Rust harness.""" + state = (state + 0x9E3779B97F4A7C15) & M64 + z = state + z = ((z ^ (z >> 30)) * 0xBF58476D1CE4E5B9) & M64 + z = ((z ^ (z >> 27)) * 0x94D049BB133111EB) & M64 + return state, z ^ (z >> 31) + + +def value(k, i): + """Element i (row-major) of dataset k, exactly as concurrent_read writes it.""" + _, noise = splitmix64(i ^ (k << 40)) + return np.float32((((i >> 6) % 16384) + k) + (noise & 0xFF) / 256.0) + + +def slab_offsets(seed, count, rows, cols, slab): + s = seed + out = [] + for _ in range(count): + s, r = splitmix64(s) + s, c = splitmix64(s) + out.append((r % (rows - slab + 1), c % (cols - slab + 1))) + return out + + +def now(): + return time.clock_gettime(time.CLOCK_MONOTONIC) + + +def work(f, mode, t, threads, m, slabs, slab, verify): + """Worker t's share of one repetition on an open h5py.File.""" + n = m["rows"] * m["cols"] + if mode == "distinct": + for k in range(t, m["datasets"], threads): + got = f[f"d{k:02d}"][...] + assert got.size == n + if verify: + flat = got.reshape(-1) + for i in (0, n // 3, n - 1): + assert flat[i] == value(k, i), f"d{k:02d}[{i}]" + else: + ds = f["d00"] + cols = m["cols"] + for r, c in slabs[t::threads]: + got = ds[r : r + slab, c : c + slab] + assert got.shape == (slab, slab) + if verify: + assert got[0, 0] == value(0, r * cols + c) + last = (r + slab - 1) * cols + c + slab - 1 + assert got[-1, -1] == value(0, last) + + +# ----- process workers ------------------------------------------------------ + +_barrier = None + + +def _init(barrier): + global _barrier + _barrier = barrier + + +def _proc_task(task): + path, mode, t, threads, m, slabs, slab = task + _barrier.wait() + start = now() + with h5py.File(path, "r") as f: + work(f, mode, t, threads, m, slabs, slab, False) + return start, now() + + +def _noop(_): + return os.getpid() + + +def run_threads(path, mode, threads, m, slabs, slab): + spans = [None] * threads + barrier = threading.Barrier(threads) + with h5py.File(path, "r") as f: + + def body(t): + barrier.wait() + start = now() + work(f, mode, t, threads, m, slabs, slab, False) + spans[t] = (start, now()) + + ts = [threading.Thread(target=body, args=(t,)) for t in range(threads)] + for th in ts: + th.start() + for th in ts: + th.join() + return max(e for _, e in spans) - min(s for s, _ in spans) + + +def run_processes(pool, path, mode, threads, m, slabs, slab): + tasks = [(path, mode, t, threads, m, slabs, slab) for t in range(threads)] + # One task per worker: each blocks in the barrier until all T have + # started, so no worker can take a second task. + spans = pool.map(_proc_task, tasks, chunksize=1) + return max(e for _, e in spans) - min(s for s, _ in spans) + + +def warm(path): + with open(path, "rb") as fh: + while fh.read(1 << 24): + pass + + +def evict(path): + fd = os.open(path, os.O_RDONLY) + try: + os.posix_fadvise(fd, 0, 0, os.POSIX_FADV_DONTNEED) + finally: + os.close(fd) + + +def main(): + ap = argparse.ArgumentParser(description=__doc__.split("\n\n")[0]) + ap.add_argument("--dir", default="concurrent-read-data") + ap.add_argument("--executor", choices=["threads", "processes"], default="threads") + ap.add_argument("--threads", default="1,2,4,8,16") + ap.add_argument("--reps", type=int, default=3) + ap.add_argument("--slab", type=int, default=256) + ap.add_argument("--slabs", type=int, default=1024) + ap.add_argument("--seed", type=int, default=42) + ap.add_argument("--cold", action="store_true") + ap.add_argument("--modes", default="distinct,same") + ap.add_argument("--layouts", default="deflate,contiguous") + ap.add_argument("--json") + a = ap.parse_args() + + # The Rust harness pins this value (splitmix64_reference). + assert splitmix64(42)[1] == 0xBDD732262FEB6E95, "splitmix64 port is wrong" + + try: + with open(os.path.join(a.dir, "manifest.json")) as fh: + m = json.load(fh) + except FileNotFoundError: + sys.exit(f"{a.dir}/manifest.json not found: generate the files with " + "`cargo run --release -p clawhdf5-bench --bin concurrent_read -- --dir ...` first") + threads_list = [int(x) for x in a.threads.split(",")] + modes = a.modes.split(",") + layouts = a.layouts.split(",") + if a.slab < 1 or a.slab > min(m["rows"], m["cols"]): + sys.exit(f"--slab must be 1..={min(m['rows'], m['cols'])}") + files = dict(m["files"]) + slabs = slab_offsets(a.seed, a.slabs, m["rows"], m["cols"], a.slab) + dataset_bytes = m["rows"] * m["cols"] * 4 + tool = f"h5py-{a.executor}" + + ctx = mp.get_context("spawn") # never fork a process holding HDF5 state + pools = {} + if a.executor == "processes": + for t in threads_list: + pool = ctx.Pool(t, initializer=_init, initargs=(ctx.Barrier(t),)) + pool.map(_noop, range(t)) # start the workers outside the timing + pools[t] = pool + + rows = [] + print("| layout | mode | threads | MB/s | efficiency | median s |") + print("|---|---|---:|---:|---:|---:|") + try: + for layout in layouts: + path = os.path.join(a.dir, files[layout]) + if not a.cold: + warm(path) + for mode in modes: + with h5py.File(path, "r") as f: # untimed, checked pass + work(f, mode, 0, 1, m, slabs, a.slab, True) + nbytes = (dataset_bytes * m["datasets"] if mode == "distinct" + else a.slab * a.slab * 4 * a.slabs) + base = None + for t in threads_list: + times = [] + for _ in range(a.reps): + if a.cold: + evict(path) + if a.executor == "threads": + times.append(run_threads(path, mode, t, m, slabs, a.slab)) + else: + times.append(run_processes(pools[t], path, mode, t, m, slabs, a.slab)) + med = sorted(times)[len(times) // 2] + mb_s = nbytes / (1 << 20) / med + if t == 1: + base = mb_s + eff = mb_s / (t * base) if base else None + print(f"| {layout} | {mode} | {t} | {mb_s:.0f} | " + f"{'-' if eff is None else f'{eff:.2f}'} | {med:.4f} |") + rows.append({ + "layout": layout, "mode": mode, "threads": t, "bytes": nbytes, + "times_s": times, "median_s": med, "mb_s": mb_s, "efficiency": eff, + }) + finally: + for pool in pools.values(): + pool.terminate() + + if a.json: + doc = { + "tool": tool, + "version": h5py.__version__, + "hdf5_version": h5py.version.hdf5_version, + "python": platform.python_version(), + "host": socket.gethostname(), + "cpus": os.cpu_count(), + "unix_time": int(time.time()), + "cache": ("cold (posix_fadvise DONTNEED before each repetition)" + if a.cold else "warm"), + "decode_threads": 1, + "params": { + "datasets": m["datasets"], "rows": m["rows"], "cols": m["cols"], + "chunk": m["chunk"], "deflate_level": m["deflate_level"], + "mib": dataset_bytes // (1 << 20), "slab": a.slab, "slabs": a.slabs, + "seed": a.seed, "reps": a.reps, "dir": a.dir, + }, + "results": rows, + } + with open(a.json, "w") as fh: + json.dump(doc, fh, indent=2) + + +if __name__ == "__main__": + main() diff --git a/crates/clawhdf5-bench/src/bin/concurrent_read.rs b/crates/clawhdf5-bench/src/bin/concurrent_read.rs new file mode 100644 index 0000000..a87ac26 --- /dev/null +++ b/crates/clawhdf5-bench/src/bin/concurrent_read.rs @@ -0,0 +1,523 @@ +//! Concurrent-read harness: how does decoded read throughput scale with the +//! number of threads reading one open file? +//! +//! libhdf5 (threadsafe build) serialises every API call under one global +//! mutex, and h5py holds it too, so threads cannot decode in parallel there. +//! A clawhdf5 [`File`] is `Send + Sync`; this harness measures what that buys. +//! `crates/clawhdf5-bench/scripts/concurrent_read_h5py.py` runs the same +//! workload on the same files with h5py (threads, and processes), and +//! `compare_concurrent_read.py` tabulates the JSON both write. +//! +//! Files (generated on first use, reused while `manifest.json` matches): +//! +//! * `/deflate.h5`: `--datasets` datasets `d00`, `d01`, ... of `f32`, +//! `--mib` MiB decoded each, shape `[mib * 256, 1024]`, chunks `256 x 256`, +//! deflate level 4. +//! * `/contiguous.h5`: the same datasets, contiguous. +//! +//! Modes, for each layout and each thread count `T` (strong scaling: the total +//! work per repetition is fixed, split among the threads): +//! +//! * `distinct`: every dataset is read in full once; thread `t` reads datasets +//! `t, t + T, t + 2T, ...`. +//! * `same`: all threads read `d00`, `--slabs` random `--slab` x `--slab` +//! hyperslabs in total (slab `j` goes to thread `j % T`). The offsets come +//! from a splitmix64 stream seeded with `--seed`, identical in the h5py +//! script. +//! +//! One `File` per layout per repetition is shared by all threads (opened +//! fresh each repetition, so no chunk cache carries over). Page cache: +//! `warm` (default) reads every file once before timing; `--cold` evicts the +//! files from the page cache with `posix_fadvise(POSIX_FADV_DONTNEED)` before +//! every repetition (no root needed; it only evicts clean, unmapped pages, so +//! it is best effort — the JSON says which was used). +//! +//! Decode inside one read is itself parallel when clawhdf5-format's `parallel` +//! feature is on (it is in this binary, via clawhdf5-agent). `--decode-threads +//! N` sizes that rayon pool; `--decode-threads 1` measures the API's own +//! thread scaling, comparable with h5py where each call decodes on the +//! calling thread. +//! +//! ```text +//! cargo run --release -p clawhdf5-bench --bin concurrent_read -- \ +//! --dir /data/concurrent-read --json clawhdf5.json +//! cargo run --release -p clawhdf5-bench --bin concurrent_read -- \ +//! --dir /tmp/cr --datasets 4 --mib 1 --threads 1,2 --slabs 16 --reps 1 # smoke +//! ``` + +use std::path::{Path, PathBuf}; +use std::sync::Barrier; +use std::time::Instant; + +use clawhdf5::{File, FileBuilder, Selection}; +use serde::{Deserialize, Serialize}; + +const COLS: u64 = 1024; +const ROWS_PER_MIB: u64 = 256; // 256 rows x 1024 cols x 4 bytes = 1 MiB +const CHUNK: u64 = 256; +const DEFLATE_LEVEL: u32 = 4; +const LAYOUTS: [&str; 2] = ["deflate", "contiguous"]; +const MANIFEST_VERSION: u32 = 1; + +/// splitmix64 — shared with the h5py script, which must produce the same +/// stream (both the data and the hyperslab offsets depend on it). +fn splitmix64(state: &mut u64) -> u64 { + *state = state.wrapping_add(0x9E37_79B9_7F4A_7C15); + let mut z = *state; + z = (z ^ (z >> 30)).wrapping_mul(0xBF58_476D_1CE4_E5B9); + z = (z ^ (z >> 27)).wrapping_mul(0x94D0_49BB_1331_11EB); + z ^ (z >> 31) +} + +/// Element `i` (row-major) of dataset `k`: a slowly varying integer part plus +/// 8 bits of noise, so deflate has real work to do (about 3.1x) and every value +/// is exact in `f32` (< 2^15 with 8 fraction bits), which lets both harnesses +/// check what they read against this formula. +fn value(k: u64, i: u64) -> f32 { + let mut s = i ^ (k << 40); + let noise = splitmix64(&mut s) & 0xff; + (((i >> 6) % 16384) + k) as f32 + noise as f32 / 256.0 +} + +#[derive(Serialize, Deserialize, PartialEq, Debug, Clone)] +struct Manifest { + version: u32, + datasets: u64, + rows: u64, + cols: u64, + chunk: [u64; 2], + deflate_level: u32, + files: Vec<(String, String)>, // (layout, file name) + writer: String, +} + +fn manifest_for(datasets: u64, mib: u64) -> Manifest { + Manifest { + version: MANIFEST_VERSION, + datasets, + rows: mib * ROWS_PER_MIB, + cols: COLS, + chunk: [CHUNK, CHUNK], + deflate_level: DEFLATE_LEVEL, + files: LAYOUTS + .iter() + .map(|l| (l.to_string(), format!("{l}.h5"))) + .collect(), + writer: format!("clawhdf5 {}", env!("CARGO_PKG_VERSION")), + } +} + +fn dataset_values(k: u64, n: u64) -> Vec { + (0..n).map(|i| value(k, i)).collect() +} + +/// Write the files unless `dir` already holds ones matching `want`. +fn ensure_files(dir: &Path, want: &Manifest) -> std::io::Result { + let manifest_path = dir.join("manifest.json"); + if let Ok(text) = std::fs::read_to_string(&manifest_path) + && let Ok(have) = serde_json::from_str::(&text) + && have.version == want.version + && have.datasets == want.datasets + && have.rows == want.rows + && have.cols == want.cols + && have.chunk == want.chunk + && have.deflate_level == want.deflate_level + && have.files == want.files + && want.files.iter().all(|(_, f)| dir.join(f).exists()) + { + return Ok(false); + } + std::fs::create_dir_all(dir)?; + // A stale manifest must not survive a half-written regeneration. + let _ = std::fs::remove_file(&manifest_path); + let n = want.rows * want.cols; + for (layout, file) in &want.files { + // One layout at a time keeps the peak memory to about twice one + // file's decoded size. + let mut b = FileBuilder::new(); + for k in 0..want.datasets { + let ds = b.create_dataset(&format!("d{k:02}")); + ds.with_f32_data(&dataset_values(k, n)) + .with_shape(&[want.rows, want.cols]); + if layout == "deflate" { + ds.with_chunks(&[CHUNK.min(want.rows), CHUNK]) + .with_deflate(DEFLATE_LEVEL); + } + } + b.write(dir.join(file)).map_err(std::io::Error::other)?; + } + std::fs::write( + &manifest_path, + serde_json::to_string_pretty(want).map_err(std::io::Error::other)?, + )?; + Ok(true) +} + +fn slab_offsets(seed: u64, count: usize, rows: u64, cols: u64, slab: u64) -> Vec<(u64, u64)> { + let mut s = seed; + (0..count) + .map(|_| { + let r = splitmix64(&mut s) % (rows - slab + 1); + let c = splitmix64(&mut s) % (cols - slab + 1); + (r, c) + }) + .collect() +} + +/// Warm the page cache by reading every byte of `path`. +fn warm(path: &Path) -> std::io::Result<()> { + let mut f = std::fs::File::open(path)?; + std::io::copy(&mut f, &mut std::io::sink())?; + Ok(()) +} + +/// Ask the kernel to drop `path`'s pages from the page cache. +fn evict(path: &Path) -> std::io::Result<()> { + use std::os::fd::AsRawFd; + let f = std::fs::File::open(path)?; + // SAFETY: plain syscall on a valid, open file descriptor. + let rc = unsafe { libc::posix_fadvise(f.as_raw_fd(), 0, 0, libc::POSIX_FADV_DONTNEED) }; + if rc != 0 { + return Err(std::io::Error::from_raw_os_error(rc)); + } + Ok(()) +} + +#[derive(Serialize)] +struct Row { + layout: String, + mode: String, + threads: usize, + /// Decoded (selected) bytes read per repetition. + bytes: u64, + times_s: Vec, + median_s: f64, + mb_s: f64, + /// `mb_s / (threads * mb_s at threads = 1)`; null without a 1-thread row. + efficiency: Option, +} + +struct Args { + dir: PathBuf, + datasets: u64, + mib: u64, + threads: Vec, + reps: usize, + slab: u64, + slabs: usize, + seed: u64, + cold: bool, + decode_threads: usize, + modes: Vec, + layouts: Vec, + json: Option, +} + +const USAGE: &str = "\ +usage: concurrent_read [--dir DIR] [--datasets N] [--mib N] [--threads 1,2,4,8,16] + [--reps N] [--slab N] [--slabs N] [--seed N] [--cold] + [--decode-threads N] [--modes distinct,same] + [--layouts deflate,contiguous] [--json FILE]"; + +fn parse_list(s: &str) -> Result, String> { + s.split(',') + .map(|x| x.trim().parse().map_err(|_| format!("bad list item {x:?}"))) + .collect() +} + +fn parse_args() -> Result { + let mut a = Args { + dir: PathBuf::from("concurrent-read-data"), + datasets: 64, + mib: 64, + threads: vec![1, 2, 4, 8, 16], + reps: 3, + slab: 256, + slabs: 1024, + seed: 42, + cold: false, + decode_threads: 0, + modes: vec!["distinct".into(), "same".into()], + layouts: LAYOUTS.iter().map(|s| s.to_string()).collect(), + json: None, + }; + let mut it = std::env::args().skip(1); + while let Some(flag) = it.next() { + if flag == "--cold" { + a.cold = true; + continue; + } + if flag == "-h" || flag == "--help" { + return Err(USAGE.into()); + } + let v = it.next().ok_or(format!("{flag} needs a value\n{USAGE}"))?; + let num = |v: &str| { + v.parse::() + .map_err(|_| format!("{flag}: bad number {v:?}")) + }; + match flag.as_str() { + "--dir" => a.dir = v.into(), + "--datasets" => a.datasets = num(&v)?, + "--mib" => a.mib = num(&v)?, + "--threads" => a.threads = parse_list(&v)?, + "--reps" => a.reps = num(&v)? as usize, + "--slab" => a.slab = num(&v)?, + "--slabs" => a.slabs = num(&v)? as usize, + "--seed" => a.seed = num(&v)?, + "--decode-threads" => a.decode_threads = num(&v)? as usize, + "--modes" => a.modes = parse_list(&v)?, + "--layouts" => a.layouts = parse_list(&v)?, + "--json" => a.json = Some(v.into()), + _ => return Err(format!("unknown flag {flag}\n{USAGE}")), + } + } + if a.datasets == 0 || a.datasets > 100 { + return Err("--datasets must be 1..=100".into()); + } + if a.mib == 0 || a.reps == 0 || a.slabs == 0 || a.threads.contains(&0) { + return Err("--mib, --reps, --slabs and every --threads value must be > 0".into()); + } + if a.slab == 0 || a.slab > COLS || a.slab > a.mib * ROWS_PER_MIB { + return Err(format!( + "--slab must be 1..={}", + COLS.min(a.mib * ROWS_PER_MIB) + )); + } + for m in &a.modes { + if m != "distinct" && m != "same" { + return Err(format!("unknown mode {m:?}")); + } + } + for l in &a.layouts { + if !LAYOUTS.contains(&l.as_str()) { + return Err(format!("unknown layout {l:?}")); + } + } + Ok(a) +} + +/// One timed repetition: `T` threads on one shared `File`. Returns seconds. +fn run_once( + path: &Path, + mode: &str, + threads: usize, + m: &Manifest, + slabs: &[(u64, u64)], + slab: u64, + verify: bool, +) -> f64 { + let file = File::open(path).expect("open"); + let barrier = Barrier::new(threads + 1); // + the spawning thread + let n = m.rows * m.cols; + // Each thread times itself from the barrier; the repetition spans the + // earliest start to the latest finish (timing on the spawning thread + // instead undercounts whenever it is scheduled after the workers ran). + let spans: Vec<(Instant, Instant)> = std::thread::scope(|s| { + let handles: Vec<_> = (0..threads) + .map(|t| { + let (file, barrier) = (&file, &barrier); + s.spawn(move || { + barrier.wait(); + let start = Instant::now(); + match mode { + "distinct" => { + for k in (t as u64..m.datasets).step_by(threads) { + let got = file.dataset(&format!("d{k:02}")).unwrap().read_f32(); + let got = got.unwrap(); + assert_eq!(got.len() as u64, n); + if verify { + for i in [0, n / 3, n - 1] { + assert_eq!(got[i as usize], value(k, i), "d{k:02}[{i}]"); + } + } + std::hint::black_box(got); + } + } + _ => { + let ds = file.dataset("d00").unwrap(); + for &(r, c) in slabs.iter().skip(t).step_by(threads) { + let sel = Selection::Hyperslab { + start: vec![r, c], + stride: vec![1, 1], + count: vec![slab, slab], + block: vec![1, 1], + }; + let got = ds.read_f32_selection(&sel).unwrap(); + assert_eq!(got.len() as u64, slab * slab); + if verify { + let last = (r + slab - 1) * m.cols + c + slab - 1; + assert_eq!(got[0], value(0, r * m.cols + c)); + assert_eq!(*got.last().unwrap(), value(0, last)); + } + std::hint::black_box(got); + } + } + } + (start, Instant::now()) + }) + }) + .collect(); + barrier.wait(); + handles.into_iter().map(|h| h.join().unwrap()).collect() + }); + let start = spans.iter().map(|s| s.0).min().unwrap(); + let end = spans.iter().map(|s| s.1).max().unwrap(); + (end - start).as_secs_f64() +} + +fn median(v: &[f64]) -> f64 { + let mut s = v.to_vec(); + s.sort_by(f64::total_cmp); + s[s.len() / 2] +} + +fn hostname() -> String { + std::fs::read_to_string("/proc/sys/kernel/hostname") + .map(|s| s.trim().to_string()) + .unwrap_or_else(|_| "unknown".into()) +} + +fn main() { + let args = match parse_args() { + Ok(a) => a, + Err(e) => { + eprintln!("{e}"); + std::process::exit(2); + } + }; + if cfg!(debug_assertions) { + eprintln!("warning: debug build — numbers are meaningless. Use --release."); + } + if args.decode_threads > 0 { + rayon::ThreadPoolBuilder::new() + .num_threads(args.decode_threads) + .build_global() + .expect("configure rayon pool"); + } + + let manifest = manifest_for(args.datasets, args.mib); + let t = Instant::now(); + match ensure_files(&args.dir, &manifest) { + Ok(true) => eprintln!( + "generated {} in {:.1} s", + args.dir.display(), + t.elapsed().as_secs_f64() + ), + Ok(false) => eprintln!("reusing {}", args.dir.display()), + Err(e) => { + eprintln!("cannot write test files in {}: {e}", args.dir.display()); + std::process::exit(1); + } + } + let path_of = |layout: &str| args.dir.join(format!("{layout}.h5")); + let slabs = slab_offsets( + args.seed, + args.slabs, + manifest.rows, + manifest.cols, + args.slab, + ); + let dataset_bytes = manifest.rows * manifest.cols * 4; + + let mut rows: Vec = Vec::new(); + println!("| layout | mode | threads | MB/s | efficiency | median s |"); + println!("|---|---|---:|---:|---:|---:|"); + for layout in &args.layouts { + let path = path_of(layout); + // Untimed pass: page cache warm (unless --cold), results checked. + if !args.cold { + warm(&path).expect("warm page cache"); + } + for mode in &args.modes { + run_once(&path, mode, 1, &manifest, &slabs, args.slab, true); + let bytes = match mode.as_str() { + "distinct" => dataset_bytes * manifest.datasets, + _ => args.slab * args.slab * 4 * args.slabs as u64, + }; + let mut base: Option = None; + for &threads in &args.threads { + let times: Vec = (0..args.reps) + .map(|_| { + if args.cold { + evict(&path).expect("posix_fadvise"); + } + run_once(&path, mode, threads, &manifest, &slabs, args.slab, false) + }) + .collect(); + let med = median(×); + let mb_s = bytes as f64 / (1 << 20) as f64 / med; + if threads == 1 { + base = Some(mb_s); + } + let efficiency = base.map(|b| mb_s / (threads as f64 * b)); + println!( + "| {layout} | {mode} | {threads} | {mb_s:.0} | {} | {med:.4} |", + efficiency.map_or("-".into(), |e| format!("{e:.2}")) + ); + rows.push(Row { + layout: layout.clone(), + mode: mode.clone(), + threads, + bytes, + times_s: times, + median_s: med, + mb_s, + efficiency, + }); + } + } + } + + if let Some(out) = &args.json { + let doc = serde_json::json!({ + "tool": "clawhdf5", + "version": env!("CARGO_PKG_VERSION"), + "host": hostname(), + "cpus": std::thread::available_parallelism().map_or(0, |n| n.get()), + "unix_time": std::time::SystemTime::now() + .duration_since(std::time::UNIX_EPOCH) + .map_or(0, |d| d.as_secs()), + "cache": if args.cold { "cold (posix_fadvise DONTNEED before each repetition)" } else { "warm" }, + "decode_threads": rayon::current_num_threads(), + "params": { + "datasets": manifest.datasets, + "mib": args.mib, + "rows": manifest.rows, + "cols": manifest.cols, + "chunk": manifest.chunk, + "deflate_level": manifest.deflate_level, + "slab": args.slab, + "slabs": args.slabs, + "seed": args.seed, + "reps": args.reps, + "dir": args.dir, + }, + "results": rows, + }); + std::fs::write(out, serde_json::to_string_pretty(&doc).unwrap()).expect("write json"); + } +} + +#[cfg(test)] +mod tests { + use super::*; + + #[test] + fn values_are_exact_in_f32() { + for k in [0, 7, 63] { + for i in [0u64, 1, 4095, 1 << 20, (1 << 24) - 1] { + let v = value(k, i); + assert_eq!(v, (v as f64) as f32); + assert!(v < 32768.0); + assert_eq!((v * 256.0).fract(), 0.0); + } + } + } + + /// The h5py script hard-codes this vector to check its splitmix64 port. + #[test] + fn splitmix64_reference() { + let mut s = 42; + assert_eq!(splitmix64(&mut s), 0xBDD7_3226_2FEB_6E95); + } +} diff --git a/crates/clawhdf5-bench/tests/concurrent_read_smoke.rs b/crates/clawhdf5-bench/tests/concurrent_read_smoke.rs new file mode 100644 index 0000000..a8ab5c7 --- /dev/null +++ b/crates/clawhdf5-bench/tests/concurrent_read_smoke.rs @@ -0,0 +1,148 @@ +//! Keeps the concurrent-read harnesses working: runs `concurrent_read`, the +//! h5py script (threads and processes) and the comparison script end to end +//! on tiny files. h5py reading the files also checks, element by element at +//! spot positions, that both harnesses generate the same data and slabs. +//! +//! The h5py half is skipped when python3 with h5py is unavailable, unless +//! `CLAWHDF5_REQUIRE_INTEROP=1`; `CLAWHDF5_PYTHON` picks the interpreter. + +use std::path::{Path, PathBuf}; +use std::process::Command; + +fn python() -> String { + std::env::var("CLAWHDF5_PYTHON").unwrap_or_else(|_| "python3".to_string()) +} + +fn interop_required() -> bool { + std::env::var("CLAWHDF5_REQUIRE_INTEROP").is_ok_and(|v| v == "1") +} + +fn python_available() -> bool { + Command::new(python()) + .args(["-c", "import h5py, numpy"]) + .output() + .map(|o| o.status.success()) + .unwrap_or(false) +} + +fn scripts() -> PathBuf { + Path::new(env!("CARGO_MANIFEST_DIR")).join("scripts") +} + +fn run(cmd: &mut Command) -> String { + let out = cmd.output().expect("spawn"); + assert!( + out.status.success(), + "{cmd:?} failed\nSTDOUT:\n{}\nSTDERR:\n{}", + String::from_utf8_lossy(&out.stdout), + String::from_utf8_lossy(&out.stderr) + ); + String::from_utf8_lossy(&out.stdout).into_owned() +} + +const SMALL: [&str; 8] = [ + "--threads", + "1,2", + "--slabs", + "8", + "--reps", + "1", + "--slab", + "64", +]; + +fn results(path: &Path) -> serde_json::Value { + serde_json::from_str(&std::fs::read_to_string(path).unwrap()).unwrap() +} + +#[test] +fn harnesses_run_end_to_end_on_tiny_files() { + let dir = tempfile::TempDir::new().unwrap(); + let data = dir.path().join("data"); + let claw = dir.path().join("claw.json"); + + let bin = env!("CARGO_BIN_EXE_concurrent_read"); + run(Command::new(bin) + .arg("--dir") + .arg(&data) + .args(["--datasets", "3", "--mib", "1"]) + .args(SMALL) + .arg("--json") + .arg(&claw)); + // Second run reuses the files (and exercises --cold). + let out = Command::new(bin) + .arg("--dir") + .arg(&data) + .args(["--datasets", "3", "--mib", "1", "--cold"]) + .args(SMALL) + .output() + .unwrap(); + assert!(out.status.success()); + assert!(String::from_utf8_lossy(&out.stderr).contains("reusing")); + + let doc = results(&claw); + assert_eq!(doc["tool"], "clawhdf5"); + // 2 layouts x 2 modes x 2 thread counts. + assert_eq!(doc["results"].as_array().unwrap().len(), 8); + for r in doc["results"].as_array().unwrap() { + assert!(r["mb_s"].as_f64().unwrap() > 0.0, "{r}"); + } + + if !python_available() { + assert!( + !interop_required(), + "CLAWHDF5_REQUIRE_INTEROP=1 but {} has no h5py", + python() + ); + eprintln!("skipping the h5py half: no h5py in {}", python()); + return; + } + let mut jsons = vec![claw]; + for executor in ["threads", "processes"] { + let out = dir.path().join(format!("h5py-{executor}.json")); + run(Command::new(python()) + .arg(scripts().join("concurrent_read_h5py.py")) + .arg("--dir") + .arg(&data) + .args(["--executor", executor]) + .args(SMALL) + .arg("--json") + .arg(&out)); + let doc = results(&out); + assert_eq!(doc["tool"], format!("h5py-{executor}")); + assert_eq!(doc["results"].as_array().unwrap().len(), 8); + jsons.push(out); + } + let table = run(Command::new(python()) + .arg(scripts().join("compare_concurrent_read.py")) + .args(&jsons)); + assert!(table.contains("| deflate | same | 2 |"), "{table}"); + assert!(table.contains("clawhdf5 / h5py-processes"), "{table}"); + + // A different workload must not be compared. + let other = dir.path().join("other.json"); + run(Command::new(python()) + .arg(scripts().join("concurrent_read_h5py.py")) + .arg("--dir") + .arg(&data) + .args([ + "--threads", + "1", + "--slabs", + "4", + "--reps", + "1", + "--slab", + "64", + ]) + .arg("--json") + .arg(&other)); + let out = Command::new(python()) + .arg(scripts().join("compare_concurrent_read.py")) + .arg(&jsons[0]) + .arg(&other) + .output() + .unwrap(); + assert!(!out.status.success()); + assert!(String::from_utf8_lossy(&out.stderr).contains("slabs")); +}