master
1//! The engines provided here should be initialized from an external source.
2//! For a thread-local cryptographically secure pseudo random number generator,
3//! use `std.crypto.random`.
4//! Be sure to use a CSPRNG when required, otherwise using a normal PRNG will
5//! be faster and use substantially less stack space.
6
7const std = @import("std.zig");
8const math = std.math;
9const mem = std.mem;
10const assert = std.debug.assert;
11const maxInt = std.math.maxInt;
12const Random = @This();
13
14/// Fast unbiased random numbers.
15pub const DefaultPrng = Xoshiro256;
16
17/// Cryptographically secure random numbers.
18pub const DefaultCsprng = ChaCha;
19
20pub const Ascon = @import("Random/Ascon.zig");
21pub const ChaCha = @import("Random/ChaCha.zig");
22
23pub const Isaac64 = @import("Random/Isaac64.zig");
24pub const Pcg = @import("Random/Pcg.zig");
25pub const Xoroshiro128 = @import("Random/Xoroshiro128.zig");
26pub const Xoshiro256 = @import("Random/Xoshiro256.zig");
27pub const Sfc64 = @import("Random/Sfc64.zig");
28pub const RomuTrio = @import("Random/RomuTrio.zig");
29pub const SplitMix64 = @import("Random/SplitMix64.zig");
30pub const ziggurat = @import("Random/ziggurat.zig");
31
32/// Any comparison of this field may result in illegal behavior, since it may be set to
33/// `undefined` in cases where the random implementation does not have any associated
34/// state.
35ptr: *anyopaque,
36fillFn: *const fn (ptr: *anyopaque, buf: []u8) void,
37
38pub fn init(pointer: anytype, comptime fillFn: fn (ptr: @TypeOf(pointer), buf: []u8) void) Random {
39 const Ptr = @TypeOf(pointer);
40 assert(@typeInfo(Ptr) == .pointer); // Must be a pointer
41 assert(@typeInfo(Ptr).pointer.size == .one); // Must be a single-item pointer
42 assert(@typeInfo(@typeInfo(Ptr).pointer.child) == .@"struct"); // Must point to a struct
43 const gen = struct {
44 fn fill(ptr: *anyopaque, buf: []u8) void {
45 const self: Ptr = @ptrCast(@alignCast(ptr));
46 fillFn(self, buf);
47 }
48 };
49
50 return .{
51 .ptr = pointer,
52 .fillFn = gen.fill,
53 };
54}
55
56/// Read random bytes into the specified buffer until full.
57pub fn bytes(r: Random, buf: []u8) void {
58 r.fillFn(r.ptr, buf);
59}
60
61pub fn array(r: Random, comptime E: type, comptime N: usize) [N]E {
62 var result: [N]E = undefined;
63 bytes(r, &result);
64 return result;
65}
66
67pub fn boolean(r: Random) bool {
68 return r.int(u1) != 0;
69}
70
71/// Returns a random value from an enum, evenly distributed.
72///
73/// Note that this will not yield consistent results across all targets
74/// due to dependence on the representation of `usize` as an index.
75/// See `enumValueWithIndex` for further commentary.
76pub inline fn enumValue(r: Random, comptime EnumType: type) EnumType {
77 return r.enumValueWithIndex(EnumType, usize);
78}
79
80/// Returns a random value from an enum, evenly distributed.
81///
82/// An index into an array of all named values is generated using the
83/// specified `Index` type to determine the return value.
84/// This allows for results to be independent of `usize` representation.
85///
86/// Prefer `enumValue` if this isn't important.
87///
88/// See `uintLessThan`, which this function uses in most cases,
89/// for commentary on the runtime of this function.
90pub fn enumValueWithIndex(r: Random, comptime EnumType: type, comptime Index: type) EnumType {
91 comptime assert(@typeInfo(EnumType) == .@"enum");
92
93 // We won't use int -> enum casting because enum elements can have
94 // arbitrary values. Instead we'll randomly pick one of the type's values.
95 const values = comptime std.enums.values(EnumType);
96 comptime assert(values.len > 0); // can't return anything
97 comptime assert(maxInt(Index) >= values.len - 1); // can't access all values
98 if (values.len == 1) return values[0];
99
100 const index = if (comptime values.len - 1 == maxInt(Index))
101 r.int(Index)
102 else
103 r.uintLessThan(Index, values.len);
104
105 const MinInt = MinArrayIndex(Index);
106 return values[@as(MinInt, @intCast(index))];
107}
108
109/// Returns a random int `i` such that `minInt(T) <= i <= maxInt(T)`.
110/// `i` is evenly distributed.
111pub fn int(r: Random, comptime T: type) T {
112 const bits = @typeInfo(T).int.bits;
113 const UnsignedT = std.meta.Int(.unsigned, bits);
114 const ceil_bytes = comptime std.math.divCeil(u16, bits, 8) catch unreachable;
115 const ByteAlignedT = std.meta.Int(.unsigned, ceil_bytes * 8);
116
117 var rand_bytes: [ceil_bytes]u8 = undefined;
118 r.bytes(&rand_bytes);
119
120 // use LE instead of native endian for better portability maybe?
121 // TODO: endian portability is pointless if the underlying prng isn't endian portable.
122 // TODO: document the endian portability of this library.
123 const byte_aligned_result = mem.readInt(ByteAlignedT, &rand_bytes, .little);
124 const unsigned_result: UnsignedT = @truncate(byte_aligned_result);
125 return @bitCast(unsigned_result);
126}
127
128/// Constant-time implementation off `uintLessThan`.
129/// The results of this function may be biased.
130pub fn uintLessThanBiased(r: Random, comptime T: type, less_than: T) T {
131 comptime assert(@typeInfo(T).int.signedness == .unsigned);
132 assert(0 < less_than);
133 return limitRangeBiased(T, r.int(T), less_than);
134}
135
136/// Returns an evenly distributed random unsigned integer `0 <= i < less_than`.
137/// This function assumes that the underlying `fillFn` produces evenly distributed values.
138/// Within this assumption, the runtime of this function is exponentially distributed.
139/// If `fillFn` were backed by a true random generator,
140/// the runtime of this function would technically be unbounded.
141/// However, if `fillFn` is backed by any evenly distributed pseudo random number generator,
142/// this function is guaranteed to return.
143/// If you need deterministic runtime bounds, use `uintLessThanBiased`.
144pub fn uintLessThan(r: Random, comptime T: type, less_than: T) T {
145 comptime assert(@typeInfo(T).int.signedness == .unsigned);
146 const bits = @typeInfo(T).int.bits;
147 assert(0 < less_than);
148
149 // adapted from:
150 // http://www.pcg-random.org/posts/bounded-rands.html
151 // "Lemire's (with an extra tweak from me)"
152 var x = r.int(T);
153 var m = math.mulWide(T, x, less_than);
154 var l: T = @truncate(m);
155 if (l < less_than) {
156 var t = -%less_than;
157
158 if (t >= less_than) {
159 t -= less_than;
160 if (t >= less_than) {
161 t %= less_than;
162 }
163 }
164 while (l < t) {
165 x = r.int(T);
166 m = math.mulWide(T, x, less_than);
167 l = @truncate(m);
168 }
169 }
170 return @intCast(m >> bits);
171}
172
173/// Constant-time implementation off `uintAtMost`.
174/// The results of this function may be biased.
175pub fn uintAtMostBiased(r: Random, comptime T: type, at_most: T) T {
176 assert(@typeInfo(T).int.signedness == .unsigned);
177 if (at_most == maxInt(T)) {
178 // have the full range
179 return r.int(T);
180 }
181 return r.uintLessThanBiased(T, at_most + 1);
182}
183
184/// Returns an evenly distributed random unsigned integer `0 <= i <= at_most`.
185/// See `uintLessThan`, which this function uses in most cases,
186/// for commentary on the runtime of this function.
187pub fn uintAtMost(r: Random, comptime T: type, at_most: T) T {
188 assert(@typeInfo(T).int.signedness == .unsigned);
189 if (at_most == maxInt(T)) {
190 // have the full range
191 return r.int(T);
192 }
193 return r.uintLessThan(T, at_most + 1);
194}
195
196/// Constant-time implementation off `intRangeLessThan`.
197/// The results of this function may be biased.
198pub fn intRangeLessThanBiased(r: Random, comptime T: type, at_least: T, less_than: T) T {
199 assert(at_least < less_than);
200 const info = @typeInfo(T).int;
201 if (info.signedness == .signed) {
202 // Two's complement makes this math pretty easy.
203 const UnsignedT = std.meta.Int(.unsigned, info.bits);
204 const lo: UnsignedT = @bitCast(at_least);
205 const hi: UnsignedT = @bitCast(less_than);
206 const result = lo +% r.uintLessThanBiased(UnsignedT, hi -% lo);
207 return @bitCast(result);
208 } else {
209 // The signed implementation would work fine, but we can use stricter arithmetic operators here.
210 return at_least + r.uintLessThanBiased(T, less_than - at_least);
211 }
212}
213
214/// Returns an evenly distributed random integer `at_least <= i < less_than`.
215/// See `uintLessThan`, which this function uses in most cases,
216/// for commentary on the runtime of this function.
217pub fn intRangeLessThan(r: Random, comptime T: type, at_least: T, less_than: T) T {
218 assert(at_least < less_than);
219 const info = @typeInfo(T).int;
220 if (info.signedness == .signed) {
221 // Two's complement makes this math pretty easy.
222 const UnsignedT = std.meta.Int(.unsigned, info.bits);
223 const lo: UnsignedT = @bitCast(at_least);
224 const hi: UnsignedT = @bitCast(less_than);
225 const result = lo +% r.uintLessThan(UnsignedT, hi -% lo);
226 return @bitCast(result);
227 } else {
228 // The signed implementation would work fine, but we can use stricter arithmetic operators here.
229 return at_least + r.uintLessThan(T, less_than - at_least);
230 }
231}
232
233/// Constant-time implementation off `intRangeAtMostBiased`.
234/// The results of this function may be biased.
235pub fn intRangeAtMostBiased(r: Random, comptime T: type, at_least: T, at_most: T) T {
236 assert(at_least <= at_most);
237 const info = @typeInfo(T).int;
238 if (info.signedness == .signed) {
239 // Two's complement makes this math pretty easy.
240 const UnsignedT = std.meta.Int(.unsigned, info.bits);
241 const lo: UnsignedT = @bitCast(at_least);
242 const hi: UnsignedT = @bitCast(at_most);
243 const result = lo +% r.uintAtMostBiased(UnsignedT, hi -% lo);
244 return @bitCast(result);
245 } else {
246 // The signed implementation would work fine, but we can use stricter arithmetic operators here.
247 return at_least + r.uintAtMostBiased(T, at_most - at_least);
248 }
249}
250
251/// Returns an evenly distributed random integer `at_least <= i <= at_most`.
252/// See `uintLessThan`, which this function uses in most cases,
253/// for commentary on the runtime of this function.
254pub fn intRangeAtMost(r: Random, comptime T: type, at_least: T, at_most: T) T {
255 assert(at_least <= at_most);
256 const info = @typeInfo(T).int;
257 if (info.signedness == .signed) {
258 // Two's complement makes this math pretty easy.
259 const UnsignedT = std.meta.Int(.unsigned, info.bits);
260 const lo: UnsignedT = @bitCast(at_least);
261 const hi: UnsignedT = @bitCast(at_most);
262 const result = lo +% r.uintAtMost(UnsignedT, hi -% lo);
263 return @bitCast(result);
264 } else {
265 // The signed implementation would work fine, but we can use stricter arithmetic operators here.
266 return at_least + r.uintAtMost(T, at_most - at_least);
267 }
268}
269
270/// Return a floating point value evenly distributed in the range [0, 1).
271pub fn float(r: Random, comptime T: type) T {
272 // Generate a uniformly random value for the mantissa.
273 // Then generate an exponentially biased random value for the exponent.
274 // This covers every possible value in the range.
275 switch (T) {
276 f32 => {
277 // Use 23 random bits for the mantissa, and the rest for the exponent.
278 // If all 41 bits are zero, generate additional random bits, until a
279 // set bit is found, or 126 bits have been generated.
280 const rand = r.int(u64);
281 var rand_lz = @clz(rand);
282 if (rand_lz >= 41) {
283 @branchHint(.unlikely);
284 rand_lz = 41 + @clz(r.int(u64));
285 if (rand_lz == 41 + 64) {
286 @branchHint(.unlikely);
287 // It is astronomically unlikely to reach this point.
288 rand_lz += @clz(r.int(u32) | 0x7FF);
289 }
290 }
291 const mantissa: u23 = @truncate(rand);
292 const exponent = @as(u32, 126 - rand_lz) << 23;
293 return @bitCast(exponent | mantissa);
294 },
295 f64 => {
296 // Use 52 random bits for the mantissa, and the rest for the exponent.
297 // If all 12 bits are zero, generate additional random bits, until a
298 // set bit is found, or 1022 bits have been generated.
299 const rand = r.int(u64);
300 var rand_lz: u64 = @clz(rand);
301 if (rand_lz >= 12) {
302 rand_lz = 12;
303 while (true) {
304 // It is astronomically unlikely for this loop to execute more than once.
305 const addl_rand_lz = @clz(r.int(u64));
306 rand_lz += addl_rand_lz;
307 if (addl_rand_lz != 64) {
308 @branchHint(.likely);
309 break;
310 }
311 if (rand_lz >= 1022) {
312 rand_lz = 1022;
313 break;
314 }
315 }
316 }
317 const mantissa = rand & 0xFFFFFFFFFFFFF;
318 const exponent = (1022 - rand_lz) << 52;
319 return @bitCast(exponent | mantissa);
320 },
321 else => @compileError("unknown floating point type"),
322 }
323}
324
325/// Return a floating point value normally distributed with mean = 0, stddev = 1.
326///
327/// To use different parameters, use: floatNorm(...) * desiredStddev + desiredMean.
328pub fn floatNorm(r: Random, comptime T: type) T {
329 const value = ziggurat.next_f64(r, ziggurat.NormDist);
330 switch (T) {
331 f32 => return @floatCast(value),
332 f64 => return value,
333 else => @compileError("unknown floating point type"),
334 }
335}
336
337/// Return an exponentially distributed float with a rate parameter of 1.
338///
339/// To use a different rate parameter, use: floatExp(...) / desiredRate.
340pub fn floatExp(r: Random, comptime T: type) T {
341 const value = ziggurat.next_f64(r, ziggurat.ExpDist);
342 switch (T) {
343 f32 => return @floatCast(value),
344 f64 => return value,
345 else => @compileError("unknown floating point type"),
346 }
347}
348
349/// Shuffle a slice into a random order.
350///
351/// Note that this will not yield consistent results across all targets
352/// due to dependence on the representation of `usize` as an index.
353/// See `shuffleWithIndex` for further commentary.
354pub inline fn shuffle(r: Random, comptime T: type, buf: []T) void {
355 r.shuffleWithIndex(T, buf, usize);
356}
357
358/// Shuffle a slice into a random order, using an index of a
359/// specified type to maintain distribution across targets.
360/// Asserts the index type can represent `buf.len`.
361///
362/// Indexes into the slice are generated using the specified `Index`
363/// type, which determines distribution properties. This allows for
364/// results to be independent of `usize` representation.
365///
366/// Prefer `shuffle` if this isn't important.
367///
368/// See `intRangeLessThan`, which this function uses,
369/// for commentary on the runtime of this function.
370pub fn shuffleWithIndex(r: Random, comptime T: type, buf: []T, comptime Index: type) void {
371 const MinInt = MinArrayIndex(Index);
372 if (buf.len < 2) {
373 return;
374 }
375
376 // `i <= j < max <= maxInt(MinInt)`
377 const max: MinInt = @intCast(buf.len);
378 var i: MinInt = 0;
379 while (i < max - 1) : (i += 1) {
380 const j: MinInt = @intCast(r.intRangeLessThan(Index, i, max));
381 mem.swap(T, &buf[i], &buf[j]);
382 }
383}
384
385/// Randomly selects an index into `proportions`, where the likelihood of each
386/// index is weighted by that proportion.
387/// It is more likely for the index of the last proportion to be returned
388/// than the index of the first proportion in the slice, and vice versa.
389///
390/// This is useful for selecting an item from a slice where weights are not equal.
391/// `T` must be a numeric type capable of holding the sum of `proportions`.
392pub fn weightedIndex(r: Random, comptime T: type, proportions: []const T) usize {
393 // This implementation works by summing the proportions and picking a
394 // random point in [0, sum). We then loop over the proportions,
395 // accumulating until our accumulator is greater than the random point.
396
397 const sum = s: {
398 var sum: T = 0;
399 for (proportions) |v| sum += v;
400 break :s sum;
401 };
402
403 const point = switch (@typeInfo(T)) {
404 .int => |int_info| switch (int_info.signedness) {
405 .signed => r.intRangeLessThan(T, 0, sum),
406 .unsigned => r.uintLessThan(T, sum),
407 },
408 // take care that imprecision doesn't lead to a value slightly greater than sum
409 .float => @min(r.float(T) * sum, sum - std.math.floatEps(T)),
410 else => @compileError("weightedIndex does not support proportions of type " ++
411 @typeName(T)),
412 };
413
414 assert(point < sum);
415
416 var accumulator: T = 0;
417 for (proportions, 0..) |p, index| {
418 accumulator += p;
419 if (point < accumulator) return index;
420 } else unreachable;
421}
422
423/// Convert a random integer 0 <= random_int <= maxValue(T),
424/// into an integer 0 <= result < less_than.
425/// This function introduces a minor bias.
426pub fn limitRangeBiased(comptime T: type, random_int: T, less_than: T) T {
427 comptime assert(@typeInfo(T).int.signedness == .unsigned);
428 const bits = @typeInfo(T).int.bits;
429
430 // adapted from:
431 // http://www.pcg-random.org/posts/bounded-rands.html
432 // "Integer Multiplication (Biased)"
433 const m = math.mulWide(T, random_int, less_than);
434 return @intCast(m >> bits);
435}
436
437/// Returns the smallest of `Index` and `usize`.
438fn MinArrayIndex(comptime Index: type) type {
439 const index_info = @typeInfo(Index).int;
440 assert(index_info.signedness == .unsigned);
441 return if (index_info.bits >= @typeInfo(usize).int.bits) usize else Index;
442}
443
444test {
445 std.testing.refAllDecls(@This());
446 _ = @import("Random/test.zig");
447}