commit 81f7bd84d0dd8bc25b9034959fb5fe0380bae4d9 Author: Owen Rummage Date: Fri Jul 17 13:34:44 2026 -0500 initial commit for alpha release and testing on MacOS diff --git a/Makefile b/Makefile new file mode 100644 index 0000000..e74f4f3 --- /dev/null +++ b/Makefile @@ -0,0 +1,152 @@ +# fossmark - multi-core CPU benchmark +# +# The assembly kernels are architecture-specific: +# src/fossmark.S AArch64 (ARM64) +# src/fossmark_x86_64.S x86-64 (AMD64) +# The C driver (src/main.c) is portable across architectures and OSes. A +# "binary that runs everywhere" is not possible - each OS/arch pair uses a +# different executable format and instruction set - so output is named per +# platform, e.g. dist/fossmark-linux-arm64, dist/fossmark-linux-amd64. +# +# Common targets: +# make build for the host arch (dist/fossmark--) +# make linux-arm64 build the Linux/ARM64 binary +# make linux-amd64 build the Linux/AMD64 binary +# make macos-arm64 build the macOS/ARM64 binary +# make macos-amd64 build the macOS/AMD64 binary +# make all build both Linux binaries +# make bench build for the host and run it +# make test build and run the kernel correctness tests (host arch) +# make clean remove dist/ +# +# Cross-compiling: linux-amd64 on an ARM64 host (or vice versa) needs the +# matching cross toolchain. The compiler for each target defaults to the host +# `cc` when the host arch already matches, and to the conventional GNU cross +# compiler otherwise. Override with CC_ARM64=... / CC_AMD64=... if your +# toolchain is named differently, e.g.: +# make linux-amd64 CC_AMD64=x86_64-linux-gnu-gcc-14 +# make linux-arm64 CC_ARM64="clang --target=aarch64-linux-gnu" +# +# On macOS, Apple Clang can build both architectures. The macOS compiler may +# be overridden for an osxcross or other cross toolchain: +# make macos-arm64 CC_MACOS_ARM64=clang +# make macos-amd64 CC_MACOS_AMD64=clang + +CC ?= cc +CFLAGS ?= -O2 -Wall -Wextra +LDLIBS ?= -lm +# The driver spreads each workload across all cores with pthreads. +PTHREAD := -pthread + +DIST := dist +DRIVER := src/main.c +ASM_ARM64 := src/fossmark.S +ASM_AMD64 := src/fossmark_x86_64.S + +# ---- host detection: normalise `uname -m` to our arch names ---- +HOST_ARCH := $(shell uname -m) +ifneq (,$(filter aarch64 arm64,$(HOST_ARCH))) + HOST_ARCHNAME := arm64 + HOST_ASM := $(ASM_ARM64) +else ifneq (,$(filter x86_64 amd64,$(HOST_ARCH))) + HOST_ARCHNAME := amd64 + HOST_ASM := $(ASM_AMD64) +else + HOST_ARCHNAME := $(HOST_ARCH) + HOST_ASM := $(ASM_ARM64) +endif + +# ---- host OS name for the native binary ---- +UNAME_S := $(shell uname -s) +ifeq ($(UNAME_S),Linux) + OSNAME := linux +else ifeq ($(UNAME_S),Darwin) + OSNAME := macos +else ifeq ($(OS),Windows_NT) + OSNAME := windows +else + OSNAME := $(shell uname -s | tr '[:upper:]' '[:lower:]') +endif + +# ---- per-target compilers: native cc if the host matches, else a cross gcc ---- +ifeq ($(HOST_ARCHNAME),arm64) + CC_ARM64 ?= $(CC) +else + CC_ARM64 ?= aarch64-linux-gnu-gcc +endif +ifeq ($(HOST_ARCHNAME),amd64) + CC_AMD64 ?= $(CC) +else + CC_AMD64 ?= x86_64-linux-gnu-gcc +endif +CC_MACOS_ARM64 ?= $(CC) +CC_MACOS_AMD64 ?= $(CC) + +NATIVE_BIN := $(DIST)/fossmark-$(OSNAME)-$(HOST_ARCHNAME) + +# `make` with no target builds the host binary, as before. +.DEFAULT_GOAL := native +.PHONY: all native linux-arm64 linux-amd64 macos-arm64 macos-amd64 bench test clean + +# `make all` builds both Linux binaries. +all: linux-arm64 linux-amd64 + +# `make native` (and bare `make`) build for whatever host you are on. +native: $(NATIVE_BIN) + +linux-arm64: $(DIST)/fossmark-linux-arm64 +linux-amd64: $(DIST)/fossmark-linux-amd64 +macos-arm64: $(DIST)/fossmark-macos-arm64 +macos-amd64: $(DIST)/fossmark-macos-amd64 + +$(DIST)/fossmark-linux-arm64: $(DRIVER) $(ASM_ARM64) | $(DIST) + $(CC_ARM64) $(CFLAGS) $(PTHREAD) -o $@ $(DRIVER) $(ASM_ARM64) $(LDLIBS) + @echo "built $@" + +$(DIST)/fossmark-linux-amd64: $(DRIVER) $(ASM_AMD64) | $(DIST) + $(CC_AMD64) $(CFLAGS) $(PTHREAD) -o $@ $(DRIVER) $(ASM_AMD64) $(LDLIBS) + @echo "built $@" + +$(DIST)/fossmark-macos-arm64: $(DRIVER) $(ASM_ARM64) | $(DIST) + $(CC_MACOS_ARM64) -arch arm64 $(CFLAGS) $(PTHREAD) -o $@ $(DRIVER) $(ASM_ARM64) $(LDLIBS) + @echo "built $@" + +$(DIST)/fossmark-macos-amd64: $(DRIVER) $(ASM_AMD64) | $(DIST) + $(CC_MACOS_AMD64) -arch x86_64 $(CFLAGS) $(PTHREAD) -o $@ $(DRIVER) $(ASM_AMD64) $(LDLIBS) + @echo "built $@" + +# When the host is Linux/ARM64 or Linux/AMD64, the native binary IS one of the +# linux-* targets above, so no separate recipe is defined (that would be a +# duplicate). Otherwise - e.g. macOS/ARM64 - provide the native recipe here. +ifeq ($(OSNAME)-$(HOST_ARCHNAME),linux-arm64) +NATIVE_HAS_RULE := yes +endif +ifeq ($(OSNAME)-$(HOST_ARCHNAME),linux-amd64) +NATIVE_HAS_RULE := yes +endif +ifeq ($(OSNAME)-$(HOST_ARCHNAME),macos-arm64) +NATIVE_HAS_RULE := yes +endif +ifeq ($(OSNAME)-$(HOST_ARCHNAME),macos-amd64) +NATIVE_HAS_RULE := yes +endif +ifneq ($(NATIVE_HAS_RULE),yes) +$(NATIVE_BIN): $(DRIVER) $(HOST_ASM) | $(DIST) + $(CC) $(CFLAGS) $(PTHREAD) -o $@ $(DRIVER) $(HOST_ASM) $(LDLIBS) + @echo "built $@" +endif + +$(DIST): + mkdir -p $(DIST) + +# Build for the host and run the benchmark. +bench: $(NATIVE_BIN) + ./$(NATIVE_BIN) + +# Build and run the kernel correctness tests for the host arch. +test: | $(DIST) + $(CC) $(CFLAGS) $(PTHREAD) -o $(DIST)/test_kernels src/test_kernels.c $(HOST_ASM) $(LDLIBS) + ./$(DIST)/test_kernels + +clean: + rm -rf $(DIST) diff --git a/README.md b/README.md new file mode 100644 index 0000000..b0e0c20 --- /dev/null +++ b/README.md @@ -0,0 +1,124 @@ +# Fossmark + +A single-threaded CPU benchmark for ARM64 (AArch64), with the numeric kernels +hand-written in assembly and a small portable C driver to run and score them. + +## What it measures + +Nine workloads, each a tight assembly kernel: + +| # | Test | What it exercises | +|---|-------------------------|--------------------------------------------------------------| +| 1 | Integer Math | 64-bit ALU: `madd`, `umulh`/`smulh`, `udiv`/`sdiv`, bit ops | +| 2 | Floating Point Math | scalar double: `fmadd`, `fdiv`, `fsqrt` | +| 3 | Prime Numbers | sieve of Eratosthenes to 2,000,000 (strided memory + ALU) | +| 4 | Extended Instructions | NEON/ASIMD: 128-bit integer, widening, table, float vectors | +| 5 | Compression | LZ77 match-finder over a 4 MiB corpus (branchy, cache probe) | +| 6 | Encryption | ChaCha20, 20 rounds, NEON, over 1 MiB | +| 7 | Physics | 512-body direct-summation gravity, double precision | +| 8 | Sorting | in-place heapsort of 1M `uint32` (branch + cache stress) | +| 9 | Single-Threaded | dependent-load pointer chase over 16 MiB (memory latency) | + +Each test auto-calibrates its iteration count until it runs long enough to be +timed reliably, then reports the best of several runs (the run least disturbed +by the OS scheduler). Every kernel returns a checksum that the driver verifies +across runs, so a miscompiled or non-deterministic kernel is caught rather than +silently mis-scored. + +## Scoring + +Each test's raw rate is normalised against a **reference machine** into a +unitless score, and the overall is a **weighted geometric mean** of those +scores: + +``` +S_i = TARGET * (rate_i / REF_i) (per-test score) +Overall = TARGET * exp( Σ w_i·ln(rate_i/REF_i) / Σ w_i ) (weighted geo. mean) +``` + +The reference rates are the tuning machine's own rates, and `TARGET` is 20000, +so that machine scores ~20000 on every test and overall. Scaling is linear in +performance: a machine half as fast scores ~10000, one 10× slower ~2000, and a +future machine twice as fast ~40000 — so there is unbounded room both below and +above the reference. + +The weights reflect each test's influence on **everyday, common-workload user +experience** — integer/general-purpose throughput and memory-latency-bound +responsiveness matter most; specialised floating-point and physics matter +least. This mirrors the weighted, integer-dominant approach of mainstream +suites such as Geekbench 6 (which splits integer/FP roughly 65/35 and combines +real-world workloads with a weighted mean). + +| Test | Weight | +|---|---:| +| Integer Math | 20% | +| Single-Threaded | 16% | +| Compression | 14% | +| Sorting | 12% | +| Extended Instructions | 11% | +| Floating Point | 9% | +| Encryption | 8% | +| Prime Numbers | 6% | +| Physics | 4% | + +Everything above is configurable via `#define`s at the top of `src/main.c`: +`FM_TARGET_SCORE`, the nine `FM_REF_*` reference rates, and the nine +`FM_WEIGHT_*` weights. Weights are relative — the code normalises by their sum, +so you can change one without rebalancing the rest. To re-baseline for a +different reference machine, set each `FM_REF_*` to that machine's measured +rate. + +Note: the pointer-chase (Single-Threaded) test measures raw memory latency and +is the noisiest to sample, so the overall typically varies ~1–2% run to run. + +## "Runs on all operating systems" + +The **assembly is** OS-independent: `src/fossmark.S` contains no system calls, +no libc calls, and no external relocations. Every routine is a pure function of +its arguments under the AAPCS64 calling convention, so the same source +assembles and runs correctly on Linux (ELF), macOS (Mach-O), Windows (COFF) and +the BSDs. It avoids `x18` (reserved on Darwin/Windows) and the `v8`–`v15` +callee-saved vector bank. + +A single *binary* that runs everywhere is not possible — Linux, macOS and +Windows use incompatible executable formats and system-call ABIs. So the +portable C driver (`src/main.c`) supplies the per-OS parts (timing, memory, +I/O), and you build one binary per platform. The Linux build is named +`fossmark-linux-arm64`. + +## Build + +```sh +make # builds dist/fossmark-- for the host +make linux-arm64 +make linux-amd64 +make macos-arm64 +make macos-amd64 +make bench # build and run the benchmark +make test # build and run the kernel correctness tests +``` + +Both macOS targets can be built on either Apple Silicon or Intel Macs; Apple +Clang selects the requested architecture with `-arch`. They produce +`dist/fossmark-macos-arm64` and `dist/fossmark-macos-amd64`, respectively. + +Or by hand: + +```sh +cc -O2 src/main.c src/fossmark.S -o dist/fossmark-linux-arm64 -lm +``` + +On macOS the same command produces a native binary (name it +`fossmark-macos-arm64`); on Windows use `clang` from the LLVM/MSVC toolchain. + +## Testing + +`src/test_kernels.c` is a standalone harness that validates each kernel against +an independent reference or invariant — the sieve against a C reference sieve, +the NEON ChaCha20 against a scalar reference anchored to the RFC 8439 +known-answer vector, the sort against `qsort`, the N-body step against +conservation of momentum, and so on. It exits non-zero if any check fails. + +```sh +make test +``` diff --git a/dist/fossmark-linux-amd64 b/dist/fossmark-linux-amd64 new file mode 100755 index 0000000..94308bd Binary files /dev/null and b/dist/fossmark-linux-amd64 differ diff --git a/dist/fossmark-linux-arm64 b/dist/fossmark-linux-arm64 new file mode 100755 index 0000000..18e54ba Binary files /dev/null and b/dist/fossmark-linux-arm64 differ diff --git a/dist/test_kernels b/dist/test_kernels new file mode 100755 index 0000000..0abdfd4 Binary files /dev/null and b/dist/test_kernels differ diff --git a/src/fossmark.S b/src/fossmark.S new file mode 100644 index 0000000..f7d01bf --- /dev/null +++ b/src/fossmark.S @@ -0,0 +1,942 @@ +/* + * fossmark.S - AArch64 CPU benchmark kernels + * + * OS-independent: contains no syscalls, no libc calls, no relocations against + * external data. Every routine is a pure function of its arguments under the + * AAPCS64 procedure call standard, so this source assembles and runs correctly + * on Linux (ELF), macOS (Mach-O), Windows (COFF), and the BSDs. + * + * x18 is never used: it is the platform register on Darwin and Windows. + * v8-v15 are never used: only their low 64 bits are callee-saved, which makes + * them a trap for 128-bit vector code. + */ + + .arch armv8-a + +#if defined(__APPLE__) +# define SYM(name) _##name +#else +# define SYM(name) name +#endif + +#if defined(__ELF__) +# define FN_BEGIN(name) .p2align 4 ; .globl SYM(name) ; .type SYM(name), %function ; SYM(name): +# define FN_END(name) .size SYM(name), . - SYM(name) +#else +# define FN_BEGIN(name) .p2align 4 ; .globl SYM(name) ; SYM(name): +# define FN_END(name) +#endif + +/* Materialise a 64-bit constant without a literal pool, so no relocations and + * no .ltorg placement worries across object formats. */ +#define MOV64(reg, val) \ + movz reg, #((val) & 0xffff) ;\ + movk reg, #(((val) >> 16) & 0xffff), lsl #16 ;\ + movk reg, #(((val) >> 32) & 0xffff), lsl #32 ;\ + movk reg, #(((val) >> 48) & 0xffff), lsl #48 + + .text + +/* =================================================================== + * uint64_t fm_int_math(uint64_t iters) + * + * Four largely independent accumulator chains to expose instruction-level + * parallelism, mixed with high-latency serialising ops (udiv/sdiv) and the + * bit-manipulation instructions. Returns a checksum so the compiler and the + * driver cannot elide the work. + * =================================================================== */ +FN_BEGIN(fm_int_math) + cbz x0, .Lim_zero + + MOV64(x1, 0x9E3779B97F4A7C15) /* a */ + MOV64(x2, 0xBF58476D1CE4E5B9) /* b */ + MOV64(x3, 0x94D049BB133111EB) /* c */ + MOV64(x4, 0x2545F4914F6CDD1D) /* d */ + MOV64(x5, 0x00000000DEADBEEF) /* odd multiplier, never zero */ + mov x6, x0 /* trip count */ + +.Lim_loop: + /* four independent multiply-accumulate chains */ + madd x1, x1, x5, x2 + madd x2, x2, x5, x3 + madd x3, x3, x5, x4 + madd x4, x4, x5, x1 + + /* cross-mix with shifts and logic ops (free shifter operands) */ + eor x1, x1, x3, lsr #29 + eor x2, x2, x4, lsl #17 + eor x3, x3, x1, ror #31 + bic x4, x4, x2, asr #7 + + /* wide multiplies: umulh/smulh are the long-latency multiplier path */ + umulh x9, x1, x3 + smulh x10, x2, x4 + add x1, x1, x9 + add x2, x2, x10 + + /* bit manipulation */ + rbit x11, x1 + clz x12, x2 + rev x13, x3 + eor x4, x4, x11 + add x4, x4, x12 + eor x1, x1, x13 + + /* division: fully serialising, ~10-20 cycle latency, not pipelined */ + orr x14, x5, #1 /* guarantee a non-zero divisor */ + udiv x15, x1, x14 + sdiv x16, x2, x14 + msub x3, x15, x14, x3 + add x4, x4, x16 + + /* bitfield ops */ + ror x2, x2, #11 + extr x1, x1, x2, #23 + + subs x6, x6, #1 + b.ne .Lim_loop + + eor x0, x1, x2 + eor x0, x0, x3 + eor x0, x0, x4 + ret + +.Lim_zero: + mov x0, xzr + ret +FN_END(fm_int_math) + + +/* =================================================================== + * uint64_t fm_fp_math(uint64_t iters) + * + * Double-precision scalar FP. Four fmadd chains cover the pipelined + * multiply-add path; fdiv and fsqrt cover the non-pipelined divide/sqrt unit, + * which is usually the real differentiator between cores. + * Returns the result bit-cast to u64. + * =================================================================== */ +FN_BEGIN(fm_fp_math) + cbz x0, .Lfp_zero + mov x6, x0 + + /* Constants come from a table rather than fmov immediates: the AArch64 + * 8-bit FP immediate can only encode a narrow set of values, and + * several of the ones we want fall outside it. */ + adr x7, .Lfp_consts + ldp d0, d1, [x7] /* a = 1.5, b = 2.5 */ + ldp d2, d3, [x7, #16] /* c = 3.5, d = 0.5 */ + ldp d4, d5, [x7, #32] /* mul, small addend */ + ldp d6, d7, [x7, #48] /* 2.0, 1.0 */ + +.Lfp_loop: + /* four independent fused multiply-add chains */ + fmadd d0, d0, d4, d5 + fmadd d1, d1, d4, d5 + fmadd d2, d2, d4, d5 + fmadd d3, d3, d4, d5 + + /* keep the accumulators bounded so they never reach inf/NaN */ + fmin d0, d0, d6 + fmin d1, d1, d6 + fmin d2, d2, d6 + fmin d3, d3, d6 + + /* square root: long latency, low throughput */ + fsqrt d16, d0 + fsqrt d17, d1 + fadd d2, d2, d16 + fadd d3, d3, d17 + + /* divide: the other long-latency unit */ + fadd d18, d2, d7 /* divisor >= 1, never zero */ + fdiv d19, d7, d18 + fadd d0, d0, d19 + + fadd d20, d3, d7 + fdiv d21, d7, d20 + fadd d1, d1, d21 + + /* abs/neg/compare-select: cheap ops to balance the mix */ + fabs d2, d2 + fneg d22, d3 + fabs d3, d22 + fmax d3, d3, d7 + + subs x6, x6, #1 + b.ne .Lfp_loop + + fadd d0, d0, d1 + fadd d2, d2, d3 + fadd d0, d0, d2 + fmov x0, d0 + ret + +.Lfp_zero: + mov x0, xzr + ret +FN_END(fm_fp_math) + + .p2align 4 +.Lfp_consts: + .double 1.5, 2.5 + .double 3.5, 0.5 + .double 1.0625, 0.0009765625 + .double 2.0, 1.0 + + +/* =================================================================== + * uint64_t fm_primes(uint64_t limit, uint8_t *sieve) + * + * Sieve of Eratosthenes over [0, limit). The caller supplies `limit` bytes of + * scratch; this routine clears it itself, so the clearing pass is part of the + * measured work (as it would be in any real use). Returns the prime count. + * + * Strided stores over a buffer larger than L1 make this a memory-hierarchy + * test as much as an arithmetic one. + * =================================================================== */ +FN_BEGIN(fm_primes) + cmp x0, #2 + b.lo .Lpr_none + + mov x2, x1 /* sieve base */ + mov x3, x0 /* limit */ + + /* zero the sieve, 32 bytes per iteration */ + movi v0.16b, #0 + mov x4, xzr + and x5, x3, #~31 /* bulk portion */ +.Lpr_clear32: + cmp x4, x5 + b.hs .Lpr_clear1 + add x6, x2, x4 + stp q0, q0, [x6] + add x4, x4, #32 + b .Lpr_clear32 +.Lpr_clear1: + cmp x4, x3 + b.hs .Lpr_clear_done + strb wzr, [x2, x4] + add x4, x4, #1 + b .Lpr_clear1 +.Lpr_clear_done: + + /* mark 0 and 1 as composite */ + mov w6, #1 + strb w6, [x2] + strb w6, [x2, #1] + + /* outer loop: i = 2; i*i < limit; i++ */ + mov x7, #2 +.Lpr_outer: + mul x9, x7, x7 + cmp x9, x3 + b.hs .Lpr_count + + ldrb w10, [x2, x7] + cbnz w10, .Lpr_outer_next /* already composite, skip */ + + /* inner loop: mark multiples starting at i*i, stride i */ + mov x11, x9 +.Lpr_inner: + cmp x11, x3 + b.hs .Lpr_outer_next + strb w6, [x2, x11] + add x11, x11, x7 + b .Lpr_inner + +.Lpr_outer_next: + add x7, x7, #1 + b .Lpr_outer + + /* count the survivors */ +.Lpr_count: + mov x0, xzr /* count */ + mov x4, #2 +.Lpr_count_loop: + cmp x4, x3 + b.hs .Lpr_done + ldrb w10, [x2, x4] + cmp w10, #0 + cinc x0, x0, eq + add x4, x4, #1 + b .Lpr_count_loop + +.Lpr_done: + ret + +.Lpr_none: + mov x0, xzr + ret +FN_END(fm_primes) + + +/* =================================================================== + * uint64_t fm_simd(uint64_t iters, void *buf) + * + * "Extended instructions": the ASIMD/NEON unit, which is architecturally + * mandatory on AArch64 and therefore safe to use without runtime feature + * detection. (AES, SHA and DotProd are all *optional* extensions; using them + * unguarded would fault on cores that lack them, so they are deliberately + * avoided here.) + * + * buf must be at least 128 bytes and 16-byte aligned. Returns a checksum. + * =================================================================== */ +FN_BEGIN(fm_simd) + cbz x0, .Lsd_zero + mov x6, x0 + + /* seed eight vectors from the scratch buffer */ + ldp q0, q1, [x1] + ldp q2, q3, [x1, #32] + ldp q4, q5, [x1, #64] + ldp q6, q7, [x1, #96] + + movi v28.4s, #3 + movi v29.4s, #7 + movi v30.16b, #0x5A + + /* byte-permute table for tbl: reverses the 16 byte lanes */ + adr x9, .Lsd_perm + ldr q31, [x9] + + /* float accumulator and operands: derived from the integer seeds by + * conversion, so the FP pipeline sees real finite values rather than + * reinterpreted integer bit patterns (which would be NaNs/denormals + * and would measure the slow path, not the common one). */ + movi v22.4s, #0 + scvtf v26.4s, v0.4s + scvtf v27.4s, v1.4s + +.Lsd_loop: + /* integer SIMD: multiply-accumulate across four independent vectors */ + mla v0.4s, v1.4s, v28.4s + mla v1.4s, v2.4s, v29.4s + mla v2.4s, v3.4s, v28.4s + mla v3.4s, v0.4s, v29.4s + + /* saturating and halving arithmetic */ + sqadd v4.4s, v4.4s, v0.4s + uhadd v5.4s, v5.4s, v1.4s + srhadd v6.4s, v6.4s, v2.4s + sqsub v7.4s, v7.4s, v3.4s + + /* widening ops: 16->32 bit lane expansion */ + umull v16.4s, v0.4h, v1.4h + umull2 v17.4s, v0.8h, v1.8h + saddw v2.4s, v2.4s, v16.4h + ssubw2 v3.4s, v3.4s, v17.8h + + /* pairwise reduction */ + uaddlp v18.2d, v4.4s + addp v19.4s, v5.4s, v6.4s + add v4.4s, v4.4s, v19.4s + + /* table lookup: the byte-permute network */ + tbl v20.16b, {v0.16b}, v31.16b + eor v1.16b, v1.16b, v20.16b + + /* shifts, logic, min/max, reverse */ + shl v21.4s, v2.4s, #3 + usra v21.4s, v2.4s, #29 + orr v2.16b, v2.16b, v21.16b + bic v3.16b, v3.16b, v30.16b + umax v5.4s, v5.4s, v0.4s + umin v6.4s, v6.4s, v1.4s + rev32 v7.16b, v7.16b + + /* single-precision float SIMD: 4-wide fmla, plus the reciprocal and + * rsqrt estimate instructions that shader-style code leans on */ + /* v16/v17 are dead after the widening ops above, so they are reused + * here as scratch rather than touching the callee-saved v8-v15 bank. */ + fmla v22.4s, v26.4s, v27.4s + frecpe v23.4s, v26.4s + frsqrte v16.4s, v22.4s + fadd v22.4s, v22.4s, v23.4s + fmul v26.4s, v26.4s, v16.4s + + /* population count and leading-zero count */ + cnt v24.16b, v0.16b + uaddlv h25, v24.8b + clz v17.4s, v1.4s + add v0.4s, v0.4s, v17.4s + + /* fold the reduction results back in so nothing is dead code */ + add v4.2d, v4.2d, v18.2d + + subs x6, x6, #1 + b.ne .Lsd_loop + + /* horizontal fold to a single 64-bit checksum */ + eor v0.16b, v0.16b, v1.16b + eor v2.16b, v2.16b, v3.16b + eor v4.16b, v4.16b, v5.16b + eor v6.16b, v6.16b, v7.16b + eor v0.16b, v0.16b, v2.16b + eor v4.16b, v4.16b, v6.16b + eor v0.16b, v0.16b, v4.16b + addv s0, v0.4s + fmov w0, s0 + ret + +.Lsd_zero: + mov x0, xzr + ret +FN_END(fm_simd) + + .p2align 4 +.Lsd_perm: + .byte 15,14,13,12,11,10,9,8,7,6,5,4,3,2,1,0 + + +/* =================================================================== + * uint64_t fm_compress(const uint8_t *src, uint64_t len, uint32_t *ht) + * + * The match-finding inner loop of an LZ77 compressor (the LZ4 fast strategy): + * hash the next 4 bytes, probe a single-entry-per-bucket table, verify, then + * extend the match. This is where real compressors spend their time - it is + * branch-heavy with a data-dependent, cache-missing table probe. + * + * ht must hold 1<<16 uint32_t (256 KiB); this routine clears it itself. + * Returns the encoded size in bytes. + * =================================================================== */ +FN_BEGIN(fm_compress) + stp x29, x30, [sp, #-96]! + mov x29, sp + stp x19, x20, [sp, #16] + stp x21, x22, [sp, #32] + stp x23, x24, [sp, #48] + stp x25, x26, [sp, #64] + stp x27, x28, [sp, #80] + + mov x19, x0 /* src */ + mov x20, x1 /* len */ + mov x21, x2 /* ht */ + + /* clear the hash table: 1<<16 entries * 4 bytes = 262144 bytes */ + movi v0.16b, #0 + mov x9, xzr + MOV64(x10, 262144) +.Lcm_clear: + add x11, x21, x9 + stp q0, q0, [x11] + stp q0, q0, [x11, #32] + add x9, x9, #64 + cmp x9, x10 + b.lo .Lcm_clear + + cmp x20, #16 + b.lo .Lcm_tiny + + mov x22, x19 /* ip */ + mov x23, x19 /* anchor */ + add x24, x19, x20 /* end */ + sub x25, x24, #12 /* mflimit */ + mov x26, xzr /* outsize */ + + MOV64(x27, 2654435761) /* Knuth multiplicative hash */ + +.Lcm_loop: + cmp x22, x25 + b.hs .Lcm_flush + + ldr w9, [x22] /* seq = load32(ip) */ + mul w10, w9, w27 + lsr w10, w10, #16 /* h = (seq * prime) >> 16 */ + + ldr w11, [x21, x10, lsl #2] /* ref_off = ht[h] */ + sub x12, x22, x19 /* cur_off = ip - src */ + str w12, [x21, x10, lsl #2] /* ht[h] = cur_off */ + + add x13, x19, x11 /* ref = src + ref_off */ + cmp x13, x22 + b.hs .Lcm_no_match /* ref must be strictly behind ip */ + + sub x14, x22, x13 /* distance */ + MOV64(x15, 65536) + cmp x14, x15 + b.hs .Lcm_no_match /* 16-bit offset window */ + + ldr w16, [x13] + cmp w16, w9 + b.ne .Lcm_no_match + + /* match confirmed: extend it byte by byte */ + mov x28, #4 /* ml */ +.Lcm_extend: + add x9, x22, x28 + cmp x9, x24 + b.hs .Lcm_emit + ldrb w10, [x22, x28] + ldrb w11, [x13, x28] + cmp w10, w11 + b.ne .Lcm_emit + add x28, x28, #1 + b .Lcm_extend + +.Lcm_emit: + /* token(1) + offset(2) + literals + varint extensions */ + sub x9, x22, x23 /* literal run length */ + add x26, x26, x9 + add x26, x26, #3 + cmp x9, #15 + cinc x26, x26, hs /* literal-length extension byte */ + cmp x28, #19 + cinc x26, x26, hs /* match-length extension byte */ + + add x22, x22, x28 + mov x23, x22 + b .Lcm_loop + +.Lcm_no_match: + add x22, x22, #1 + b .Lcm_loop + +.Lcm_flush: + /* trailing literals */ + sub x9, x24, x23 + add x26, x26, x9 + add x26, x26, #1 + mov x0, x26 + b .Lcm_ret + +.Lcm_tiny: + add x0, x20, #1 + +.Lcm_ret: + ldp x27, x28, [sp, #80] + ldp x25, x26, [sp, #64] + ldp x23, x24, [sp, #48] + ldp x21, x22, [sp, #32] + ldp x19, x20, [sp, #16] + ldp x29, x30, [sp], #96 + ret +FN_END(fm_compress) + + +/* =================================================================== + * uint64_t fm_chacha20(uint8_t *buf, uint64_t len, const uint8_t key[32], + * uint64_t rounds) + * + * ChaCha20 stream cipher, NEON, four 128-bit state rows. Chosen over AES + * deliberately: the ARMv8 AES extension is optional, so an AES-instruction + * benchmark would fault on cores without it. ChaCha20 needs only baseline + * ASIMD and is a real, widely deployed cipher (TLS, WireGuard, SSH). + * + * len is rounded down to a multiple of 64. `rounds` = number of passes over + * the buffer. Returns a checksum of the keystream output. + * =================================================================== */ + +/* rotate each 32-bit lane left by n, via shl + shift-right-and-insert */ +#define VROTL(vd, vs, vt, n) \ + shl vt##.4s, vs##.4s, #(n) ;\ + sri vt##.4s, vs##.4s, #(32 - (n)) ;\ + mov vd##.16b, vt##.16b + +/* one ChaCha quarter-round over rows a,b,c,d using v24 as scratch */ +#define QROUND(a, b, c, d) \ + add a##.4s, a##.4s, b##.4s ;\ + eor d##.16b, d##.16b, a##.16b ;\ + rev32 d##.8h, d##.8h ;\ + add c##.4s, c##.4s, d##.4s ;\ + eor b##.16b, b##.16b, c##.16b ;\ + VROTL(b, b, v24, 12) ;\ + add a##.4s, a##.4s, b##.4s ;\ + eor d##.16b, d##.16b, a##.16b ;\ + VROTL(d, d, v24, 8) ;\ + add c##.4s, c##.4s, d##.4s ;\ + eor b##.16b, b##.16b, c##.16b ;\ + VROTL(b, b, v24, 7) + +FN_BEGIN(fm_chacha20) + and x1, x1, #~63 /* whole 64-byte blocks only */ + cbz x1, .Lcc_zero + cbz x3, .Lcc_zero + + stp x29, x30, [sp, #-32]! + mov x29, sp + stp x19, x20, [sp, #16] + + mov x19, x0 /* buf */ + mov x20, x1 /* len */ + + /* v4..v7 hold the base state */ + adr x9, .Lcc_sigma + ldr q4, [x9] /* "expand 32-byte k" */ + ldp q5, q6, [x2] /* key[0..31] */ + movi v7.4s, #0 /* counter || nonce */ + + movi v25.16b, #0 /* running checksum */ + mov x10, xzr /* counter value */ + +.Lcc_pass: + mov x11, xzr /* byte offset into buf */ + +.Lcc_block: + /* working state = base state, with the block counter in lane 0 of v7 */ + mov v0.16b, v4.16b + mov v1.16b, v5.16b + mov v2.16b, v6.16b + mov v3.16b, v7.16b + mov v3.s[0], w10 + + /* keep originals for the final feed-forward add */ + mov v16.16b, v0.16b + mov v17.16b, v1.16b + mov v18.16b, v2.16b + mov v19.16b, v3.16b + + mov w12, #10 /* 10 double rounds = 20 rounds */ +.Lcc_rounds: + /* column round */ + QROUND(v0, v1, v2, v3) + + /* rotate lanes to form the diagonals */ + ext v1.16b, v1.16b, v1.16b, #4 + ext v2.16b, v2.16b, v2.16b, #8 + ext v3.16b, v3.16b, v3.16b, #12 + + /* diagonal round */ + QROUND(v0, v1, v2, v3) + + /* undo the lane rotation */ + ext v1.16b, v1.16b, v1.16b, #12 + ext v2.16b, v2.16b, v2.16b, #8 + ext v3.16b, v3.16b, v3.16b, #4 + + subs w12, w12, #1 + b.ne .Lcc_rounds + + /* feed-forward: keystream = working + original */ + add v0.4s, v0.4s, v16.4s + add v1.4s, v1.4s, v17.4s + add v2.4s, v2.4s, v18.4s + add v3.4s, v3.4s, v19.4s + + /* XOR the keystream into the buffer */ + add x13, x19, x11 + ldp q20, q21, [x13] + ldp q22, q23, [x13, #32] + eor v20.16b, v20.16b, v0.16b + eor v21.16b, v21.16b, v1.16b + eor v22.16b, v22.16b, v2.16b + eor v23.16b, v23.16b, v3.16b + stp q20, q21, [x13] + stp q22, q23, [x13, #32] + + /* accumulate a checksum of the keystream */ + eor v25.16b, v25.16b, v0.16b + eor v25.16b, v25.16b, v3.16b + + add x10, x10, #1 + add x11, x11, #64 + cmp x11, x20 + b.lo .Lcc_block + + subs x3, x3, #1 + b.ne .Lcc_pass + + addv s25, v25.4s + fmov w0, s25 + + ldp x19, x20, [sp, #16] + ldp x29, x30, [sp], #32 + ret + +.Lcc_zero: + mov x0, xzr + ret +FN_END(fm_chacha20) + + .p2align 4 +.Lcc_sigma: + .word 0x61707865, 0x3320646e, 0x79622d32, 0x6b206574 + + +/* =================================================================== + * uint64_t fm_physics(double *bodies, uint64_t n, uint64_t steps) + * + * Direct-summation N-body gravity, O(n^2) per step, double precision. + * Layout per body, 8 doubles (64 bytes, one cache line): + * [0]=x [1]=y [2]=z [3]=mass [4]=vx [5]=vy [6]=vz [7]=pad + * + * The 1/sqrt is done with a real fsqrt+fdiv rather than the frsqrte estimate, + * so this exercises the divide/sqrt unit the way physics code actually does. + * Returns a checksum bit-cast from the final velocity sum. + * =================================================================== */ +FN_BEGIN(fm_physics) + cbz x1, .Lph_zero + cbz x2, .Lph_zero + + stp x29, x30, [sp, #-64]! + mov x29, sp + stp x19, x20, [sp, #16] + stp x21, x22, [sp, #32] + stp x23, x24, [sp, #48] + + mov x19, x0 /* bodies */ + mov x20, x1 /* n */ + mov x21, x2 /* steps */ + + adr x9, .Lph_consts + ldp d28, d29, [x9] /* dt, eps^2 */ + ldr d30, [x9, #16] /* 1.0 */ + +.Lph_step: + mov x22, xzr /* i */ + +.Lph_body_i: + lsl x9, x22, #6 /* i * 64 */ + add x23, x19, x9 /* &bodies[i] */ + + ldp d0, d1, [x23] /* xi, yi */ + ldr d2, [x23, #16] /* zi */ + + movi d16, #0 /* ax */ + movi d17, #0 /* ay */ + movi d18, #0 /* az */ + + mov x24, xzr /* j */ + mov x10, x19 /* &bodies[j] */ + +.Lph_body_j: + ldp d3, d4, [x10] /* xj, yj */ + ldp d5, d6, [x10, #16] /* zj, mj */ + + fsub d3, d3, d0 /* dx */ + fsub d4, d4, d1 /* dy */ + fsub d5, d5, d2 /* dz */ + + /* d2 = dx*dx + dy*dy + dz*dz + eps^2 (always >= eps^2, never zero, + * so the i==j self-term is finite and contributes exactly 0 below) */ + fmul d7, d3, d3 + fmadd d7, d4, d4, d7 + fmadd d7, d5, d5, d7 + fadd d7, d7, d29 + + fsqrt d19, d7 /* r */ + fdiv d20, d30, d19 /* 1/r */ + fmul d21, d20, d20 /* 1/r^2 */ + fmul d21, d21, d20 /* 1/r^3 */ + fmul d21, d21, d6 /* m/r^3 */ + + fmadd d16, d3, d21, d16 /* ax += dx * m/r^3 */ + fmadd d17, d4, d21, d17 + fmadd d18, d5, d21, d18 + + add x10, x10, #64 + add x24, x24, #1 + cmp x24, x20 + b.lo .Lph_body_j + + /* v += a * dt */ + ldp d22, d23, [x23, #32] + ldr d24, [x23, #48] + fmadd d22, d16, d28, d22 + fmadd d23, d17, d28, d23 + fmadd d24, d18, d28, d24 + stp d22, d23, [x23, #32] + str d24, [x23, #48] + + add x22, x22, #1 + cmp x22, x20 + b.lo .Lph_body_i + + /* second pass: x += v * dt (positions updated only after all forces) */ + mov x22, xzr + mov x10, x19 +.Lph_integrate: + ldp d0, d1, [x10] + ldr d2, [x10, #16] + ldp d22, d23, [x10, #32] + ldr d24, [x10, #48] + fmadd d0, d22, d28, d0 + fmadd d1, d23, d28, d1 + fmadd d2, d24, d28, d2 + stp d0, d1, [x10] + str d2, [x10, #16] + add x10, x10, #64 + add x22, x22, #1 + cmp x22, x20 + b.lo .Lph_integrate + + subs x21, x21, #1 + b.ne .Lph_step + + /* checksum: sum of all velocity components */ + movi d0, #0 + mov x22, xzr + mov x10, x19 +.Lph_sum: + ldp d22, d23, [x10, #32] + ldr d24, [x10, #48] + fadd d0, d0, d22 + fadd d0, d0, d23 + fadd d0, d0, d24 + add x10, x10, #64 + add x22, x22, #1 + cmp x22, x20 + b.lo .Lph_sum + + fmov x0, d0 + + ldp x23, x24, [sp, #48] + ldp x21, x22, [sp, #32] + ldp x19, x20, [sp, #16] + ldp x29, x30, [sp], #64 + ret + +.Lph_zero: + mov x0, xzr + ret +FN_END(fm_physics) + + .p2align 4 +.Lph_consts: + .double 0.0078125 /* dt */ + .double 0.0625 /* eps^2 */ + .double 1.0 + + +/* =================================================================== + * uint64_t fm_sort(uint32_t *a, uint64_t n) + * + * In-place heapsort. Chosen over quicksort because it needs no recursion or + * explicit stack, yet is aggressively branch-unpredictable and touches memory + * in a scattered pattern - it stresses the branch predictor and the cache + * hierarchy, which is what a sort benchmark should measure. + * + * Returns an order-sensitive checksum, which also verifies the sort. + * =================================================================== */ +FN_BEGIN(fm_sort) + cmp x1, #2 + b.lo .Lst_trivial + + stp x29, x30, [sp, #-48]! + mov x29, sp + stp x19, x20, [sp, #16] + stp x21, x22, [sp, #32] + + mov x19, x0 /* a */ + mov x20, x1 /* n */ + + /* ---- build the max-heap: for i = n/2 - 1 down to 0 ---- */ + lsr x21, x20, #1 /* i = n/2 */ +.Lst_build: + cbz x21, .Lst_extract + sub x21, x21, #1 /* i-- */ + mov x0, x21 /* root */ + mov x1, x20 /* end */ + bl .Lst_siftdown + cbnz x21, .Lst_build + + /* ---- extract: for end = n-1 down to 1 ---- */ +.Lst_extract: + sub x22, x20, #1 /* end = n-1 */ +.Lst_extract_loop: + cbz x22, .Lst_checksum + + /* swap a[0] and a[end] */ + ldr w9, [x19] + ldr w10, [x19, x22, lsl #2] + str w10, [x19] + str w9, [x19, x22, lsl #2] + + mov x0, xzr /* root = 0 */ + mov x1, x22 /* end = end */ + bl .Lst_siftdown + + sub x22, x22, #1 + b .Lst_extract_loop + + /* ---- order-sensitive checksum ---- */ +.Lst_checksum: + mov x0, xzr + mov x9, xzr +.Lst_cksum_loop: + ldr w10, [x19, x9, lsl #2] + eor x0, x0, x10 + ror x0, x0, #7 + add x0, x0, x10 + add x9, x9, #1 + cmp x9, x20 + b.lo .Lst_cksum_loop + + ldp x21, x22, [sp, #32] + ldp x19, x20, [sp, #16] + ldp x29, x30, [sp], #48 + ret + +.Lst_trivial: + mov x0, xzr + cbz x1, .Lst_trivial_ret + ldr w0, [x0] +.Lst_trivial_ret: + ret + + /* ---- local helper: siftdown(root = x0, end = x1) + * clobbers x9-x15 only; x19 (base) is live across the call. ---- */ +.Lst_siftdown: + mov x11, x0 /* root */ +.Lst_sift_loop: + lsl x12, x11, #1 + add x12, x12, #1 /* child = 2*root + 1 */ + cmp x12, x1 + b.hs .Lst_sift_done /* no children */ + + /* pick the larger of the two children */ + add x13, x12, #1 /* child + 1 */ + cmp x13, x1 + b.hs .Lst_sift_have_child + ldr w14, [x19, x12, lsl #2] + ldr w15, [x19, x13, lsl #2] + cmp w15, w14 + csel x12, x13, x12, hi + +.Lst_sift_have_child: + ldr w14, [x19, x11, lsl #2] /* a[root] */ + ldr w15, [x19, x12, lsl #2] /* a[child] */ + cmp w14, w15 + b.hs .Lst_sift_done /* heap property holds */ + + /* swap and descend */ + str w15, [x19, x11, lsl #2] + str w14, [x19, x12, lsl #2] + mov x11, x12 + b .Lst_sift_loop + +.Lst_sift_done: + ret +FN_END(fm_sort) + + +/* =================================================================== + * uint64_t fm_chase(void **ptrs, uint64_t steps) + * + * Pointer chase around a randomised cycle. Every load depends on the previous + * one, so nothing can be prefetched, overlapped or reordered - this measures + * the pure serial latency of the memory hierarchy, which is the single + * hardest thing for a wide out-of-order core to hide. It is the truest + * "single-threaded" test in the suite. + * =================================================================== */ +FN_BEGIN(fm_chase) + cbz x1, .Lch_zero + mov x2, x0 /* p = ptrs */ + mov x3, x1 + +.Lch_loop: + ldr x2, [x2] + subs x3, x3, #1 + b.ne .Lch_loop + + sub x0, x2, x0 /* final offset, keeps p live */ + ret + +.Lch_zero: + mov x0, xzr + ret +FN_END(fm_chase) + + +#if defined(__ELF__) + .section .note.GNU-stack, "", %progbits +#endif diff --git a/src/fossmark_x86_64.S b/src/fossmark_x86_64.S new file mode 100644 index 0000000..4332ecf --- /dev/null +++ b/src/fossmark_x86_64.S @@ -0,0 +1,987 @@ +/* + * fossmark_x86_64.S - x86-64 (AMD64) CPU benchmark kernels + * + * The AMD64 counterpart to fossmark.S. Same nine routines, same contract: each + * is a pure function of its arguments under the System V AMD64 ABI, contains no + * syscalls, no libc calls and no external data relocations, so it assembles and + * runs on Linux (ELF), macOS (Mach-O) and the BSDs. The portable C driver in + * main.c is shared unchanged between this file and the AArch64 one. + * + * Only baseline instructions are used: general-purpose AMD64 plus SSE2, which + * is architecturally mandatory on x86-64. The optional extensions (SSE4, AVX, + * FMA, AES-NI, POPCNT/BMI) are deliberately avoided - using them unguarded + * would fault (#UD) on cores that lack them - so, exactly as the NEON file + * sticks to mandatory ASIMD and shuns the optional AES/SHA/DotProd, this file + * sticks to mandatory SSE2 and shuns everything above it. FMA in particular is + * optional here, so every fused multiply-add is written as a separate multiply + * and add. + * + * Register conventions (System V AMD64): + * integer args rdi, rsi, rdx, rcx, r8, r9 (return in rax) + * callee-saved rbx, rbp, r12, r13, r14, r15 (saved when used) + * all of xmm0-15 are caller-saved, so no vector register need be preserved. + */ + + .intel_syntax noprefix + +#if defined(__APPLE__) +# define SYM(name) _##name +#else +# define SYM(name) name +#endif + +#if defined(__ELF__) +# define FN_BEGIN(name) .p2align 4 ; .globl SYM(name) ; .type SYM(name), @function ; SYM(name): +# define FN_END(name) .size SYM(name), . - SYM(name) +#else +# define FN_BEGIN(name) .p2align 4 ; .globl SYM(name) ; SYM(name): +# define FN_END(name) +#endif + + .text + +/* =================================================================== + * uint64_t fm_int_math(uint64_t iters) [rdi = iters] + * + * Four independent multiply-accumulate chains for instruction-level + * parallelism, mixed with the long-latency serialising ops (mul/div) and + * bit-manipulation. Returns a checksum so nothing can be elided. + * =================================================================== */ +FN_BEGIN(fm_int_math) + test rdi, rdi + jz .Lim_zero + + mov r8, 0x9E3779B97F4A7C15 /* a */ + mov r9, 0xBF58476D1CE4E5B9 /* b */ + mov r10, 0x94D049BB133111EB /* c */ + mov r11, 0x2545F4914F6CDD1D /* d */ + mov rsi, 0x00000000DEADBEEF /* odd multiplier, never zero */ + +.Lim_loop: + /* four independent multiply-accumulate chains */ + imul r8, rsi + add r8, r9 + imul r9, rsi + add r9, r10 + imul r10, rsi + add r10, r11 + imul r11, rsi + add r11, r8 + + /* cross-mix with shifts and logic ops */ + mov rax, r10 + shr rax, 29 + xor r8, rax + mov rax, r11 + shl rax, 17 + xor r9, rax + mov rax, r8 + ror rax, 31 + xor r10, rax + mov rax, r9 + sar rax, 7 + not rax + and r11, rax /* r11 &= ~(r9 >> 7 arith) */ + + /* wide multiplies: the long-latency 128-bit multiplier path */ + mov rax, r8 + mul r10 /* rdx:rax = a*c, high in rdx */ + add r8, rdx + mov rax, r9 + imul r11 /* rdx:rax = b*d signed, high rdx */ + add r9, rdx + + /* bit manipulation: byte reverse */ + mov rax, r10 + bswap rax + xor r8, rax + mov rax, r11 + bswap rax + xor r9, rax + + /* division: fully serialising, not pipelined */ + mov rcx, rsi + or rcx, 1 /* guarantee a non-zero divisor */ + mov rax, r8 + xor edx, edx + div rcx /* rax = a / rcx (unsigned) */ + imul rax, rcx + sub r10, rax /* c -= (a/div)*div */ + mov rax, r9 + cqo + idiv rcx /* rax = b / rcx (signed) */ + add r11, rax + + /* bitfield ops */ + ror r9, 11 + shld r8, r9, 23 + + dec rdi + jnz .Lim_loop + + mov rax, r8 + xor rax, r9 + xor rax, r10 + xor rax, r11 + ret + +.Lim_zero: + xor eax, eax + ret +FN_END(fm_int_math) + + +/* =================================================================== + * uint64_t fm_fp_math(uint64_t iters) [rdi = iters] + * + * Double-precision scalar FP. Four multiply-add chains for the pipelined + * path; sqrtsd and divsd for the non-pipelined divide/sqrt unit that usually + * separates cores. Returns the result bit-cast to u64. + * =================================================================== */ +FN_BEGIN(fm_fp_math) + test rdi, rdi + jz .Lfp_zero + + movsd xmm0, [rip + .Lfp_consts + 0] /* a = 1.5 */ + movsd xmm1, [rip + .Lfp_consts + 8] /* b = 2.5 */ + movsd xmm2, [rip + .Lfp_consts + 16] /* c = 3.5 */ + movsd xmm3, [rip + .Lfp_consts + 24] /* d = 0.5 */ + movsd xmm4, [rip + .Lfp_consts + 32] /* mul */ + movsd xmm5, [rip + .Lfp_consts + 40] /* addend */ + movsd xmm6, [rip + .Lfp_consts + 48] /* 2.0 */ + movsd xmm7, [rip + .Lfp_consts + 56] /* 1.0 */ + +.Lfp_loop: + /* four independent multiply-add chains (no baseline FMA) */ + mulsd xmm0, xmm4 + addsd xmm0, xmm5 + mulsd xmm1, xmm4 + addsd xmm1, xmm5 + mulsd xmm2, xmm4 + addsd xmm2, xmm5 + mulsd xmm3, xmm4 + addsd xmm3, xmm5 + + /* keep the accumulators bounded so they never reach inf/NaN */ + minsd xmm0, xmm6 + minsd xmm1, xmm6 + minsd xmm2, xmm6 + minsd xmm3, xmm6 + + /* square root: long latency, low throughput */ + sqrtsd xmm8, xmm0 + sqrtsd xmm9, xmm1 + addsd xmm2, xmm8 + addsd xmm3, xmm9 + + /* divide: 1/(c+1), divisor >= 1 so never zero */ + movapd xmm10, xmm2 + addsd xmm10, xmm7 + movapd xmm11, xmm7 + divsd xmm11, xmm10 + addsd xmm0, xmm11 + + movapd xmm10, xmm3 + addsd xmm10, xmm7 + movapd xmm11, xmm7 + divsd xmm11, xmm10 + addsd xmm1, xmm11 + + /* abs/neg/max: cheap ops to balance the mix */ + andpd xmm2, [rip + .Lfp_absmask] /* fabs(c) */ + xorpd xmm3, [rip + .Lfp_signmask] /* fneg(d) */ + andpd xmm3, [rip + .Lfp_absmask] /* fabs() */ + maxsd xmm3, xmm7 + + dec rdi + jnz .Lfp_loop + + addsd xmm0, xmm1 + addsd xmm2, xmm3 + addsd xmm0, xmm2 + movq rax, xmm0 + ret + +.Lfp_zero: + xor eax, eax + ret +FN_END(fm_fp_math) + + .p2align 4 +.Lfp_consts: + .double 1.5, 2.5 + .double 3.5, 0.5 + .double 1.0625, 0.0009765625 + .double 2.0, 1.0 + .p2align 4 +.Lfp_absmask: + .quad 0x7fffffffffffffff, 0x7fffffffffffffff +.Lfp_signmask: + .quad 0x8000000000000000, 0x8000000000000000 + + +/* =================================================================== + * uint64_t fm_primes(uint64_t limit, uint8_t *sieve) [rdi, rsi] + * + * Sieve of Eratosthenes over [0, limit). The routine clears the caller's + * scratch itself, so the clearing pass counts as measured work. Strided stores + * over a buffer larger than L1 make this a memory-hierarchy test too. Returns + * the prime count. + * =================================================================== */ +FN_BEGIN(fm_primes) + cmp rdi, 2 + jb .Lpr_none + + /* zero the sieve, 32 bytes per iteration */ + pxor xmm0, xmm0 + xor rax, rax /* index */ + mov rcx, rdi + and rcx, -32 /* bulk portion */ +.Lpr_clear32: + cmp rax, rcx + jae .Lpr_clear1 + movdqu [rsi + rax], xmm0 + movdqu [rsi + rax + 16], xmm0 + add rax, 32 + jmp .Lpr_clear32 +.Lpr_clear1: + cmp rax, rdi + jae .Lpr_clear_done + mov byte ptr [rsi + rax], 0 + inc rax + jmp .Lpr_clear1 +.Lpr_clear_done: + + /* mark 0 and 1 as composite */ + mov byte ptr [rsi], 1 + mov byte ptr [rsi + 1], 1 + + /* outer loop: i = 2; i*i < limit; i++ */ + mov r8, 2 +.Lpr_outer: + mov rax, r8 + imul rax, r8 /* i*i */ + cmp rax, rdi + jae .Lpr_count + + movzx edx, byte ptr [rsi + r8] + test dl, dl + jnz .Lpr_outer_next /* already composite, skip */ + + /* inner loop: mark multiples starting at i*i, stride i */ + mov r9, rax /* j = i*i */ +.Lpr_inner: + cmp r9, rdi + jae .Lpr_outer_next + mov byte ptr [rsi + r9], 1 + add r9, r8 + jmp .Lpr_inner + +.Lpr_outer_next: + inc r8 + jmp .Lpr_outer + + /* count the survivors */ +.Lpr_count: + xor eax, eax /* count */ + mov r9, 2 +.Lpr_count_loop: + cmp r9, rdi + jae .Lpr_done + cmp byte ptr [rsi + r9], 0 + jne .Lpr_count_next + inc rax +.Lpr_count_next: + inc r9 + jmp .Lpr_count_loop + +.Lpr_done: + ret + +.Lpr_none: + xor eax, eax + ret +FN_END(fm_primes) + + +/* =================================================================== + * uint64_t fm_simd(uint64_t iters, void *buf) [rdi = iters, rsi = buf] + * + * "Extended instructions": the SSE2 unit, which is architecturally mandatory + * on x86-64 and therefore safe without runtime feature detection. Packed + * 16-bit integer arithmetic is used (SSE2's widest integer multiply is 16-bit; + * 32-bit packed multiply, pmulld, is an SSE4.1 extension and is avoided), plus + * saturating/averaging ops, widening multiply-add, shuffles and the packed + * single-precision float path including the reciprocal/rsqrt estimates. + * + * buf must be at least 128 bytes. Returns a checksum. + * =================================================================== */ +FN_BEGIN(fm_simd) + test rdi, rdi + jz .Lsd_zero + + /* seed eight vectors from the scratch buffer */ + movdqu xmm0, [rsi] + movdqu xmm1, [rsi + 16] + movdqu xmm2, [rsi + 32] + movdqu xmm3, [rsi + 48] + movdqu xmm4, [rsi + 64] + movdqu xmm5, [rsi + 80] + movdqu xmm6, [rsi + 96] + movdqu xmm7, [rsi + 112] + + /* float operands: convert the integer seeds to finite floats rather + * than reinterpreting bit patterns (which would be NaNs/denormals) */ + cvtdq2ps xmm12, xmm0 + cvtdq2ps xmm13, xmm1 + pxor xmm14, xmm14 /* float accumulator */ + +.Lsd_loop: + /* 16-bit integer multiply-accumulate across independent vectors */ + pmullw xmm0, xmm1 + paddw xmm0, xmm2 + pmullw xmm1, xmm2 + paddw xmm1, xmm3 + pmullw xmm2, xmm3 + paddw xmm2, xmm0 + + /* saturating and averaging arithmetic */ + paddsw xmm4, xmm0 + paddusw xmm5, xmm1 + psubsw xmm6, xmm2 + psubusw xmm7, xmm3 + + /* widening multiply-add: 16->32 bit lanes */ + movdqa xmm8, xmm0 + pmaddwd xmm8, xmm1 + paddd xmm3, xmm8 + + /* high-half multiply */ + movdqa xmm9, xmm0 + pmulhw xmm9, xmm1 + pxor xmm2, xmm9 + + /* shifts and logic */ + movdqa xmm10, xmm2 + pslld xmm10, 3 + psrld xmm2, 29 + por xmm2, xmm10 + pand xmm3, xmm4 + + /* min/max */ + pmaxsw xmm5, xmm0 + pminsw xmm6, xmm1 + + /* byte average and sum-of-absolute-differences */ + pavgb xmm7, xmm4 + movdqa xmm11, xmm0 + psadbw xmm11, xmm5 + paddw xmm4, xmm11 + + /* lane shuffle: reverse the four 32-bit lanes */ + pshufd xmm0, xmm0, 0x1B + pxor xmm1, xmm0 + + /* single-precision float SIMD: multiply-add plus the reciprocal and + * rsqrt estimates that shader-style code leans on */ + movaps xmm15, xmm12 + mulps xmm15, xmm13 + addps xmm14, xmm15 + rcpps xmm8, xmm12 + rsqrtps xmm9, xmm14 + addps xmm14, xmm8 + mulps xmm12, xmm9 + + dec rdi + jnz .Lsd_loop + + /* fold the eight integer vectors together */ + pxor xmm0, xmm1 + pxor xmm2, xmm3 + pxor xmm4, xmm5 + pxor xmm6, xmm7 + pxor xmm0, xmm2 + pxor xmm4, xmm6 + pxor xmm0, xmm4 + + /* fold in the float accumulator (truncate to int lanes) */ + cvttps2dq xmm14, xmm14 + pxor xmm0, xmm14 + + /* horizontal add of the four 32-bit lanes -> single checksum */ + pshufd xmm1, xmm0, 0x4E + paddd xmm0, xmm1 + pshufd xmm1, xmm0, 0xB1 + paddd xmm0, xmm1 + movd eax, xmm0 + ret + +.Lsd_zero: + xor eax, eax + ret +FN_END(fm_simd) + + +/* =================================================================== + * uint64_t fm_compress(const uint8_t *src, uint64_t len, uint32_t *ht) + * [rdi, rsi, rdx] + * + * The match-finding inner loop of an LZ77 compressor (the LZ4 fast strategy): + * hash the next 4 bytes, probe a single-entry-per-bucket table, verify, then + * extend. Branch-heavy with a data-dependent, cache-missing table probe. + * + * ht must hold 1<<16 uint32_t (256 KiB); this routine clears it itself. + * Returns the encoded size in bytes. + * =================================================================== */ +FN_BEGIN(fm_compress) + push rbp + push rbx + push r12 + push r13 + push r14 + push r15 + + mov r12, rdi /* src */ + mov r13, rdx /* ht */ + + /* clear the hash table: 1<<16 entries * 4 bytes = 262144 bytes */ + pxor xmm0, xmm0 + xor rax, rax + mov ecx, 262144 +.Lcm_clear: + movdqu [r13 + rax], xmm0 + movdqu [r13 + rax + 16], xmm0 + movdqu [r13 + rax + 32], xmm0 + movdqu [r13 + rax + 48], xmm0 + add rax, 64 + cmp rax, rcx + jb .Lcm_clear + + cmp rsi, 16 + jb .Lcm_tiny + + mov r14, r12 /* ip */ + mov r15, r12 /* anchor */ + lea rbx, [r12 + rsi] /* end */ + lea r10, [rbx - 12] /* mflimit = end - 12 */ + xor ebp, ebp /* outsize */ + +.Lcm_loop: + cmp r14, r10 + jae .Lcm_flush + + mov eax, [r14] /* seq = load32(ip) */ + imul eax, eax, 0x9E3779B1 /* * Knuth prime 2654435761 */ + shr eax, 16 /* h = (seq*prime) >> 16 */ + + mov ecx, [r13 + rax*4] /* ref_off = ht[h] */ + mov rdx, r14 + sub rdx, r12 /* cur_off = ip - src */ + mov [r13 + rax*4], edx /* ht[h] = cur_off */ + + lea rsi, [r12 + rcx] /* ref = src + ref_off */ + cmp rsi, r14 + jae .Lcm_no_match /* ref must be strictly behind ip */ + + mov rax, r14 + sub rax, rsi /* distance */ + cmp rax, 65536 + jae .Lcm_no_match /* 16-bit offset window */ + + mov eax, [rsi] + cmp eax, [r14] /* verify the 4-byte match */ + jne .Lcm_no_match + + /* match confirmed: extend it byte by byte */ + mov r9, 4 /* ml */ +.Lcm_extend: + lea rax, [r14 + r9] + cmp rax, rbx /* ip + ml vs end */ + jae .Lcm_emit + mov cl, [r14 + r9] + cmp cl, [rsi + r9] + jne .Lcm_emit + inc r9 + jmp .Lcm_extend + +.Lcm_emit: + /* token(1) + offset(2) + literals + varint extensions */ + mov rax, r14 + sub rax, r15 /* literal run length */ + add rbp, rax + add rbp, 3 + cmp rax, 15 + jb .Lcm_no_lit_ext + inc rbp /* literal-length extension byte */ +.Lcm_no_lit_ext: + cmp r9, 19 + jb .Lcm_no_ml_ext + inc rbp /* match-length extension byte */ +.Lcm_no_ml_ext: + add r14, r9 /* ip += ml */ + mov r15, r14 /* anchor = ip */ + jmp .Lcm_loop + +.Lcm_no_match: + inc r14 + jmp .Lcm_loop + +.Lcm_flush: + /* trailing literals */ + mov rax, rbx + sub rax, r15 /* end - anchor */ + add rbp, rax + add rbp, 1 + mov rax, rbp + jmp .Lcm_ret + +.Lcm_tiny: + lea rax, [rsi + 1] /* len + 1 */ + +.Lcm_ret: + pop r15 + pop r14 + pop r13 + pop r12 + pop rbx + pop rbp + ret +FN_END(fm_compress) + + +/* =================================================================== + * uint64_t fm_chacha20(uint8_t *buf, uint64_t len, + * const uint8_t key[32], uint64_t rounds) + * [rdi, rsi, rdx, rcx] + * + * ChaCha20 stream cipher, SSE2, four 128-bit state rows. Chosen over AES for + * the same reason the NEON file chose it: the AES-NI extension is optional, so + * an AES-instruction benchmark would fault (#UD) on cores without it. ChaCha20 + * needs only baseline SSE2 and is a real, widely deployed cipher. + * + * len is rounded down to a multiple of 64. `rounds` = passes over the buffer. + * Returns a checksum of the keystream output. + * =================================================================== */ + +/* rotate each 32-bit lane left by n, via shift-left + shift-right + or */ +#define ROL32(v, n) \ + movdqa xmm14, v ;\ + pslld v, n ;\ + psrld xmm14, (32 - (n)) ;\ + por v, xmm14 + +/* one ChaCha quarter-round over rows a,b,c,d (xmm14 is scratch, via ROL32) */ +#define QROUND(a, b, c, d) \ + paddd a, b ;\ + pxor d, a ;\ + ROL32(d, 16) ;\ + paddd c, d ;\ + pxor b, c ;\ + ROL32(b, 12) ;\ + paddd a, b ;\ + pxor d, a ;\ + ROL32(d, 8) ;\ + paddd c, d ;\ + pxor b, c ;\ + ROL32(b, 7) + +FN_BEGIN(fm_chacha20) + and rsi, -64 /* whole 64-byte blocks only */ + jz .Lcc_zero + test rcx, rcx + jz .Lcc_zero + + /* xmm4..7 hold the base state */ + movdqa xmm4, [rip + .Lcc_sigma] /* "expand 32-byte k" */ + movdqu xmm5, [rdx] /* key[0..15] */ + movdqu xmm6, [rdx + 16] /* key[16..31] */ + pxor xmm7, xmm7 /* counter || nonce = 0 */ + + pxor xmm13, xmm13 /* running checksum */ + xor r8, r8 /* block counter value */ + +.Lcc_pass: + xor r9, r9 /* byte offset into buf */ + +.Lcc_block: + /* working state = base state, block counter in lane 0 of row 3. + * Row 3 is all-zero (nonce and counter), so a plain movd both sets + * the counter lane and clears the nonce lanes. */ + movdqa xmm0, xmm4 + movdqa xmm1, xmm5 + movdqa xmm2, xmm6 + movd xmm3, r8d + + /* keep originals for the final feed-forward add */ + movdqa xmm8, xmm0 + movdqa xmm9, xmm1 + movdqa xmm10, xmm2 + movdqa xmm11, xmm3 + + mov r10d, 10 /* 10 double rounds = 20 rounds */ +.Lcc_rounds: + /* column round */ + QROUND(xmm0, xmm1, xmm2, xmm3) + + /* rotate lanes to form the diagonals */ + pshufd xmm1, xmm1, 0x39 /* <<< 1 lane */ + pshufd xmm2, xmm2, 0x4E /* <<< 2 lanes */ + pshufd xmm3, xmm3, 0x93 /* <<< 3 lanes */ + + /* diagonal round */ + QROUND(xmm0, xmm1, xmm2, xmm3) + + /* undo the lane rotation */ + pshufd xmm1, xmm1, 0x93 + pshufd xmm2, xmm2, 0x4E + pshufd xmm3, xmm3, 0x39 + + dec r10d + jnz .Lcc_rounds + + /* feed-forward: keystream = working + original */ + paddd xmm0, xmm8 + paddd xmm1, xmm9 + paddd xmm2, xmm10 + paddd xmm3, xmm11 + + /* XOR the keystream into the buffer */ + lea rax, [rdi + r9] + movdqu xmm12, [rax] + pxor xmm12, xmm0 + movdqu [rax], xmm12 + movdqu xmm12, [rax + 16] + pxor xmm12, xmm1 + movdqu [rax + 16], xmm12 + movdqu xmm12, [rax + 32] + pxor xmm12, xmm2 + movdqu [rax + 32], xmm12 + movdqu xmm12, [rax + 48] + pxor xmm12, xmm3 + movdqu [rax + 48], xmm12 + + /* accumulate a checksum of the keystream */ + pxor xmm13, xmm0 + pxor xmm13, xmm3 + + inc r8 /* counter++ */ + add r9, 64 + cmp r9, rsi + jb .Lcc_block + + dec rcx + jnz .Lcc_pass + + /* horizontal add of the checksum lanes */ + pshufd xmm0, xmm13, 0x4E + paddd xmm13, xmm0 + pshufd xmm0, xmm13, 0xB1 + paddd xmm13, xmm0 + movd eax, xmm13 + ret + +.Lcc_zero: + xor eax, eax + ret +FN_END(fm_chacha20) + + .p2align 4 +.Lcc_sigma: + .long 0x61707865, 0x3320646e, 0x79622d32, 0x6b206574 + + +/* =================================================================== + * uint64_t fm_physics(double *bodies, uint64_t n, uint64_t steps) + * [rdi, rsi, rdx] + * + * Direct-summation N-body gravity, O(n^2) per step, double precision. + * Layout per body, 8 doubles (64 bytes): [x y z mass vx vy vz pad]. + * The 1/sqrt is a real sqrtsd+divsd (not the rsqrt estimate), exercising the + * divide/sqrt unit the way physics code does. Returns a velocity checksum. + * =================================================================== */ +FN_BEGIN(fm_physics) + test rsi, rsi + jz .Lph_zero + test rdx, rdx + jz .Lph_zero + + movsd xmm13, [rip + .Lph_dt] /* dt */ + movsd xmm14, [rip + .Lph_eps2] /* eps^2 */ + movsd xmm15, [rip + .Lph_one] /* 1.0 */ + +.Lph_step: + xor r8, r8 /* i */ + +.Lph_body_i: + mov rax, r8 + shl rax, 6 /* i * 64 */ + lea r9, [rdi + rax] /* &bodies[i] */ + + movsd xmm0, [r9] /* xi */ + movsd xmm1, [r9 + 8] /* yi */ + movsd xmm2, [r9 + 16] /* zi */ + + xorpd xmm3, xmm3 /* ax */ + xorpd xmm4, xmm4 /* ay */ + xorpd xmm5, xmm5 /* az */ + + xor r10, r10 /* j */ + mov r11, rdi /* &bodies[j] */ + +.Lph_body_j: + movsd xmm6, [r11] /* xj */ + movsd xmm7, [r11 + 8] /* yj */ + movsd xmm8, [r11 + 16] /* zj */ + movsd xmm9, [r11 + 24] /* mj */ + + subsd xmm6, xmm0 /* dx */ + subsd xmm7, xmm1 /* dy */ + subsd xmm8, xmm2 /* dz */ + + /* r2 = dx*dx + dy*dy + dz*dz + eps^2 (>= eps^2, so the i==j self-term + * is finite and contributes exactly 0 below) */ + movsd xmm10, xmm6 + mulsd xmm10, xmm6 + movsd xmm11, xmm7 + mulsd xmm11, xmm7 + addsd xmm10, xmm11 + movsd xmm11, xmm8 + mulsd xmm11, xmm8 + addsd xmm10, xmm11 + addsd xmm10, xmm14 /* + eps^2 */ + + sqrtsd xmm10, xmm10 /* r */ + movsd xmm11, xmm15 + divsd xmm11, xmm10 /* 1/r */ + movsd xmm12, xmm11 + mulsd xmm12, xmm11 /* 1/r^2 */ + mulsd xmm12, xmm11 /* 1/r^3 */ + mulsd xmm12, xmm9 /* m/r^3 */ + + mulsd xmm6, xmm12 /* dx * m/r^3 */ + addsd xmm3, xmm6 + mulsd xmm7, xmm12 + addsd xmm4, xmm7 + mulsd xmm8, xmm12 + addsd xmm5, xmm8 + + add r11, 64 + inc r10 + cmp r10, rsi + jb .Lph_body_j + + /* v += a * dt */ + movsd xmm6, [r9 + 32] /* vx */ + movsd xmm7, [r9 + 40] /* vy */ + movsd xmm8, [r9 + 48] /* vz */ + mulsd xmm3, xmm13 + addsd xmm6, xmm3 + mulsd xmm4, xmm13 + addsd xmm7, xmm4 + mulsd xmm5, xmm13 + addsd xmm8, xmm5 + movsd [r9 + 32], xmm6 + movsd [r9 + 40], xmm7 + movsd [r9 + 48], xmm8 + + inc r8 + cmp r8, rsi + jb .Lph_body_i + + /* second pass: x += v * dt (positions move only after all forces) */ + xor r8, r8 + mov r11, rdi +.Lph_integrate: + movsd xmm0, [r11] + movsd xmm1, [r11 + 8] + movsd xmm2, [r11 + 16] + movsd xmm6, [r11 + 32] + movsd xmm7, [r11 + 40] + movsd xmm8, [r11 + 48] + mulsd xmm6, xmm13 + addsd xmm0, xmm6 + mulsd xmm7, xmm13 + addsd xmm1, xmm7 + mulsd xmm8, xmm13 + addsd xmm2, xmm8 + movsd [r11], xmm0 + movsd [r11 + 8], xmm1 + movsd [r11 + 16], xmm2 + add r11, 64 + inc r8 + cmp r8, rsi + jb .Lph_integrate + + dec rdx + jnz .Lph_step + + /* checksum: sum of all velocity components */ + xorpd xmm0, xmm0 + xor r8, r8 + mov r11, rdi +.Lph_sum: + movsd xmm6, [r11 + 32] + movsd xmm7, [r11 + 40] + movsd xmm8, [r11 + 48] + addsd xmm0, xmm6 + addsd xmm0, xmm7 + addsd xmm0, xmm8 + add r11, 64 + inc r8 + cmp r8, rsi + jb .Lph_sum + + movq rax, xmm0 + ret + +.Lph_zero: + xor eax, eax + ret +FN_END(fm_physics) + + .p2align 4 +.Lph_dt: + .double 0.0078125 /* dt */ +.Lph_eps2: + .double 0.0625 /* eps^2 */ +.Lph_one: + .double 1.0 + + +/* =================================================================== + * uint64_t fm_sort(uint32_t *a, uint64_t n) [rdi = a, rsi = n] + * + * In-place heapsort: no recursion or explicit stack, aggressively + * branch-unpredictable, with scattered memory access - it stresses the branch + * predictor and the cache hierarchy. Returns an order-sensitive checksum, + * which also verifies the sort. + * + * Uses only caller-saved registers, so no prologue is needed; the internal + * siftdown is reached with `call` (contract below). + * =================================================================== */ +FN_BEGIN(fm_sort) + cmp rsi, 2 + jb .Lst_trivial + + /* ---- build the max-heap: for i = n/2 - 1 down to 0 ---- */ + mov r8, rsi + shr r8, 1 /* i = n/2 */ +.Lst_build: + test r8, r8 + jz .Lst_extract + dec r8 /* i-- */ + mov rcx, r8 /* root */ + mov rdx, rsi /* end */ + call .Lst_siftdown + test r8, r8 + jnz .Lst_build + + /* ---- extract: for end = n-1 down to 1 ---- */ +.Lst_extract: + mov r9, rsi + dec r9 /* end = n-1 */ +.Lst_extract_loop: + test r9, r9 + jz .Lst_checksum + + /* swap a[0] and a[end] */ + mov eax, [rdi] + mov r10d, [rdi + r9*4] + mov [rdi], r10d + mov [rdi + r9*4], eax + + xor ecx, ecx /* root = 0 */ + mov rdx, r9 /* end = end */ + call .Lst_siftdown + + dec r9 + jmp .Lst_extract_loop + + /* ---- order-sensitive checksum ---- */ +.Lst_checksum: + xor eax, eax + xor rcx, rcx /* index */ +.Lst_cksum_loop: + mov r10d, [rdi + rcx*4] + xor rax, r10 + ror rax, 7 + add rax, r10 + inc rcx + cmp rcx, rsi + jb .Lst_cksum_loop + ret + +.Lst_trivial: + xor eax, eax + test rsi, rsi + jz .Lst_trivial_ret + mov eax, [rdi] +.Lst_trivial_ret: + ret + + /* ---- local helper: siftdown(root = rcx, end = rdx) + * base a = rdi; clobbers rax, rcx, r10, r11 only. r8 (build i), + * r9 (extract end), rsi (n), rdx (end) all survive. ---- */ +.Lst_siftdown: + mov r11, rcx /* root */ +.Lst_sift_loop: + lea r10, [r11 + r11 + 1] /* child = 2*root + 1 */ + cmp r10, rdx + jae .Lst_sift_done /* no children */ + + lea rax, [r10 + 1] /* right = child + 1 */ + cmp rax, rdx + jae .Lst_sift_have_child /* no right child */ + mov ecx, [rdi + rax*4] /* a[right] */ + cmp ecx, [rdi + r10*4] /* a[right] vs a[child] */ + jbe .Lst_sift_have_child /* keep child if a[right] <= it */ + mov r10, rax /* else child = right */ + +.Lst_sift_have_child: + mov eax, [rdi + r11*4] /* a[root] */ + mov ecx, [rdi + r10*4] /* a[child] */ + cmp eax, ecx + jae .Lst_sift_done /* heap property holds */ + + /* swap and descend */ + mov [rdi + r11*4], ecx + mov [rdi + r10*4], eax + mov r11, r10 + jmp .Lst_sift_loop + +.Lst_sift_done: + ret +FN_END(fm_sort) + + +/* =================================================================== + * uint64_t fm_chase(void **ptrs, uint64_t steps) [rdi = ptrs, rsi = steps] + * + * Pointer chase around a randomised cycle. Every load depends on the previous + * one, so nothing can be prefetched, overlapped or reordered - this measures + * the pure serial latency of the memory hierarchy. The truest single-threaded + * test in the suite. + * =================================================================== */ +FN_BEGIN(fm_chase) + test rsi, rsi + jz .Lch_zero + mov rax, rdi /* p = ptrs */ + mov rcx, rsi + +.Lch_loop: + mov rax, [rax] + dec rcx + jnz .Lch_loop + + sub rax, rdi /* final offset, keeps p live */ + ret + +.Lch_zero: + xor eax, eax + ret +FN_END(fm_chase) + + +#if defined(__ELF__) + .section .note.GNU-stack, "", @progbits +#endif diff --git a/src/main.c b/src/main.c new file mode 100644 index 0000000..ee8124c --- /dev/null +++ b/src/main.c @@ -0,0 +1,709 @@ +/* + * fossmark - a multi-core AArch64 CPU benchmark + * + * This file is the portable driver: it owns everything the assembly kernels + * deliberately do not (timing, memory, I/O, scoring). The kernels in + * fossmark.S are pure computation and identical on every OS; only this file + * knows what an operating system is. + * + * Every workload is run twice: once on a single core, and once on all available + * cores at once - one identical copy of the kernel per core, each with its own + * private buffers, so the machine is driven to 100%% and the rate is whole-machine + * throughput. From these two passes fossmark reports two composite scores, a + * SINGLECORE and a MULTICORE, from the same tests and the same weights. + * + * Build: cc -O2 -pthread main.c fossmark.S -o fossmark -lm + */ + +#include +#include +#include +#include +#include +#include +#include + +/* ---------- platform identification (for the banner only) ---------- */ + +#if defined(_WIN32) +# define FM_OS "Windows" +#elif defined(__APPLE__) +# define FM_OS "macOS" +#elif defined(__linux__) +# define FM_OS "Linux" +#else +# define FM_OS "POSIX" +#endif + +#if defined(__aarch64__) || defined(_M_ARM64) +# define FM_ARCH "ARM64" +# define D_INT "64-bit ALU: madd, umulh, udiv, bitops" +# define D_FP "double: fmadd, fdiv, fsqrt" +# define D_SIMD "NEON ASIMD: 128-bit integer + float" +#elif defined(__x86_64__) || defined(_M_X64) +# define FM_ARCH "x86-64" +# define D_INT "64-bit ALU: imul, mul, div, bitops" +# define D_FP "double: mulsd/addsd, divsd, sqrtsd" +# define D_SIMD "SSE2: 128-bit integer + float" +#else +# define FM_ARCH "unknown" +# define D_INT "64-bit integer ALU" +# define D_FP "double-precision FP" +# define D_SIMD "128-bit SIMD: integer + float" +#endif + +/* ---------- portable monotonic clock ---------- */ + +#if defined(_WIN32) +# define WIN32_LEAN_AND_MEAN +# include +static double now_seconds(void) +{ + LARGE_INTEGER f, t; + QueryPerformanceFrequency(&f); + QueryPerformanceCounter(&t); + return (double)t.QuadPart / (double)f.QuadPart; +} +#else +# include +static double now_seconds(void) +{ + struct timespec ts; + clock_gettime(CLOCK_MONOTONIC, &ts); + return (double)ts.tv_sec + (double)ts.tv_nsec * 1e-9; +} +#endif + +/* ---------- the assembly kernels ---------- */ + +extern uint64_t fm_int_math(uint64_t iters); +extern uint64_t fm_fp_math(uint64_t iters); +extern uint64_t fm_primes(uint64_t limit, uint8_t *sieve); +extern uint64_t fm_simd(uint64_t iters, void *buf); +extern uint64_t fm_compress(const uint8_t *src, uint64_t len, uint32_t *ht); +extern uint64_t fm_chacha20(uint8_t *buf, uint64_t len, + const uint8_t key[32], uint64_t rounds); +extern uint64_t fm_physics(double *bodies, uint64_t n, uint64_t steps); +extern uint64_t fm_sort(uint32_t *a, uint64_t n); +extern uint64_t fm_chase(void **ptrs, uint64_t steps); + +/* ---------- tuning ---------- */ + +#define PRIME_LIMIT (2u * 1000u * 1000u) /* sieve span */ +#define COMPRESS_LEN (4u * 1024u * 1024u) /* corpus size */ +#define HT_ENTRIES (1u << 16) /* LZ77 hash buckets */ +#define CIPHER_LEN (1u * 1024u * 1024u) /* plaintext size */ +#define SIMD_BUF 256 /* NEON scratch */ +#define NBODY_N 512 /* bodies */ +#define SORT_N (1u << 20) /* elements to sort */ +#define CHASE_NODES (1u << 21) /* 16 MiB cycle, > any L2 */ + +#define MIN_SECONDS 2.0 /* per-test measured floor */ +#define REPEATS 3 /* best-of, to reject noise */ + +/* ---------- scoring configuration ---------- + * + * The overall score is a WEIGHTED geometric mean of each test's rate expressed + * relative to a reference machine. Two knobs per test: + * + * FM_REF_* the reference rate (this machine's measured rate). A machine + * matching the reference scores FM_TARGET_SCORE on that test. + * FM_WEIGHT_* how much that test counts toward the overall, by its + * influence on everyday user experience. Weights are relative: + * only their ratios matter, so they need not sum to anything - + * the code normalises by their sum. (They happen to sum to 100 + * here, so each reads as a percent.) + * + * Per-test score: S_i = FM_TARGET_SCORE * (rate_i / FM_REF_i) + * Overall score: Overall = FM_TARGET_SCORE * + * exp( Sum(w_i * ln(rate_i/FM_REF_i)) / Sum(w_i) ) + * + * On the reference machine every ratio is 1, so every S_i and the overall come + * out to exactly FM_TARGET_SCORE, regardless of the weights. Scaling is linear + * in performance, so far slower machines fall well below (half as fast -> half + * the score) and faster future machines rise above. + */ + +#define FM_TARGET_SCORE 10000.0 /* reference-machine overall */ + +/* Reference rates: this machine, in each test's native unit (see tests[]). */ +#define FM_REF_INT 3086.0 /* Mops/s */ +#define FM_REF_FP 1682.0 /* Mops/s */ +#define FM_REF_PRIMES 812.0 /* Mcand/s */ +#define FM_REF_SIMD 6576.0 /* Mops/s */ +#define FM_REF_COMPRESS 674.0 /* MB/s */ +#define FM_REF_CRYPTO 406.0 /* MB/s */ +#define FM_REF_PHYSICS 631.0 /* Mpair/s */ +#define FM_REF_SORT 363.0 /* Mkey-cmp/s*/ +#define FM_REF_CHASE 79.0 /* Mhop/s (scoring); shown as ns/access */ + +/* Weights: influence on day-to-day, common-workload user experience. + * Rationale: integer/general-purpose code and memory-latency-bound + * responsiveness dominate everyday use; specialised FP/physics matter least. + * Roughly an 80/20 integer-vs-FP split, in the spirit of Geekbench 6's + * weighted, integer-dominant methodology. Retune freely. */ +#define FM_WEIGHT_INT 20.0 /* general-purpose ALU: everything */ +#define FM_WEIGHT_CHASE 16.0 /* memory latency: responsiveness */ +#define FM_WEIGHT_COMPRESS 14.0 /* web, storage, RAM compression */ +#define FM_WEIGHT_SORT 12.0 /* general data-structure work */ +#define FM_WEIGHT_SIMD 11.0 /* codecs, mem/string ops, parsing */ +#define FM_WEIGHT_FP 9.0 /* spreadsheets, app/media math */ +#define FM_WEIGHT_CRYPTO 8.0 /* TLS, disk encryption (small frac) */ +#define FM_WEIGHT_PRIMES 6.0 /* synthetic ALU+memory proxy */ +#define FM_WEIGHT_PHYSICS 4.0 /* niche simulation/games */ + +/* ---------- deterministic PRNG (splitmix64) ---------- */ + +static uint64_t rng_state = 0x853c49e6748fea9bULL; + +static uint64_t rng_next(void) +{ + uint64_t z = (rng_state += 0x9e3779b97f4a7c15ULL); + z = (z ^ (z >> 30)) * 0xbf58476d1ce4e5b9ULL; + z = (z ^ (z >> 27)) * 0x94d049bb133111ebULL; + return z ^ (z >> 31); +} + +static void rng_reset(void) { rng_state = 0x853c49e6748fea9bULL; } + +/* ---------- aligned allocation ---------- */ + +static void *xalloc(size_t n) +{ + void *p = NULL; +#if defined(_WIN32) + p = _aligned_malloc(n, 64); +#else + if (posix_memalign(&p, 64, n) != 0) + p = NULL; +#endif + if (!p) { + fprintf(stderr, "fossmark: out of memory (%zu bytes)\n", n); + exit(1); + } + return p; +} + +static void xfree(void *p) +{ +#if defined(_WIN32) + _aligned_free(p); +#else + free(p); +#endif +} + +/* ---------- workload state ---------- + * + * Because every core runs the same kernel simultaneously, each core needs its + * OWN mutable buffers - sharing them would be a data race and would corrupt + * both the results and the determinism check. Per-core scratch lives in a + * `workspace`, one per thread. Read-only inputs (the corpus, the pristine + * physics/sort seeds, the key, the chase graph) are genuinely shared. + */ +struct workspace { + uint8_t *sieve; /* prime sieve scratch */ + uint32_t *ht; /* LZ77 hash table scratch */ + uint8_t *cipher_buf; /* ChaCha20 buffer, encrypted in place */ + uint8_t *simd_buf; /* NEON scratch */ + double *bodies; /* n-body integration buffer */ + uint32_t *sort_work; /* the buffer we actually sort */ + void **chase; /* private 16 MiB pointer-chase cycle */ +}; + +static long g_ncores = 1; /* active online cores */ +static struct workspace *g_ws; /* g_ncores per-thread workspaces */ + +static uint8_t *g_corpus; /* shared, read-only compression input */ +static uint8_t g_key[32]; /* shared, read-only cipher key */ +static uint8_t *g_cipher_src; /* pristine plaintext, copied per-core */ +static uint8_t *g_simd_src; /* pristine NEON seed, copied per-core */ +static double *g_bodies_src; /* pristine initial conditions */ +static uint32_t *g_sort_src; /* pristine unsorted data */ + +/* + * Synthesise a compressible corpus. Random bytes would be incompressible and + * would make the match-finder trivially miss every probe, measuring nothing + * interesting. This builds text-like data with realistic repetition instead. + */ +static void build_corpus(uint8_t *buf, size_t len) +{ + static const char *words[] = { + "the", "quick", "brown", "fox", "jumps", "over", "lazy", + "dog", "benchmark", "processor", "assembly", "vector", + "memory", "cache", "pipeline", "instruction", "compress", + "data", "system", "performance", "register", "kernel" + }; + const size_t nwords = sizeof(words) / sizeof(words[0]); + size_t pos = 0; + + while (pos < len) { + const char *w = words[rng_next() % nwords]; + size_t wl = strlen(w); + + if (pos + wl + 1 > len) + break; + memcpy(buf + pos, w, wl); + pos += wl; + buf[pos++] = (rng_next() % 8 == 0) ? '\n' : ' '; + } + while (pos < len) + buf[pos++] = ' '; +} + +/* Build a single random cycle through the node array (Sattolo's algorithm), + * guaranteeing one cycle of exactly CHASE_NODES steps with no early closure. */ +static void build_chase(void **nodes, size_t n) +{ + size_t *perm = xalloc(n * sizeof(size_t)); + size_t i; + + for (i = 0; i < n; i++) + perm[i] = i; + for (i = n - 1; i > 0; i--) { + size_t j = (size_t)(rng_next() % i); /* strictly j < i */ + size_t t = perm[i]; + perm[i] = perm[j]; + perm[j] = t; + } + for (i = 0; i < n; i++) + nodes[perm[i]] = (void *)&nodes[perm[(i + 1) % n]]; + + xfree(perm); +} + +static void setup(void) +{ + size_t i; + long t; + + rng_reset(); + + /* shared read-only inputs and pristine per-core seeds */ + g_corpus = xalloc(COMPRESS_LEN); + g_cipher_src = xalloc(CIPHER_LEN); + g_simd_src = xalloc(SIMD_BUF); + g_bodies_src = xalloc(NBODY_N * 8 * sizeof(double)); + g_sort_src = xalloc(SORT_N * sizeof(uint32_t)); + + build_corpus(g_corpus, COMPRESS_LEN); + + for (i = 0; i < CIPHER_LEN; i++) + g_cipher_src[i] = (uint8_t)rng_next(); + for (i = 0; i < 32; i++) + g_key[i] = (uint8_t)rng_next(); + for (i = 0; i < SIMD_BUF; i++) + g_simd_src[i] = (uint8_t)rng_next(); + for (i = 0; i < SORT_N; i++) + g_sort_src[i] = (uint32_t)rng_next(); + + /* bodies: [x y z mass vx vy vz pad], positions in a unit-ish cube */ + for (i = 0; i < NBODY_N; i++) { + double *b = &g_bodies_src[i * 8]; + b[0] = (double)(rng_next() % 2000) / 1000.0 - 1.0; + b[1] = (double)(rng_next() % 2000) / 1000.0 - 1.0; + b[2] = (double)(rng_next() % 2000) / 1000.0 - 1.0; + b[3] = (double)(rng_next() % 900) / 1000.0 + 0.1; /* mass > 0 */ + b[4] = b[5] = b[6] = 0.0; + b[7] = 0.0; + } + + /* one private workspace per core, every copy seeded identically so all + * cores compute the same deterministic result. Each core also gets its + * own pointer-chase cycle: sharing one would collapse the multi-core + * latency test into a shared-cache test instead of a memory test. */ + g_ws = xalloc((size_t)g_ncores * sizeof *g_ws); + for (t = 0; t < g_ncores; t++) { + struct workspace *w = &g_ws[t]; + + w->sieve = xalloc(PRIME_LIMIT); + w->ht = xalloc(HT_ENTRIES * sizeof(uint32_t)); + w->cipher_buf = xalloc(CIPHER_LEN); + w->simd_buf = xalloc(SIMD_BUF); + w->bodies = xalloc(NBODY_N * 8 * sizeof(double)); + w->sort_work = xalloc(SORT_N * sizeof(uint32_t)); + w->chase = xalloc(CHASE_NODES * sizeof(void *)); + + memcpy(w->cipher_buf, g_cipher_src, CIPHER_LEN); + memcpy(w->simd_buf, g_simd_src, SIMD_BUF); + build_chase(w->chase, CHASE_NODES); + } +} + +static void teardown(void) +{ + long t; + + for (t = 0; t < g_ncores; t++) { + struct workspace *w = &g_ws[t]; + + xfree(w->sieve); xfree(w->ht); xfree(w->cipher_buf); + xfree(w->simd_buf); xfree(w->bodies); xfree(w->sort_work); + xfree(w->chase); + } + xfree(g_ws); + + xfree(g_corpus); xfree(g_cipher_src); xfree(g_simd_src); + xfree(g_bodies_src); xfree(g_sort_src); +} + +/* ---------- the test harness ---------- */ + +/* + * Each test runs a kernel `n` times against a per-core workspace and returns a + * checksum. The harness auto-calibrates `n` upward until the run exceeds + * MIN_SECONDS, so the result is insensitive to clock granularity and to how + * fast the machine is. + */ +typedef uint64_t (*run_fn)(uint64_t n, struct workspace *ws); + +struct test { + const char *name; + const char *detail; + run_fn run; + uint64_t start_n; + double work_per_n; /* abstract work units, for scoring */ + const char *unit; + double ref_rate; /* reference-machine rate, in `unit` */ + double weight; /* relative weight in the overall score */ +}; + +static uint64_t run_int(uint64_t n, struct workspace *ws) +{ + (void)ws; + return fm_int_math(n * 100000); +} +static uint64_t run_fp(uint64_t n, struct workspace *ws) +{ + (void)ws; + return fm_fp_math(n * 100000); +} +static uint64_t run_primes(uint64_t n, struct workspace *ws) +{ + uint64_t c = 0; + for (uint64_t i = 0; i < n; i++) + c += fm_primes(PRIME_LIMIT, ws->sieve); + return c; +} +static uint64_t run_simd(uint64_t n, struct workspace *ws) +{ + return fm_simd(n * 100000, ws->simd_buf); +} +static uint64_t run_compress(uint64_t n, struct workspace *ws) +{ + uint64_t c = 0; + for (uint64_t i = 0; i < n; i++) + c += fm_compress(g_corpus, COMPRESS_LEN, ws->ht); + return c; +} +static uint64_t run_crypto(uint64_t n, struct workspace *ws) +{ + return fm_chacha20(ws->cipher_buf, CIPHER_LEN, g_key, n); +} +static uint64_t run_physics(uint64_t n, struct workspace *ws) +{ + /* restore initial conditions: the integrator mutates the bodies, so + * a re-run must start from the same state to be reproducible */ + memcpy(ws->bodies, g_bodies_src, NBODY_N * 8 * sizeof(double)); + return fm_physics(ws->bodies, NBODY_N, n); +} +static uint64_t run_sort(uint64_t n, struct workspace *ws) +{ + uint64_t c = 0; + for (uint64_t i = 0; i < n; i++) { + /* restore the pristine data: sorting an already-sorted array + * would measure the best case, not the real one */ + memcpy(ws->sort_work, g_sort_src, SORT_N * sizeof(uint32_t)); + c ^= fm_sort(ws->sort_work, SORT_N); + } + return c; +} +static uint64_t run_chase(uint64_t n, struct workspace *ws) +{ + return fm_chase(ws->chase, n * 1000000); +} + +static const struct test tests[] = { + { "Integer Math", D_INT, + run_int, 20, 100000.0 * 24, "Mops/s", + FM_REF_INT, FM_WEIGHT_INT }, + { "Floating Point Math", D_FP, + run_fp, 20, 100000.0 * 20, "Mops/s", + FM_REF_FP, FM_WEIGHT_FP }, + { "Prime Numbers", "sieve of Eratosthenes to 2M", + run_primes, 1, (double)PRIME_LIMIT, "Mcand/s", + FM_REF_PRIMES, FM_WEIGHT_PRIMES }, + { "Extended Instructions",D_SIMD, + run_simd, 10, 100000.0 * 32, "Mops/s", + FM_REF_SIMD, FM_WEIGHT_SIMD }, + { "Compression", "LZ77 match finder, 4 MiB corpus", + run_compress, 1, (double)COMPRESS_LEN, "MB/s", + FM_REF_COMPRESS, FM_WEIGHT_COMPRESS }, + { "Encryption", "ChaCha20, 20 rounds, 1 MiB", + run_crypto, 4, (double)CIPHER_LEN, "MB/s", + FM_REF_CRYPTO, FM_WEIGHT_CRYPTO }, + { "Physics", "512-body direct-sum gravity", + run_physics, 4, (double)NBODY_N * NBODY_N, "Mpair/s", + FM_REF_PHYSICS, FM_WEIGHT_PHYSICS }, + { "Sorting", "heapsort, 1M uint32", + run_sort, 1, (double)SORT_N * 20, "Mkey-cmp/s", + FM_REF_SORT, FM_WEIGHT_SORT }, + { "Memory Latency", "dependent-load pointer chase, 16 MiB", + run_chase, 1, 1000000.0, "ns/access", + FM_REF_CHASE, FM_WEIGHT_CHASE }, +}; + +#define NTESTS (sizeof(tests) / sizeof(tests[0])) + +struct result { + double rate; /* work units per second */ + double score; + uint64_t checksum; + double seconds; + uint64_t iters; + int threads; /* cores this test was spread across */ +}; + +/* + * One unit of parallel work: run `run(n, ws)` on a private workspace. Every + * core executes the identical kernel on identically-seeded data, so all cores + * return the same checksum; the harness sums them into one aggregate that stays + * deterministic across repeats. + */ +struct job { + run_fn run; + uint64_t n; + struct workspace *ws; + uint64_t result; +}; + +static void *job_entry(void *arg) +{ + struct job *j = arg; + j->result = j->run(j->n, j->ws); + return NULL; +} + +/* + * Run the kernel on `threads` cores at once and return the summed checksum. + * The calling thread runs job 0 itself; threads 1..N-1 run on spawned workers. + * A thread that fails to spawn simply runs inline, so the benchmark still + * completes (with less parallelism) rather than aborting. + */ +static uint64_t dispatch(run_fn run, uint64_t n, int threads) +{ + struct job *jobs = xalloc((size_t)threads * sizeof *jobs); + pthread_t *tids = threads > 1 + ? xalloc((size_t)(threads - 1) * sizeof *tids) : NULL; + int i, spawned = 0; + uint64_t agg = 0; + + for (i = 0; i < threads; i++) { + jobs[i].run = run; + jobs[i].n = n; + jobs[i].ws = &g_ws[i]; + } + for (i = 1; i < threads; i++) { + if (pthread_create(&tids[spawned], NULL, job_entry, &jobs[i]) == 0) + spawned++; + else + job_entry(&jobs[i]); /* fall back to inline */ + } + + job_entry(&jobs[0]); /* this thread runs job 0 */ + + for (i = 0; i < spawned; i++) + pthread_join(tids[i], NULL); + for (i = 0; i < threads; i++) + agg += jobs[i].result; + + xfree(jobs); + xfree(tids); + return agg; +} + +static struct result run_test(const struct test *t, int threads) +{ + struct result r; + uint64_t n = t->start_n; + double elapsed = 0.0, best = 0.0; + uint64_t checksum = 0; + int i; + + /* calibrate: grow n until a single run clears the noise floor */ + for (;;) { + double t0 = now_seconds(); + checksum = dispatch(t->run, n, threads); + elapsed = now_seconds() - t0; + + if (elapsed >= MIN_SECONDS) + break; + if (elapsed < 0.001) { + n *= 8; /* far too fast to measure */ + } else { + double scale = (MIN_SECONDS * 1.3) / elapsed; + if (scale < 1.5) + scale = 1.5; + if (scale > 8.0) + scale = 8.0; + n = (uint64_t)((double)n * scale) + 1; + } + } + + /* best-of: the fastest run is the one least disturbed by the OS */ + best = elapsed; + for (i = 1; i < REPEATS; i++) { + double t0 = now_seconds(); + uint64_t c = dispatch(t->run, n, threads); + double e = now_seconds() - t0; + + if (c != checksum) { + fprintf(stderr, + "fossmark: %s is non-deterministic " + "(checksum %llu != %llu)\n", t->name, + (unsigned long long)c, + (unsigned long long)checksum); + exit(2); + } + if (e < best) + best = e; + } + + r.seconds = best; + r.iters = n; + r.checksum = checksum; + r.threads = threads; + /* aggregate throughput: `threads` cores each did n*work_per_n of work in + * the same wall-clock window, so the machine's rate is their sum */ + r.rate = ((double)threads * (double)n * t->work_per_n) / best / 1e6; + /* normalise against the reference machine: this is the per-test score */ + r.score = FM_TARGET_SCORE * (r.rate / t->ref_rate); + return r; +} + +/* + * The number shown in the RATE column. Most tests report throughput in their + * `unit`. The memory-latency test is different: throughput (hops/s) is not what + * anyone reasons about for memory, so we report the actual per-access latency + * in nanoseconds instead. That is a per-core property -- the time for one + * dependent load in the chain -- so it is derived from a single core's hop + * count and is independent of how many cores ran, unlike the aggregate `rate`. + */ +static double display_metric(const struct test *t, const struct result *r) +{ + if (t->run == run_chase) { + double hops = (double)r->iters * t->work_per_n; /* per core */ + return r->seconds / hops * 1e9; /* ns/access */ + } + return r->rate; +} + +/* ---------- output ---------- */ + +static void print_header(void) +{ + printf("\n"); + printf(" fossmark 1.0 - multi-core CPU benchmark\n"); + printf(" ------------------------------------------------------------------\n"); + printf(" platform: %s/%s\n", FM_OS, FM_ARCH); + printf(" cores: %ld (each test is run once on 1 core, once on all %ld)\n", + g_ncores, g_ncores); + printf("\n"); + printf(" RATE/TIME/SCORE below are the all-core (multi-core) pass.\n"); + printf("\n"); + printf(" %-24s %12s %-11s %8s %9s\n", + "TEST", "RATE", "UNIT", "TIME", "SCORE"); + printf(" --------------------------------------------------------------------------\n"); + fflush(stdout); +} + +int main(int argc, char **argv) +{ + struct result multi[NTESTS], single[NTESTS]; + double multi_log_sum = 0.0, single_log_sum = 0.0; + double weight_sum = 0.0; + int verbose = 0; + size_t i; + + for (i = 1; i < (size_t)argc; i++) { + if (strcmp(argv[i], "-v") == 0 || + strcmp(argv[i], "--verbose") == 0) { + verbose = 1; + } else if (strcmp(argv[i], "-h") == 0 || + strcmp(argv[i], "--help") == 0) { + printf("usage: %s [-v|--verbose]\n", argv[0]); + return 0; + } else { + fprintf(stderr, "fossmark: unknown option '%s'\n", + argv[i]); + return 1; + } + } + + { + long n = sysconf(_SC_NPROCESSORS_ONLN); + g_ncores = n > 0 ? n : 1; + } + + printf("\n preparing workloads..."); + fflush(stdout); + setup(); + printf(" done\n"); + + print_header(); + + for (i = 0; i < NTESTS; i++) { + double sm, ss; + + printf(" %-24s", tests[i].name); + fflush(stdout); + + /* each test runs twice: the all-core pass (shown) and the + * single-core pass (folded into the SINGLECORE score) */ + multi[i] = run_test(&tests[i], (int)g_ncores); + single[i] = run_test(&tests[i], 1); + + printf(" %12.1f %-11s %7.2fs %9.0f\n", + display_metric(&tests[i], &multi[i]), tests[i].unit, + multi[i].seconds, multi[i].score); + if (verbose) + printf(" %-24s %s\n" + " %-24s weight=%.0f%% 1-core: %.1f %s / %.0f %ld-core: %.1f %s / %.0f\n", + "", tests[i].detail, "", tests[i].weight, + display_metric(&tests[i], &single[i]), tests[i].unit, + single[i].score, g_ncores, + display_metric(&tests[i], &multi[i]), tests[i].unit, + multi[i].score); + fflush(stdout); + + /* accumulate the weighted geometric mean of BOTH passes, same + * weights, so the two composite scores are directly comparable */ + sm = multi[i].score > 0.0 ? multi[i].score : 1e-9; + ss = single[i].score > 0.0 ? single[i].score : 1e-9; + multi_log_sum += tests[i].weight * log(sm); + single_log_sum += tests[i].weight * log(ss); + weight_sum += tests[i].weight; + } + + printf(" --------------------------------------------------------------------------\n"); + + /* + * Two composite scores, each the WEIGHTED geometric mean of the per-test + * scores from one pass. Per-test scores are already normalised so the + * single-thread reference machine reads FM_TARGET_SCORE. Geometric rather + * than arithmetic so no single test dominates; weighted so tests count in + * proportion to their influence on everyday use (the FM_WEIGHT_* config). + * The two passes share tests and weights, so MULTICORE / SINGLECORE is a + * clean read of how much the machine gains from all its cores. + */ + printf(" %-24s %44.0f\n", "MULTICORE SCORE", + exp(multi_log_sum / weight_sum)); + printf(" %-24s %44.0f\n", "SINGLECORE SCORE", + exp(single_log_sum / weight_sum)); + printf(" %-24s (weighted geometric means, reference machine = %.0f)\n", + "", FM_TARGET_SCORE); + printf("\n"); + + teardown(); + return 0; +} diff --git a/src/test_kernels.c b/src/test_kernels.c new file mode 100644 index 0000000..9ae3595 --- /dev/null +++ b/src/test_kernels.c @@ -0,0 +1,462 @@ +/* + * test_kernels.c - correctness checks for the fossmark assembly kernels + * + * The benchmark's own best-of-N run guards against non-determinism, but a + * kernel can be perfectly deterministic and still wrong. This file is the + * "single C file to poke at and test with": it validates each kernel against + * an independent reference or an invariant, so a mistake in the assembly is + * caught here rather than silently skewing a score. + * + * Every check (except the single-threaded pointer-chase) is run concurrently + * on all available cores. The kernels take their buffers as arguments and hold + * no shared state, so a correct kernel must give identical, correct results no + * matter how many copies run at once; a hidden global or a reentrancy bug would + * survive a single-threaded run but fail here. + * + * Build: cc -O2 -pthread test_kernels.c fossmark.S -o test_kernels -lm + * Exit status is 0 iff every check passes. + */ + +#include +#include +#include +#include +#include +#include +#include +#include + +extern uint64_t fm_int_math(uint64_t iters); +extern uint64_t fm_fp_math(uint64_t iters); +extern uint64_t fm_primes(uint64_t limit, uint8_t *sieve); +extern uint64_t fm_simd(uint64_t iters, void *buf); +extern uint64_t fm_compress(const uint8_t *src, uint64_t len, uint32_t *ht); +extern uint64_t fm_chacha20(uint8_t *buf, uint64_t len, + const uint8_t key[32], uint64_t rounds); +extern uint64_t fm_physics(double *bodies, uint64_t n, uint64_t steps); +extern uint64_t fm_sort(uint32_t *a, uint64_t n); +extern uint64_t fm_chase(void **ptrs, uint64_t steps); + +static int failures = 0; +static int checks = 0; + +/* + * Concurrency plumbing. Each check runs on every core at once; the counters and + * stdout are shared, so ok()/note() serialise on this lock. `fm_primary` is set + * on exactly one thread per check (the one running on the main thread): it owns + * the human-readable output so the "[ ok ]" lines and diagnostics appear once, + * not once per core. Every thread still evaluates every assertion, so a failure + * on any core - even a silent secondary - is reported and counted. + */ +static pthread_mutex_t io_lock = PTHREAD_MUTEX_INITIALIZER; +static __thread int fm_primary = 1; +static long fm_ncores = 1; + +static void ok(const char *what, int cond) +{ + pthread_mutex_lock(&io_lock); + if (fm_primary) { + checks++; + if (cond) { + printf(" [ ok ] %s\n", what); + } else { + printf(" [FAIL] %s\n", what); + failures++; + } + } else if (!cond) { + /* a secondary core disagrees: surface it explicitly */ + printf(" [FAIL] %s (concurrent core)\n", what); + failures++; + } + pthread_mutex_unlock(&io_lock); +} + +/* Diagnostic output that should appear once per check, not once per core. */ +static void note(const char *fmt, ...) +{ + va_list ap; + + if (!fm_primary) + return; + pthread_mutex_lock(&io_lock); + va_start(ap, fmt); + vprintf(fmt, ap); + va_end(ap); + pthread_mutex_unlock(&io_lock); +} + +/* Run `check` on every core simultaneously. The main thread is the primary; + * fm_ncores-1 workers run the same check as silent secondaries. */ +static void *fm_worker(void *arg) +{ + void (*check)(void) = *(void (**)(void))arg; + + fm_primary = 0; + check(); + return NULL; +} + +static void parallel(void (*check)(void)) +{ + long extra = fm_ncores - 1; + pthread_t *th = NULL; + long i, spawned = 0; + + if (extra > 0) { + th = calloc((size_t)extra, sizeof *th); + if (th) { + for (i = 0; i < extra; i++) + if (pthread_create(&th[spawned], NULL, + fm_worker, &check) == 0) + spawned++; + } + } + + check(); /* primary runs on this thread */ + + for (i = 0; i < spawned; i++) + pthread_join(th[i], NULL); + free(th); +} + +/* ---------- reference implementations ---------- */ + +static uint64_t ref_prime_count(uint64_t limit) +{ + uint8_t *s = calloc(limit, 1); + uint64_t count = 0, i, j; + + for (i = 2; i * i < limit; i++) + if (!s[i]) + for (j = i * i; j < limit; j += i) + s[j] = 1; + for (i = 2; i < limit; i++) + if (!s[i]) + count++; + free(s); + return count; +} + +/* A textbook scalar ChaCha20 block function, used both to anchor against the + * RFC 8439 known-answer vector and to validate the NEON kernel block-for-block. + * `out` receives 64 keystream bytes for the given counter and 12-byte nonce. */ +#define ROTL32(x, n) (((x) << (n)) | ((x) >> (32 - (n)))) + +static void ref_chacha_block(uint32_t out_words[16], const uint8_t key[32], + uint32_t counter, const uint8_t nonce[12]) +{ + static const uint32_t c[4] = { + 0x61707865, 0x3320646e, 0x79622d32, 0x6b206574 + }; + uint32_t s[16], x[16]; + int i; + + for (i = 0; i < 4; i++) + s[i] = c[i]; + for (i = 0; i < 8; i++) + s[4 + i] = (uint32_t)key[4 * i] | (uint32_t)key[4 * i + 1] << 8 | + (uint32_t)key[4 * i + 2] << 16 | + (uint32_t)key[4 * i + 3] << 24; + s[12] = counter; + for (i = 0; i < 3; i++) + s[13 + i] = (uint32_t)nonce[4 * i] | (uint32_t)nonce[4 * i + 1] << 8 | + (uint32_t)nonce[4 * i + 2] << 16 | + (uint32_t)nonce[4 * i + 3] << 24; + + memcpy(x, s, sizeof x); +#define QR(a, b, cc, d) \ + x[a] += x[b]; x[d] ^= x[a]; x[d] = ROTL32(x[d], 16); \ + x[cc] += x[d]; x[b] ^= x[cc]; x[b] = ROTL32(x[b], 12); \ + x[a] += x[b]; x[d] ^= x[a]; x[d] = ROTL32(x[d], 8); \ + x[cc] += x[d]; x[b] ^= x[cc]; x[b] = ROTL32(x[b], 7) + for (i = 0; i < 10; i++) { + QR(0, 4, 8, 12); QR(1, 5, 9, 13); + QR(2, 6, 10, 14); QR(3, 7, 11, 15); + QR(0, 5, 10, 15); QR(1, 6, 11, 12); + QR(2, 7, 8, 13); QR(3, 4, 9, 14); + } +#undef QR + for (i = 0; i < 16; i++) + out_words[i] = x[i] + s[i]; +} + +/* ---------- checks ---------- */ + +static void check_int(void) +{ + /* determinism and non-triviality: the checksum must be stable and + * must actually change with the iteration count */ + uint64_t a = fm_int_math(1000); + uint64_t b = fm_int_math(1000); + uint64_t c = fm_int_math(2000); + + ok("int_math is deterministic", a == b); + ok("int_math depends on iters", a != c); + ok("int_math(0) is zero", fm_int_math(0) == 0); +} + +static void check_fp(void) +{ + uint64_t a = fm_fp_math(1000); + uint64_t b = fm_fp_math(1000); + double da; + + memcpy(&da, &a, sizeof da); + ok("fp_math is deterministic", a == b); + ok("fp_math result is finite", isfinite(da)); + ok("fp_math(0) is zero", fm_fp_math(0) == 0); +} + +static void check_primes(void) +{ + enum { LIM = 1000000 }; + uint8_t *sieve = malloc(LIM); + uint64_t got = fm_primes(LIM, sieve); + uint64_t ref = ref_prime_count(LIM); + + note(" primes < %d: got %llu, expected %llu\n", + LIM, (unsigned long long)got, (unsigned long long)ref); + ok("primes matches reference sieve", got == ref); + ok("primes < 10 == 4", fm_primes(10, sieve) == 4); /* 2,3,5,7 */ + ok("primes < 2 == 0", fm_primes(2, sieve) == 0); + free(sieve); +} + +static void check_simd(void) +{ + uint8_t *buf = aligned_alloc(16, 256); + uint64_t a, b; + + memset(buf, 0xA5, 256); + a = fm_simd(500, buf); + memset(buf, 0xA5, 256); + b = fm_simd(500, buf); + ok("simd is deterministic", a == b); + ok("simd(0) is zero", fm_simd(0, buf) == 0); + free(buf); +} + +static void check_compress(void) +{ + enum { N = 65536 }; + uint8_t *src = malloc(N); + uint32_t *ht = malloc((1 << 16) * sizeof(uint32_t)); + uint64_t incompressible, compressible; + size_t i; + + /* genuinely incompressible data (splitmix64 output): with no matches + * to exploit, an LZ coder's output must be at least the input size */ + { + uint64_t st = 0x1234567890abcdefULL; + for (i = 0; i < N; i++) { + uint64_t z = (st += 0x9e3779b97f4a7c15ULL); + z = (z ^ (z >> 30)) * 0xbf58476d1ce4e5b9ULL; + z = (z ^ (z >> 27)) * 0x94d049bb133111ebULL; + src[i] = (uint8_t)(z ^ (z >> 31)); + } + } + incompressible = fm_compress(src, N, ht); + + /* all-zero data is maximally compressible: it must shrink hugely */ + memset(src, 0, N); + compressible = fm_compress(src, N, ht); + + note(" 64KiB random -> %llu bytes, 64KiB zeros -> %llu bytes\n", + (unsigned long long)incompressible, + (unsigned long long)compressible); + ok("compress expands random data", incompressible >= N); + ok("compress shrinks constant data", compressible < N / 10); + ok("compress is deterministic", fm_compress(src, N, ht) == compressible); + free(src); + free(ht); +} + +static void check_crypto(void) +{ + uint8_t key[32]; + size_t i; + + /* (1) anchor the scalar reference to the RFC 8439 s.2.3.2 vector: + * key = 00,01,...,1f; counter = 1; nonce = 00,00,00,09,...,4a,... + * serialised keystream begins 10 f1 e7 e4. */ + { + uint32_t w[16]; + uint8_t rnonce[12] = {0,0,0,9, 0,0,0,0x4a, 0,0,0,0}; + uint8_t ks0[4]; + for (i = 0; i < 32; i++) + key[i] = (uint8_t)i; + ref_chacha_block(w, key, 1, rnonce); + for (i = 0; i < 4; i++) + ks0[i] = (uint8_t)(w[0] >> (8 * i)); + note(" ref keystream[0..3] = %02x %02x %02x %02x " + "(RFC 8439 expects 10 f1 e7 e4)\n", + ks0[0], ks0[1], ks0[2], ks0[3]); + ok("scalar ChaCha20 matches RFC 8439 vector", + ks0[0] == 0x10 && ks0[1] == 0xf1 && + ks0[2] == 0xe7 && ks0[3] == 0xe4); + } + + /* (2) validate the NEON kernel against that reference. The kernel + * hardwires nonce = 0 and starts the block counter at 0, so we + * compare its keystream to the reference block-for-block. */ + { + uint8_t buf[128]; + uint8_t zero_nonce[12] = {0}; + uint32_t ref0[16], ref1[16]; + int match = 1; + + for (i = 0; i < 32; i++) + key[i] = (uint8_t)(i * 5 + 1); + memset(buf, 0, sizeof buf); /* zeros -> raw keystream */ + fm_chacha20(buf, sizeof buf, key, 1); + + ref_chacha_block(ref0, key, 0, zero_nonce); + ref_chacha_block(ref1, key, 1, zero_nonce); + for (i = 0; i < 16; i++) { + uint32_t k0 = (uint32_t)buf[4 * i] | + (uint32_t)buf[4 * i + 1] << 8 | + (uint32_t)buf[4 * i + 2] << 16 | + (uint32_t)buf[4 * i + 3] << 24; + uint32_t k1 = (uint32_t)buf[64 + 4 * i] | + (uint32_t)buf[64 + 4 * i + 1] << 8 | + (uint32_t)buf[64 + 4 * i + 2] << 16 | + (uint32_t)buf[64 + 4 * i + 3] << 24; + if (k0 != ref0[i] || k1 != ref1[i]) + match = 0; + } + ok("NEON ChaCha20 matches scalar reference (2 blocks)", match); + } + + /* (3) the cipher is a real XOR stream: applying it twice is identity */ + { + uint8_t plain[128], work[128], k2[32]; + for (i = 0; i < 128; i++) + plain[i] = (uint8_t)(i * 7 + 1); + for (i = 0; i < 32; i++) + k2[i] = (uint8_t)(i * 3); + memcpy(work, plain, 128); + fm_chacha20(work, 128, k2, 1); + ok("chacha20 actually changes data", memcmp(work, plain, 128) != 0); + fm_chacha20(work, 128, k2, 1); + ok("chacha20 round-trips (XOR is involutive)", + memcmp(work, plain, 128) == 0); + } +} + +static void check_physics(void) +{ + /* two equal masses released from rest must accelerate toward each + * other: symmetric, momentum-conserving, and bounded. */ + double bodies[2 * 8] = {0}; + double total_p; + + bodies[0] = -1.0; bodies[3] = 1.0; /* body 0 at x=-1, mass 1 */ + bodies[8] = 1.0; bodies[11] = 1.0; /* body 1 at x=+1, mass 1 */ + + fm_physics(bodies, 2, 200); + + /* velocities must be equal and opposite (Newton's third law) */ + total_p = bodies[4] + bodies[12]; /* vx0 + vx1 */ + note(" 2-body: vx0=%.6f vx1=%.6f (sum should be ~0)\n", + bodies[4], bodies[12]); + ok("physics conserves momentum", fabs(total_p) < 1e-9); + ok("physics: bodies attract", bodies[4] > 0.0 && bodies[12] < 0.0); + ok("physics values stay finite", isfinite(bodies[0]) && isfinite(bodies[8])); +} + +static int cmp_u32(const void *p, const void *q) +{ + uint32_t x = *(const uint32_t *)p, y = *(const uint32_t *)q; + return (x > y) - (x < y); +} + +static int is_sorted(const uint32_t *a, size_t n) +{ + for (size_t i = 1; i < n; i++) + if (a[i - 1] > a[i]) + return 0; + return 1; +} + +static void check_sort(void) +{ + enum { N = 10000 }; + uint32_t *a = malloc(N * sizeof(uint32_t)); + uint32_t *b = malloc(N * sizeof(uint32_t)); + uint64_t s; + size_t i; + uint32_t r = 12345; + + for (i = 0; i < N; i++) { + r = r * 1103515245u + 12345u; + a[i] = r; + } + memcpy(b, a, N * sizeof(uint32_t)); + + s = fm_sort(a, N); + ok("sort produces sorted output", is_sorted(a, N)); + + /* multiset is preserved: sort the reference with the C library and + * compare element by element */ + qsort(b, N, sizeof(uint32_t), cmp_u32); + ok("sort is a permutation of the input", + memcmp(a, b, N * sizeof(uint32_t)) == 0); + + /* already-sorted input stays sorted and gives the same checksum */ + { + uint64_t s2 = fm_sort(a, N); + ok("sort is idempotent on sorted data", + is_sorted(a, N) && s2 == s); + } + + ok("sort of empty array is zero", fm_sort(a, 0) == 0); + free(a); + free(b); +} + +static void check_chase(void) +{ + /* build a tiny 4-node cycle by hand and confirm the walk returns to + * the start after exactly `n` steps (offset 0 relative to entry) */ + void *nodes[4]; + + nodes[0] = &nodes[1]; + nodes[1] = &nodes[2]; + nodes[2] = &nodes[3]; + nodes[3] = &nodes[0]; + + /* 4 hops from &nodes[0] returns to &nodes[0]; fm_chase returns the + * final pointer minus the starting pointer, so a full loop gives 0 */ + ok("chase completes a full cycle", fm_chase(nodes, 4) == 0); + ok("chase(0) is zero", fm_chase(nodes, 0) == 0); + /* one hop lands on &nodes[1], i.e. one pointer-width past the start */ + ok("chase single hop offset", + fm_chase(nodes, 1) == (uint64_t)((char *)&nodes[1] - (char *)&nodes[0])); +} + +int main(void) +{ + long n = sysconf(_SC_NPROCESSORS_ONLN); + + fm_ncores = n > 0 ? n : 1; + + printf("\nfossmark kernel correctness tests\n"); + printf("=================================\n"); + printf("running each check on %ld core%s in parallel\n\n", + fm_ncores, fm_ncores == 1 ? "" : "s"); + + printf("Integer Math:\n"); parallel(check_int); + printf("Floating Point Math:\n"); parallel(check_fp); + printf("Prime Numbers:\n"); parallel(check_primes); + printf("Extended Instructions:\n"); parallel(check_simd); + printf("Compression:\n"); parallel(check_compress); + printf("Encryption:\n"); parallel(check_crypto); + printf("Physics:\n"); parallel(check_physics); + printf("Sorting:\n"); parallel(check_sort); + /* the pointer chase is the single-threaded test: run it on one core */ + printf("Single-Threaded (chase):\n"); check_chase(); + + printf("\n=================================\n"); + printf("%d checks, %d failures\n\n", checks, failures); + return failures ? 1 : 0; +}