lib/bench/src/stats/compute.zig
daab053ee43316e1809a84551d573ddd1e5bf3d2
1 const std = @import("std");
2 const model = @import("model.zig");
3 const storage_mod = @import("storage.zig");
4
5 pub fn sample(
6 storage: *storage_mod.Storage,
7 samples: []const u64,
8 ) model.ComputeError!model.SampleStats {
9 return sampleWithBootstrap(storage, samples, .{});
10 }
11
12 pub fn sampleWithBootstrap(
13 storage: *storage_mod.Storage,
14 samples: []const u64,
15 bootstrap: model.BootstrapConfig,
16 ) model.ComputeError!model.SampleStats {
17 const regions = try storage.acquire(samples.len, bootstrap);
18 defer storage.reset();
19 @memcpy(regions.sorted, samples);
20 std.mem.sort(u64, regions.sorted, {}, lessThanU64);
21
22 var total: u128 = 0;
23 for (samples) |value| total += value;
24 return .{
25 .min_ns = regions.sorted[0],
26 .max_ns = regions.sorted[regions.sorted.len - 1],
27 .mean_ns = @as(f64, @floatFromInt(total)) /
28 @as(f64, @floatFromInt(samples.len)),
29 .median_ns = percentile(regions.sorted, 50),
30 .p75_ns = percentile(regions.sorted, 75),
31 .p95_ns = percentile(regions.sorted, 95),
32 .p99_ns = percentile(regions.sorted, 99),
33 .total_ns = if (total > std.math.maxInt(u64)) std.math.maxInt(u64) else @intCast(total),
34 .confidence_intervals = bootstrapConfidenceIntervals(regions, bootstrap),
35 };
36 }
37
38 pub fn percentile(sorted: []const u64, pct: u8) u64 {
39 if (sorted.len == 0) return 0;
40 return sorted[scaledIndex(sorted.len, pct, 100)];
41 }
42
43 fn lessThanU64(_: void, left: u64, right: u64) bool {
44 return left < right;
45 }
46
47 fn bootstrapConfidenceIntervals(
48 regions: storage_mod.Regions,
49 config: model.BootstrapConfig,
50 ) ?model.BootstrapConfidenceIntervals {
51 if (config.iterations == 0) return null;
52 var prng = std.Random.DefaultPrng.init(config.seed);
53 const random = prng.random();
54 for (0..config.iterations) |index| {
55 @memset(regions.counts, 0);
56 var total: u128 = 0;
57 for (0..regions.sorted.len) |_| {
58 const sample_index = random.intRangeLessThan(usize, 0, regions.sorted.len);
59 regions.counts[sample_index] += 1;
60 total += regions.sorted[sample_index];
61 }
62 regions.means[index] = @as(f64, @floatFromInt(total)) /
63 @as(f64, @floatFromInt(regions.sorted.len));
64 regions.medians[index] = @floatFromInt(weightedPercentile(
65 regions.sorted,
66 regions.counts,
67 50,
68 ));
69 regions.p75s[index] = @floatFromInt(weightedPercentile(
70 regions.sorted,
71 regions.counts,
72 75,
73 ));
74 regions.p95s[index] = @floatFromInt(weightedPercentile(
75 regions.sorted,
76 regions.counts,
77 95,
78 ));
79 regions.p99s[index] = @floatFromInt(weightedPercentile(
80 regions.sorted,
81 regions.counts,
82 99,
83 ));
84 }
85
86 std.mem.sort(f64, regions.means, {}, std.sort.asc(f64));
87 std.mem.sort(f64, regions.medians, {}, std.sort.asc(f64));
88 std.mem.sort(f64, regions.p75s, {}, std.sort.asc(f64));
89 std.mem.sort(f64, regions.p95s, {}, std.sort.asc(f64));
90 std.mem.sort(f64, regions.p99s, {}, std.sort.asc(f64));
91 return .{
92 .confidence = @as(f64, @floatFromInt(config.confidence_per_mille)) / 1000.0,
93 .iterations = config.iterations,
94 .seed = config.seed,
95 .mean_ns = interval(regions.means, config.confidence_per_mille),
96 .median_ns = interval(regions.medians, config.confidence_per_mille),
97 .p75_ns = interval(regions.p75s, config.confidence_per_mille),
98 .p95_ns = interval(regions.p95s, config.confidence_per_mille),
99 .p99_ns = interval(regions.p99s, config.confidence_per_mille),
100 };
101 }
102
103 fn interval(sorted: []const f64, confidence_per_mille: u16) model.ConfidenceInterval {
104 const tail = (1000 - @as(usize, confidence_per_mille)) / 2;
105 return .{
106 .low_ns = quantilePerMille(sorted, @intCast(tail)),
107 .high_ns = quantilePerMille(sorted, @intCast(1000 - tail)),
108 };
109 }
110
111 fn quantilePerMille(sorted: []const f64, per_mille: u16) f64 {
112 if (sorted.len == 0) return 0;
113 return sorted[scaledIndex(sorted.len, per_mille, 1000)];
114 }
115
116 fn weightedPercentile(sorted: []const u64, counts: []const u32, pct: u8) u64 {
117 std.debug.assert(sorted.len == counts.len);
118 const target = scaledIndex(sorted.len, pct, 100);
119 var seen: usize = 0;
120 for (sorted, counts) |value, count| {
121 seen += count;
122 if (seen > target) return value;
123 }
124 return sorted[sorted.len - 1];
125 }
126
127 fn scaledIndex(len: usize, numerator: usize, denominator: usize) usize {
128 std.debug.assert(len > 0);
129 std.debug.assert(denominator > 0);
130 const product = @as(u128, len - 1) * numerator;
131 return @intCast((product + denominator - 1) / denominator);
132 }