authorgravatar for andrew@ziglang.orgAndrew Kelley <andrew@ziglang.org> 2023-04-10 13:01:28-04:00
committergravatar for noreply@github.comGitHub <noreply@github.com> 2023-04-10 13:01:28-04:00
log311d50f9dc96c3d397081a1afcc3b519ae1e973f
tree6b580b02f960e41f7055ebd5003bb6cf37f48f5a
parenta2793f8ab831b4a9ffba24540e40a862f32b9d90
parentd6ad69d40a3ba32e9f10a12618dbeece2e62a051
signaturebadge-question-mark Signed by PGP key 4AEE18F83AFDEB23

Merge pull request #14827 from adriandelgado/log10-int

std.math: port `int_log10` from Rust

2 files changed, 178 insertions(+), 1 deletions(-)

lib/std/math.zig+1
......@@ -250,6 +250,7 @@ pub const ln = @import("math/ln.zig").ln;
250250pub const log = @import("math/log.zig").log;
251251pub const log2 = @import("math/log2.zig").log2;
252252pub const log10 = @import("math/log10.zig").log10;
253pub const log10_int = @import("math/log10.zig").log10_int;
253254pub const log1p = @import("math/log1p.zig").log1p;
254255pub const asinh = @import("math/asinh.zig").asinh;
255256pub const acosh = @import("math/acosh.zig").acosh;
lib/std/math/log10.zig+177-1
......@@ -1,7 +1,10 @@
11const std = @import("../std.zig");
2const builtin = @import("builtin");
23const math = std.math;
34const testing = std.testing;
45const maxInt = std.math.maxInt;
6const assert = std.debug.assert;
7const Log2Int = std.math.Log2Int;
58
69/// Returns the base-10 logarithm of x.
710///
......@@ -22,8 +25,181 @@ pub fn log10(x: anytype) @TypeOf(x) {
2225 },
2326 .Int => |IntType| switch (IntType.signedness) {
2427 .signed => @compileError("log10 not implemented for signed integers"),
25 .unsigned => return @floatToInt(T, @floor(@log10(@intToFloat(f64, x)))),
28 .unsigned => return log10_int(x),
2629 },
2730 else => @compileError("log10 not implemented for " ++ @typeName(T)),
2831 }
2932}
33
34// Based on Rust, which is licensed under the MIT license.
35// https://github.com/rust-lang/rust/blob/f63ccaf25f74151a5d8ce057904cd944074b01d2/LICENSE-MIT
36//
37// https://github.com/rust-lang/rust/blob/f63ccaf25f74151a5d8ce057904cd944074b01d2/library/core/src/num/int_log10.rs
38
39/// Return the log base 10 of integer value x, rounding down to the
40/// nearest integer.
41pub fn log10_int(x: anytype) Log2Int(@TypeOf(x)) {
42 const T = @TypeOf(x);
43 const OutT = Log2Int(T);
44 if (@typeInfo(T) != .Int or @typeInfo(T).Int.signedness != .unsigned)
45 @compileError("log10_int requires an unsigned integer, found " ++ @typeName(T));
46
47 assert(x != 0);
48
49 const bit_size = @typeInfo(T).Int.bits;
50
51 if (bit_size <= 8) {
52 return @intCast(OutT, log10_int_u8(x));
53 } else if (bit_size <= 16) {
54 return @intCast(OutT, less_than_5(x));
55 }
56
57 var val = x;
58 var log: u32 = 0;
59
60 inline for (0..11) |i| {
61 // Unnecesary branches should be removed by the compiler
62 if (bit_size > (1 << (11 - i)) * 5 * @log2(10.0) and val >= pow10((1 << (11 - i)) * 5)) {
63 const num_digits = (1 << (11 - i)) * 5;
64 val /= pow10(num_digits);
65 log += num_digits;
66 }
67 }
68
69 if (val >= pow10(5)) {
70 val /= pow10(5);
71 log += 5;
72 }
73
74 return @intCast(OutT, log + less_than_5(@intCast(u32, val)));
75}
76
77fn pow10(comptime y: comptime_int) comptime_int {
78 if (y == 0) return 1;
79
80 var squaring = 0;
81 var s = 1;
82
83 while (s <= y) : (s <<= 1) {
84 squaring += 1;
85 }
86
87 squaring -= 1;
88
89 var result = 10;
90
91 for (0..squaring) |_| {
92 result *= result;
93 }
94
95 const rest_exp = y - (1 << squaring);
96
97 return result * pow10(rest_exp);
98}
99
100inline fn log10_int_u8(x: u8) u32 {
101 // For better performance, avoid branches by assembling the solution
102 // in the bits above the low 8 bits.
103
104 // Adding c1 to val gives 10 in the top bits for val < 10, 11 for val >= 10
105 const C1: u32 = 0b11_00000000 - 10; // 758
106 // Adding c2 to val gives 01 in the top bits for val < 100, 10 for val >= 100
107 const C2: u32 = 0b10_00000000 - 100; // 412
108
109 // Value of top bits:
110 // +c1 +c2 1&2
111 // 0..=9 10 01 00 = 0
112 // 10..=99 11 01 01 = 1
113 // 100..=255 11 10 10 = 2
114 return ((x + C1) & (x + C2)) >> 8;
115}
116
117inline fn less_than_5(x: u32) u32 {
118 // Similar to log10u8, when adding one of these constants to val,
119 // we get two possible bit patterns above the low 17 bits,
120 // depending on whether val is below or above the threshold.
121 const C1: u32 = 0b011_00000000000000000 - 10; // 393206
122 const C2: u32 = 0b100_00000000000000000 - 100; // 524188
123 const C3: u32 = 0b111_00000000000000000 - 1000; // 916504
124 const C4: u32 = 0b100_00000000000000000 - 10000; // 514288
125
126 // Value of top bits:
127 // +c1 +c2 1&2 +c3 +c4 3&4 ^
128 // 0..=9 010 011 010 110 011 010 000 = 0
129 // 10..=99 011 011 011 110 011 010 001 = 1
130 // 100..=999 011 100 000 110 011 010 010 = 2
131 // 1000..=9999 011 100 000 111 011 011 011 = 3
132 // 10000..=99999 011 100 000 111 100 100 100 = 4
133 return (((x + C1) & (x + C2)) ^ ((x + C3) & (x + C4))) >> 17;
134}
135
136fn oldlog10(x: anytype) u8 {
137 return @floatToInt(u8, @log10(@intToFloat(f64, x)));
138}
139
140test "oldlog10 doesn't work" {
141 try testing.expect(14 != oldlog10(pow10(15) - 1));
142
143 // log10(10**15 -1) should indeed be 14
144 try testing.expect(14 == log10_int(@as(u64, pow10(15) - 1)));
145}
146
147test "log10_int vs old implementation" {
148 if (builtin.zig_backend == .stage2_wasm) return error.SkipZigTest; // TODO
149 if (builtin.zig_backend == .stage2_c) return error.SkipZigTest; // TODO
150 if (builtin.zig_backend == .stage2_x86_64) return error.SkipZigTest; // TODO
151 if (builtin.zig_backend == .stage2_aarch64) return error.SkipZigTest; // TODO
152 if (builtin.zig_backend == .stage2_sparc64) return error.SkipZigTest; // TODO
153 if (builtin.zig_backend == .stage2_arm) return error.SkipZigTest; // TODO
154 if (builtin.zig_backend == .stage2_llvm and comptime builtin.target.isWasm()) return error.SkipZigTest; // TODO
155
156 const int_types = .{ u8, u16, u32, u64, u128 };
157
158 inline for (int_types) |T| {
159 const last = @min(maxInt(T), 100_000);
160 for (1..last) |i| {
161 const x = @intCast(T, i);
162 try testing.expectEqual(oldlog10(x), log10_int(x));
163 }
164
165 const max_int: T = maxInt(T);
166 try testing.expectEqual(oldlog10(max_int), log10_int(max_int));
167 }
168}
169
170test "log10_int close to powers of 10" {
171 if (builtin.zig_backend == .stage2_wasm) return error.SkipZigTest; // TODO
172 if (builtin.zig_backend == .stage2_c) return error.SkipZigTest; // TODO
173 if (builtin.zig_backend == .stage2_x86_64) return error.SkipZigTest; // TODO
174 if (builtin.zig_backend == .stage2_aarch64) return error.SkipZigTest; // TODO
175 if (builtin.zig_backend == .stage2_sparc64) return error.SkipZigTest; // TODO
176 if (builtin.zig_backend == .stage2_arm) return error.SkipZigTest; // TODO
177 if (builtin.zig_backend == .stage2_llvm and comptime builtin.target.isWasm()) return error.SkipZigTest; // TODO
178
179 const int_types = .{ u8, u16, u32, u64, u128, u256, u512 };
180 const max_log_values: [7]usize = .{ 2, 4, 9, 19, 38, 77, 154 };
181
182 inline for (int_types, max_log_values) |T, expected_max_ilog| {
183 const max_val: T = maxInt(T);
184
185 try testing.expectEqual(expected_max_ilog, log10_int(max_val));
186
187 for (0..(expected_max_ilog + 1)) |idx| {
188 const i = @intCast(T, idx);
189 const p: T = try math.powi(T, 10, i);
190
191 const b = @intCast(Log2Int(T), i);
192
193 if (p >= 10) {
194 try testing.expectEqual(b - 1, log10_int(p - 9));
195 try testing.expectEqual(b - 1, log10_int(p - 1));
196 }
197
198 try testing.expectEqual(b, log10_int(p));
199 try testing.expectEqual(b, log10_int(p + 1));
200 if (p >= 10) {
201 try testing.expectEqual(b, log10_int(p + 9));
202 }
203 }
204 }
205}