1// Microbenchmark: batched realized accumulators (fused), rolling sum of squares, regression sufficient statistics2// and batched GARCH likelihoods — RAID CPU 1 thread vs CPU all threads vs Metal — over series counts and lengths.3// Protocol (charter §8): warm-ups, randomized order, R repetitions, medians + bootstrap CIs computed later from the raw4// records; cold-vs-warm shader compilation measured separately; records written as JSON lines with the manifest.5#include <cstdio>6#include <algorithm>7#include <random>8#include <fstream>9#include <chrono>10#include "raid/backend/backend.hpp"11#include "raid/core/rng.hpp"12#include "raid/core/manifest.hpp"13#include "raid/core/thread_pool.hpp"14#include "raid/estimators/realized.hpp"15using namespace raid; using json = nlohmann::json;1617namespace {18struct Case { std::string kernel; usize S, n; std::string backend; unsigned threads; };19f64 now() { return std::chrono::duration<f64>(std::chrono::steady_clock::now().time_since_epoch()).count(); }20}2122int main(int argc, char** argv) {23 std::string out = argc > 1 ? argv[1] : "experiments/results/raw/bench_primitives.jsonl"; int reps = argc > 2 ? std::atoi(argv[2]) : 7; const u64 seed = 20261102;24 std::vector<usize> Ss = {1, 10, 100, 1000, 10000}, ns = {390, 3900, 39000, 390000};25 const usize max_elems = 200'000'000; // memory bound (1.6 GB FP64 + 0.8 GB FP32)26 std::vector<Case> cases;27 for (auto S : Ss) for (auto n : ns) { if (S * n > max_elems) continue; for (auto k : {"realized_fused", "rolling_sumsq", "regression_stats"}) { cases.push_back({k, S, n, "cpu", 1}); cases.push_back({k, S, n, "cpu", 0}); cases.push_back({k, S, n, "metal", 0}); } }28 for (auto S : {1ul, 10ul, 100ul, 1000ul, 10000ul}) for (auto n : {390ul, 3900ul}) for (auto P : {1ul, 64ul}) { if (S * P > 200000) continue; cases.push_back({"garch_nll_P" + std::to_string(P), S, n, "cpu", 1}); cases.push_back({"garch_nll_P" + std::to_string(P), S, n, "cpu", 0}); cases.push_back({"garch_nll_P" + std::to_string(P), S, n, "metal", 0}); }29 std::mt19937_64 order(seed); std::shuffle(cases.begin(), cases.end(), order);30 ExperimentManifest man; man.experiment = "bench_primitives"; man.created_at_utc = now_utc_iso(); man.raid_version = RAID_VERSION; man.host = host_info(); man.seed = seed; man.parameters = {{"reps", std::to_string(reps)}, {"order", "shuffled"}, {"chunk_max", "4096"}};31 std::ofstream f(out, std::ios::app); if (!f.good()) { std::fprintf(stderr, "cannot write %s\n", out.c_str()); return 1; }32 auto cpu1 = backend::make_cpu_backend(1), cpuN = backend::make_cpu_backend(0);33 std::shared_ptr<backend::Backend> gpu;34#ifdef RAID_ENABLE_METAL35 { auto t0 = now(); gpu = backend::make_metal_backend(); if (gpu) { backend::KernelStats st; std::vector<f64> v(390, 1e-3); std::vector<i64> off = {0, 390}; backend::Batch b{v, off}; std::vector<f64> a, b2, c, d, e; backend::realized_fused(*gpu, b, a, b2, c, d, e, &st);36 json rec = {{"kernel", "shader_library_load"}, {"backend", "metal"}, {"compiled_now", st.compiled_now}, {"seconds", backend::metal_compile_seconds(*gpu)}, {"wall_first_call", st.wall_seconds}, {"note", st.note}, {"t", now() - t0}}; f << rec.dump() << "\n"; } }37#endif38 usize done = 0;39 for (auto& c : cases) {40 if (c.backend == "metal" && !gpu) continue;41 auto be = c.backend == "metal" ? gpu : (c.threads == 1 ? cpu1 : cpuN);42 Rng r(seed, c.S * 1000003 + c.n); std::vector<f64> v(c.S * c.n); for (auto& x : v) x = 1e-3 * r.student_t(5.0); std::vector<i64> off(c.S + 1); for (usize s = 0; s <= c.S; ++s) off[s] = static_cast<i64>(s * c.n);43 backend::Batch b{v, off};44 std::vector<f64> X, y; usize k = 4; usize P = 1; std::vector<f64> params;45 if (c.kernel == "regression_stats") { X.resize(c.S * c.n * k); y.resize(c.S * c.n); for (usize i = 0; i < c.S * c.n; ++i) { X[i * k] = 1; for (usize j = 1; j < k; ++j) X[i * k + j] = r.normal(); y[i] = r.normal(); } }46 if (c.kernel.rfind("garch_nll_P", 0) == 0) { P = std::stoul(c.kernel.substr(11)); params.resize(c.S * P * 4); for (usize s = 0; s < c.S; ++s) for (usize q = 0; q < P; ++q) { f64* th = params.data() + (s * P + q) * 4; th[0] = 0; th[1] = 1e-7; th[2] = 0.05 + 0.001 * q; th[3] = 0.9; } }47 auto run = [&](backend::KernelStats& st) {48 if (c.kernel == "realized_fused") { std::vector<f64> a, b2, cc, d, e;49#ifdef RAID_ENABLE_METAL50 if (c.backend == "metal") { backend::realized_fused(*be, b, a, b2, cc, d, e, &st); return; }51#endif52 backend::KernelStats s1, s2, s3, s4; a = be->sum_squares(b, &s1); b2 = be->sum_abs_adjacent(b, &s2); be->semivariances(b, cc, d, &s3); e = be->sum_fourth(b, &s4); st = s1; st.wall_seconds = s1.wall_seconds + s2.wall_seconds + s3.wall_seconds + s4.wall_seconds; }53 else if (c.kernel == "rolling_sumsq") { auto o = be->rolling_sum_squares(b, 22, &st); }54 else if (c.kernel == "regression_stats") { std::vector<f64> a, bb, cc; be->regression_stats(X, y, off, k, a, bb, cc, &st); }55 else { auto o = be->garch_nll_batch(b, params, P, 4, 0, 0.0, &st); }56 };57 backend::KernelStats warm; for (int w = 0; w < 2; ++w) run(warm);58 for (int rep = 0; rep < reps; ++rep) { backend::KernelStats st; const u64 rss0 = peak_rss_bytes(); const f64 t0 = now(); run(st); const f64 wall = now() - t0;59 json rec = {{"kernel", c.kernel}, {"S", c.S}, {"n", c.n}, {"elements", c.S * c.n}, {"backend", c.backend}, {"threads", c.backend == "metal" ? 0u : (c.threads == 1 ? 1u : ThreadPool::global().size())}, {"rep", rep}, {"wall_seconds", wall}, {"kernel_wall_seconds", st.wall_seconds}, {"device_seconds", st.device_seconds}, {"bytes_in", st.bytes_in}, {"buffer_allocations", st.buffer_allocations}, {"precision", std::string(to_string(st.precision))}, {"note", st.note}, {"peak_rss_before", rss0}, {"peak_rss_after", peak_rss_bytes()}, {"P", P}, {"k", k}, {"window", 22}};60 f << rec.dump() << "\n"; }61 f.flush(); ++done; std::fprintf(stderr, "[%zu/%zu] %s S=%zu n=%zu %s t=%u\n", done, cases.size(), c.kernel.c_str(), c.S, c.n, c.backend.c_str(), c.threads);62 }63 man.save(out + ".manifest.json"); std::fprintf(stderr, "wrote %s (+ manifest)\n", out.c_str()); return 0;64}65