zig/lib/std/fmt/parse_float/FloatInfo.zig
Marc Tiehuis 2085a4af56 add new float-parser based on eisel-lemire algorithm
The previous float-parsing method was lacking in a lot of areas. This
commit introduces a state-of-the art implementation that is both
accurate and fast to std.

Code is derived from working repo https://github.com/tiehuis/zig-parsefloat.
This includes more test-cases and performance numbers that are present
in this commit.

* Accuracy

The primary testing regime has been using test-data found at
https://github.com/tiehuis/parse-number-fxx-test-data. This is a fork of
upstream with support for f128 test-cases added. This data has been
verified against other independent implementations and represents
accurate round-to-even IEEE-754 floating point semantics.

* Performance

Compared to the existing parseFloat implementation there is ~5-10x
performance improvement using the above corpus. (f128 parsing excluded
in below measurements).

** Old

    $ time ./test_all_fxx_data
    3520298/5296694 succeeded (1776396 fail)

    ________________________________________________________
    Executed in   28.68 secs    fish           external
       usr time   28.48 secs    0.00 micros   28.48 secs
       sys time    0.08 secs  694.00 micros    0.08 secs

** This Implementation

    $ time ./test_all_fxx_data
    5296693/5296694 succeeded (1 fail)

    ________________________________________________________
    Executed in    4.54 secs    fish           external
       usr time    4.37 secs  515.00 micros    4.37 secs
       sys time    0.10 secs  171.00 micros    0.10 secs

Further performance numbers can be seen using the
https://github.com/tiehuis/simple_fastfloat_benchmark/ repository, which
compares against some other well-known string-to-float conversion
functions. A breakdown can be found here:

0d9f020f1a/PERFORMANCE.md (commit-b15406a0d2e18b50a4b62fceb5a6a3bb60ca5706)

In summary, we are within 20% of the C++ reference implementation and
have about ~600-700MB/s throughput on a Intel I5-6500 3.5Ghz.

* F128 Support

Finally, f128 is now completely supported with full accuracy. This does
use a slower path which is possible to improve in future.

* Behavioural Changes

There are a few behavioural changes to note.

 - `parseHexFloat` is now redundant and these are now supported directly
   in `parseFloat`.
 - We implement round-to-even in all parsing routines. This is as
   specified by IEEE-754. Previous code used different rounding
   mechanisms (standard was round-to-zero, hex-parsing looked to use
   round-up) so there may be subtle differences.

Closes #2207.
Fixes #11169.
2022-05-03 16:46:40 +12:00

132 lines
5.1 KiB
Zig
Raw Blame History

This file contains ambiguous Unicode characters

This file contains Unicode characters that might be confused with other characters. If you think that this is intentional, you can safely ignore this warning. Use the Escape button to reveal them.

