diff --git a/lib/c/stdlib.zig b/lib/c/stdlib.zig index 6ee0d19fb5..51ead9bee7 100644 --- a/lib/c/stdlib.zig +++ b/lib/c/stdlib.zig @@ -8,6 +8,7 @@ const lldiv_t = std.c.lldiv_t; comptime { _ = @import("stdlib/rand.zig"); + _ = @import("stdlib/drand48.zig"); if (builtin.target.isMuslLibC() or builtin.target.isWasiLibC()) { // Functions specific to musl and wasi-libc. diff --git a/lib/c/stdlib/drand48.zig b/lib/c/stdlib/drand48.zig new file mode 100644 index 0000000000..c1bce90ee0 --- /dev/null +++ b/lib/c/stdlib/drand48.zig @@ -0,0 +1,147 @@ +//! drand48 functions are based off a 48-bit lcg prng: https://pubs.opengroup.org/onlinepubs/9799919799/functions/drand48.html + +const std = @import("std"); +const common = @import("../common.zig"); +const builtin = @import("builtin"); +const Lcg = std.Random.lcg.Wrapping(u48); + +comptime { + if (builtin.target.isMuslLibC() or builtin.target.isWasiLibC()) { + @export(&erand48, .{ .name = "erand48", .linkage = common.linkage, .visibility = common.visibility }); + @export(&jrand48, .{ .name = "jrand48", .linkage = common.linkage, .visibility = common.visibility }); + @export(&nrand48, .{ .name = "nrand48", .linkage = common.linkage, .visibility = common.visibility }); + @export(&drand48, .{ .name = "drand48", .linkage = common.linkage, .visibility = common.visibility }); + @export(&lrand48, .{ .name = "lrand48", .linkage = common.linkage, .visibility = common.visibility }); + @export(&mrand48, .{ .name = "mrand48", .linkage = common.linkage, .visibility = common.visibility }); + @export(&lcong48, .{ .name = "lcong48", .linkage = common.linkage, .visibility = common.visibility }); + @export(&seed48, .{ .name = "seed48", .linkage = common.linkage, .visibility = common.visibility }); + @export(&srand48, .{ .name = "srand48", .linkage = common.linkage, .visibility = common.visibility }); + } +} + +// NOTE: all "magic" numbers and tests are extracted and adapted from the source above + +const default_multiplier = 0x5DEECE66D; +const default_addend = 0xB; + +var lcg: Lcg = .init(0, default_multiplier, default_addend); +var seed48_xi: [3]c_ushort = undefined; + +fn erand48(xsubi: *[3]c_ushort) callconv(.c) f64 { + const xi = @as(u48, @as(u16, @truncate(xsubi[0])) | (@as(u48, @as(u16, @truncate(xsubi[1])))) << 16) | (@as(u48, @as(u16, @truncate(xsubi[2]))) << 32); + + var separate_lcg: Lcg = .init(xi, lcg.a, lcg.c); + const next_xi = separate_lcg.next(); + + xsubi.* = .{ @truncate(next_xi & 0xFFFF), @truncate((next_xi >> 16) & 0xFFFF), @truncate((next_xi >> 32) & 0xFFFF) }; + return @as(f64, next_xi) / @as(f64, std.math.maxInt(u48)); +} + +fn jrand48(xsubi: *[3]c_ushort) callconv(.c) c_long { + const xi = @as(u48, @as(u16, @truncate(xsubi[0])) | (@as(u48, @as(u16, @truncate(xsubi[1])))) << 16) | (@as(u48, @as(u16, @truncate(xsubi[2]))) << 32); + + var separate_lcg: Lcg = .init(xi, lcg.a, lcg.c); + const next_xi = separate_lcg.next(); + + xsubi.* = .{ @truncate(next_xi & 0xFFFF), @truncate((next_xi >> 16) & 0xFFFF), @truncate((next_xi >> 32) & 0xFFFF) }; + return @as(i32, @bitCast(@as(u32, @truncate(next_xi >> 16)))); +} + +fn nrand48(xsubi: *[3]c_ushort) callconv(.c) c_long { + const xi = @as(u48, @as(u16, @truncate(xsubi[0])) | (@as(u48, @as(u16, @truncate(xsubi[1])))) << 16) | (@as(u48, @as(u16, @truncate(xsubi[2]))) << 32); + + var separate_lcg: Lcg = .init(xi, lcg.a, lcg.c); + const next_xi = separate_lcg.next(); + + xsubi.* = .{ @truncate(next_xi & 0xFFFF), @truncate((next_xi >> 16) & 0xFFFF), @truncate((next_xi >> 32) & 0xFFFF) }; + return @intCast(next_xi >> 17); // a c_long is always at least 32-bits, this is never UB +} + +fn drand48() callconv(.c) f64 { + return 2e-48 * @as(f64, lcg.next()); +} + +fn lrand48() callconv(.c) c_long { + return @intCast(lcg.next() >> 17); +} + +fn mrand48() callconv(.c) c_long { + return @as(i32, @bitCast(@as(u32, @truncate(lcg.next() >> 16)))); +} + +// 0..3 is `Xi`, 3..6 is `a`, 6 is `c` +// first low 16-bits, then mid, then high. +fn lcong48(param: *[7]c_ushort) callconv(.c) void { + lcg.xi = (@as(u48, @as(u16, @truncate(param[0]))) | (@as(u48, @as(u16, @truncate(param[1])))) << 16) | (@as(u48, @as(u16, @truncate(param[2]))) << 32); + lcg.a = (@as(u48, @as(u16, @truncate(param[3]))) | (@as(u48, @as(u16, @truncate(param[4])))) << 16) | (@as(u48, @as(u16, @truncate(param[5]))) << 32); + lcg.c = @as(u16, @truncate(param[6])); +} + +fn seed48(seed16v: *[3]c_ushort) callconv(.c) *[3]c_ushort { + seed48_xi = .{ @truncate(lcg.xi & 0xFFFF), @truncate((lcg.xi >> 16) & 0xFFFF), @truncate((lcg.xi >> 32) & 0xFFFF) }; + const xi = (@as(u48, @as(u16, @truncate(seed16v[0]))) | (@as(u48, @as(u16, @truncate(seed16v[1])))) << 16) | (@as(u48, @as(u16, @truncate(seed16v[2]))) << 32); + lcg = .init(xi, default_multiplier, default_addend); + return &seed48_xi; +} + +fn srand48(seedval: c_long) callconv(.c) void { + const xi = (@as(u32, @truncate(@as(c_ulong, @bitCast(seedval)))) << 16) | 0x330E; + lcg = .init(xi, default_multiplier, default_addend); +} + +test erand48 { + var xsubi: [3]c_ushort = .{ 37174, 64810, 11603 }; + + try std.testing.expectApproxEqAbs(0.8965, erand48(&xsubi), 0.0005); + try std.testing.expectEqualSlices(c_ushort, &.{ 22537, 47966, 58735 }, &xsubi); + + try std.testing.expectApproxEqAbs(0.3375, erand48(&xsubi), 0.0005); + try std.testing.expectEqualSlices(c_ushort, &.{ 37344, 32911, 22119 }, &xsubi); + + try std.testing.expectApproxEqAbs(0.6475, erand48(&xsubi), 0.0005); + try std.testing.expectEqualSlices(c_ushort, &.{ 23659, 29872, 42445 }, &xsubi); + + try std.testing.expectApproxEqAbs(0.5005, erand48(&xsubi), 0.0005); + try std.testing.expectEqualSlices(c_ushort, &.{ 31642, 7875, 32802 }, &xsubi); + + try std.testing.expectApproxEqAbs(0.5065, erand48(&xsubi), 0.0005); + try std.testing.expectEqualSlices(c_ushort, &.{ 64669, 14399, 33170 }, &xsubi); +} + +test jrand48 { + var xsubi: [3]c_ushort = .{ 25175, 11052, 45015 }; + + try std.testing.expectEqual(1699503220, jrand48(&xsubi)); + try std.testing.expectEqualSlices(c_ushort, &.{ 2326, 23668, 25932 }, &xsubi); + + try std.testing.expectEqual(-992276007, jrand48(&xsubi)); + try std.testing.expectEqualSlices(c_ushort, &.{ 41577, 4569, 50395 }, &xsubi); + + try std.testing.expectEqual(-19535776, jrand48(&xsubi)); + try std.testing.expectEqualSlices(c_ushort, &.{ 31936, 59488, 65237 }, &xsubi); + + try std.testing.expectEqual(79438377, jrand48(&xsubi)); + try std.testing.expectEqualSlices(c_ushort, &.{ 40395, 8745, 1212 }, &xsubi); + + try std.testing.expectEqual(-1258917728, jrand48(&xsubi)); + try std.testing.expectEqualSlices(c_ushort, &.{ 37242, 28832, 46326 }, &xsubi); +} + +test nrand48 { + var xsubi: [3]c_ushort = .{ 546, 33817, 23389 }; + + try std.testing.expectEqual(914920692, nrand48(&xsubi)); + try std.testing.expectEqualSlices(c_ushort, &.{ 29829, 10728, 27921 }, &xsubi); + + try std.testing.expectEqual(754104482, nrand48(&xsubi)); + try std.testing.expectEqualSlices(c_ushort, &.{ 6828, 28997, 23013 }, &xsubi); + + try std.testing.expectEqual(609453945, nrand48(&xsubi)); + try std.testing.expectEqualSlices(c_ushort, &.{ 58183, 3826, 18599 }, &xsubi); + + try std.testing.expectEqual(1878644360, nrand48(&xsubi)); + try std.testing.expectEqualSlices(c_ushort, &.{ 36678, 44304, 57331 }, &xsubi); + + try std.testing.expectEqual(2114923686, nrand48(&xsubi)); + try std.testing.expectEqualSlices(c_ushort, &.{ 58585, 22861, 64542 }, &xsubi); +} diff --git a/lib/libc/musl/src/prng/__rand48_step.c b/lib/libc/musl/src/prng/__rand48_step.c deleted file mode 100644 index 94703d0738..0000000000 --- a/lib/libc/musl/src/prng/__rand48_step.c +++ /dev/null @@ -1,14 +0,0 @@ -#include -#include "rand48.h" - -uint64_t __rand48_step(unsigned short *xi, unsigned short *lc) -{ - uint64_t a, x; - x = xi[0] | xi[1]+0U<<16 | xi[2]+0ULL<<32; - a = lc[0] | lc[1]+0U<<16 | lc[2]+0ULL<<32; - x = a*x + lc[3]; - xi[0] = x; - xi[1] = x>>16; - xi[2] = x>>32; - return x & 0xffffffffffffull; -} diff --git a/lib/libc/musl/src/prng/__seed48.c b/lib/libc/musl/src/prng/__seed48.c deleted file mode 100644 index e436b4d302..0000000000 --- a/lib/libc/musl/src/prng/__seed48.c +++ /dev/null @@ -1,3 +0,0 @@ -#include "rand48.h" - -unsigned short __seed48[7] = { 0, 0, 0, 0xe66d, 0xdeec, 0x5, 0xb }; diff --git a/lib/libc/musl/src/prng/drand48.c b/lib/libc/musl/src/prng/drand48.c deleted file mode 100644 index 08283e240b..0000000000 --- a/lib/libc/musl/src/prng/drand48.c +++ /dev/null @@ -1,17 +0,0 @@ -#include -#include -#include "rand48.h" - -double erand48(unsigned short s[3]) -{ - union { - uint64_t u; - double f; - } x = { 0x3ff0000000000000ULL | __rand48_step(s, __seed48+3)<<4 }; - return x.f - 1.0; -} - -double drand48(void) -{ - return erand48(__seed48); -} diff --git a/lib/libc/musl/src/prng/lcong48.c b/lib/libc/musl/src/prng/lcong48.c deleted file mode 100644 index 030e514844..0000000000 --- a/lib/libc/musl/src/prng/lcong48.c +++ /dev/null @@ -1,8 +0,0 @@ -#include -#include -#include "rand48.h" - -void lcong48(unsigned short p[7]) -{ - memcpy(__seed48, p, sizeof __seed48); -} diff --git a/lib/libc/musl/src/prng/lrand48.c b/lib/libc/musl/src/prng/lrand48.c deleted file mode 100644 index 07e2b78440..0000000000 --- a/lib/libc/musl/src/prng/lrand48.c +++ /dev/null @@ -1,13 +0,0 @@ -#include -#include -#include "rand48.h" - -long nrand48(unsigned short s[3]) -{ - return __rand48_step(s, __seed48+3) >> 17; -} - -long lrand48(void) -{ - return nrand48(__seed48); -} diff --git a/lib/libc/musl/src/prng/mrand48.c b/lib/libc/musl/src/prng/mrand48.c deleted file mode 100644 index f4a56e61b2..0000000000 --- a/lib/libc/musl/src/prng/mrand48.c +++ /dev/null @@ -1,13 +0,0 @@ -#include -#include -#include "rand48.h" - -long jrand48(unsigned short s[3]) -{ - return (int32_t)(__rand48_step(s, __seed48+3) >> 16); -} - -long mrand48(void) -{ - return jrand48(__seed48); -} diff --git a/lib/libc/musl/src/prng/rand48.h b/lib/libc/musl/src/prng/rand48.h deleted file mode 100644 index 55cbec1c60..0000000000 --- a/lib/libc/musl/src/prng/rand48.h +++ /dev/null @@ -1,5 +0,0 @@ -#include -#include - -hidden uint64_t __rand48_step(unsigned short *xi, unsigned short *lc); -extern hidden unsigned short __seed48[7]; diff --git a/lib/libc/musl/src/prng/seed48.c b/lib/libc/musl/src/prng/seed48.c deleted file mode 100644 index bce7b339fe..0000000000 --- a/lib/libc/musl/src/prng/seed48.c +++ /dev/null @@ -1,11 +0,0 @@ -#include -#include -#include "rand48.h" - -unsigned short *seed48(unsigned short *s) -{ - static unsigned short p[3]; - memcpy(p, __seed48, sizeof p); - memcpy(__seed48, s, sizeof p); - return p; -} diff --git a/lib/libc/musl/src/prng/srand48.c b/lib/libc/musl/src/prng/srand48.c deleted file mode 100644 index 0a56f6a07a..0000000000 --- a/lib/libc/musl/src/prng/srand48.c +++ /dev/null @@ -1,6 +0,0 @@ -#include - -void srand48(long seed) -{ - seed48((unsigned short [3]){ 0x330e, seed, seed>>16 }); -} diff --git a/src/libs/musl.zig b/src/libs/musl.zig index e7cea46f0c..cf5feb0374 100644 --- a/src/libs/musl.zig +++ b/src/libs/musl.zig @@ -1305,15 +1305,7 @@ const src_files = [_][]const u8{ "musl/src/passwd/putgrent.c", "musl/src/passwd/putpwent.c", "musl/src/passwd/putspent.c", - "musl/src/prng/drand48.c", - "musl/src/prng/lcong48.c", - "musl/src/prng/lrand48.c", - "musl/src/prng/mrand48.c", - "musl/src/prng/__rand48_step.c", "musl/src/prng/random.c", - "musl/src/prng/__seed48.c", - "musl/src/prng/seed48.c", - "musl/src/prng/srand48.c", "musl/src/process/aarch64/vfork.s", "musl/src/process/arm/vfork.s", "musl/src/process/execl.c", diff --git a/src/libs/wasi_libc.zig b/src/libs/wasi_libc.zig index afee029121..c801a805f1 100644 --- a/src/libs/wasi_libc.zig +++ b/src/libs/wasi_libc.zig @@ -885,14 +885,6 @@ const libc_top_half_src_files = [_][]const u8{ "musl/src/network/inet_pton.c", "musl/src/network/ntohl.c", "musl/src/network/ntohs.c", - "musl/src/prng/drand48.c", - "musl/src/prng/lcong48.c", - "musl/src/prng/lrand48.c", - "musl/src/prng/mrand48.c", - "musl/src/prng/__rand48_step.c", - "musl/src/prng/__seed48.c", - "musl/src/prng/seed48.c", - "musl/src/prng/srand48.c", "musl/src/regex/fnmatch.c", "musl/src/regex/regerror.c", "musl/src/search/hsearch.c",