const std = @import("std");
const Self = @This();
// Minimum exponent that for a fast path case, or `-⌊(MANTISSA_EXPLICIT_BITS+1)/log2(5)⌋`
min_exponent_fast_path: comptime_int,
// Maximum exponent that for a fast path case, or `⌊(MANTISSA_EXPLICIT_BITS+1)/log2(5)⌋`
max_exponent_fast_path: comptime_int,
// Maximum exponent that can be represented for a disguised-fast path case.
// This is `MAX_EXPONENT_FAST_PATH + ⌊(MANTISSA_EXPLICIT_BITS+1)/log2(10)⌋`
max_exponent_fast_path_disguised: comptime_int,
// Maximum mantissa for the fast-path (`1 << 53` for f64).
max_mantissa_fast_path: comptime_int,
// Smallest decimal exponent for a non-zero value. Including subnormals.
smallest_power_of_ten: comptime_int,
// Largest decimal exponent for a non-infinite value.
largest_power_of_ten: comptime_int,
// The number of bits in the significand, *excluding* the hidden bit.
mantissa_explicit_bits: comptime_int,
// Minimum exponent value `-(1 << (EXP_BITS - 1)) + 1`.
minimum_exponent: comptime_int,
// Round-to-even only happens for negative values of q
// when q ≥ 4 in the 64-bit case and when q ≥ 17 in
// the 32-bitcase.
//
// When q ≥ 0,we have that 5^q ≤ 2m+1. In the 64-bit case,we
// have 5^q ≤ 2m+1 ≤ 2^54 or q ≤ 23. In the 32-bit case,we have
// 5^q ≤ 2m+1 ≤ 2^25 or q ≤ 10.
//
// When q < 0, we have w ≥ (2m+1)×5^q. We must have that w < 2^64
// so (2m+1)×5^q < 2^64. We have that 2m+1 > 2^53 (64-bit case)
// or 2m+1 > 2^24 (32-bit case). Hence,we must have 2^53×5^q < 2^64
// (64-bit) and 2^24×5^q < 2^64 (32-bit). Hence we have 5^q < 2^11
// or q ≥ 4 (64-bit case) and 5^q < 2^40 or q ≥ 17 (32-bitcase).
//
// Thus we have that we only need to round ties to even when
// we have that q ∈ [4,23](in the 64-bit case) or q∈[17,10]
// (in the 32-bit case). In both cases,the power of five(5^|q|)
// fits in a 64-bit word.
min_exponent_round_to_even: comptime_int,
max_exponent_round_to_even: comptime_int,
// Largest exponent value `(1 << EXP_BITS) - 1`.
infinite_power: comptime_int,
// Following should compute based on derived calculations where possible.
pub fn from(comptime T: type) Self {
return switch (T) {
f16 => .{
// Fast-Path
.min_exponent_fast_path = -4,
.max_exponent_fast_path = 4,
.max_exponent_fast_path_disguised = 7,
.max_mantissa_fast_path = 2 << std.math.floatMantissaBits(T),
// Slow + Eisel-Lemire
.mantissa_explicit_bits = std.math.floatMantissaBits(T),
.infinite_power = 0x1f,
// Eisel-Lemire
.smallest_power_of_ten = -26, // TODO: refine, fails one test
.largest_power_of_ten = 4,
.minimum_exponent = -15,
// w >= (2m+1) * 5^-q and w < 2^64
// => 2m+1 > 2^11
// => 2^11*5^-q < 2^64
// => 5^-q < 2^53
// => q >= -23
.min_exponent_round_to_even = -22,
.max_exponent_round_to_even = 5,
},
f32 => .{
// Fast-Path
.min_exponent_fast_path = -10,
.max_exponent_fast_path = 10,
.max_exponent_fast_path_disguised = 17,
.max_mantissa_fast_path = 2 << std.math.floatMantissaBits(T),
// Slow + Eisel-Lemire
.mantissa_explicit_bits = std.math.floatMantissaBits(T),
.infinite_power = 0xff,
// Eisel-Lemire
.smallest_power_of_ten = -65,
.largest_power_of_ten = 38,
.minimum_exponent = -127,
.min_exponent_round_to_even = -17,
.max_exponent_round_to_even = 10,
},
f64 => .{
// Fast-Path
.min_exponent_fast_path = -22,
.max_exponent_fast_path = 22,
.max_exponent_fast_path_disguised = 37,
.max_mantissa_fast_path = 2 << std.math.floatMantissaBits(T),
// Slow + Eisel-Lemire
.mantissa_explicit_bits = std.math.floatMantissaBits(T),
.infinite_power = 0x7ff,
// Eisel-Lemire
.smallest_power_of_ten = -342,
.largest_power_of_ten = 308,
.minimum_exponent = -1023,
.min_exponent_round_to_even = -4,
.max_exponent_round_to_even = 23,
},
f128 => .{
// Fast-Path
.min_exponent_fast_path = -48,
.max_exponent_fast_path = 48,
.max_exponent_fast_path_disguised = 82,
.max_mantissa_fast_path = 2 << std.math.floatMantissaBits(T),
// Slow + Eisel-Lemire
.mantissa_explicit_bits = std.math.floatMantissaBits(T),
.infinite_power = 0x7fff,
// Eisel-Lemire.
// NOTE: Not yet tested (no f128 eisel-lemire implementation)
.smallest_power_of_ten = -4966,
.largest_power_of_ten = 4932,
.minimum_exponent = -16382,
// 2^113 * 5^-q < 2^128
// 5^-q < 2^15
// => q >= -6
.min_exponent_round_to_even = -6,
.max_exponent_round_to_even = 49,
},
else => unreachable,
};
}