py/formatfloat: Improve accuracy of float formatting code.

Following discussions in PR #16666, this commit updates the float
formatting code to improve the `repr` reversibility, i.e. the percentage of
valid floating point numbers that do parse back to the same number when
formatted by `repr` (in CPython it's 100%).

This new code offers a choice of 3 float conversion methods, depending on
the desired tradeoff between code size and conversion precision:

- BASIC method is the smallest code footprint

- APPROX method uses an iterative method to approximate the exact
  representation, which is a bit slower but but does not have a big impact
  on code size.  It provides `repr` reversibility on >99.8% of the cases in
  double precision, and on >98.5% in single precision (except with REPR_C,
  where reversibility is 100% as the last two bits are not taken into
  account).

- EXACT method uses higher-precision floats during conversion, which
  provides perfect results but has a higher impact on code size.  It is
  faster than APPROX method, and faster than the CPython equivalent
  implementation.  It is however not available on all compilers when using
  FLOAT_IMPL_DOUBLE.

Here is the table comparing the impact of the three conversion methods on
code footprint on PYBV10 (using single-precision floats) and reversibility
rate for both single-precision and double-precision floats.  The table
includes current situation as a baseline for the comparison:

              PYBV10  REPR_C   FLOAT  DOUBLE
    current = 364688   12.9%   27.6%   37.9%
    basic   = 364812   85.6%   60.5%   85.7%
    approx  = 365080  100.0%   98.5%   99.8%
    exact   = 366408  100.0%  100.0%  100.0%

Signed-off-by: Yoctopuce dev <dev@yoctopuce.com>

tests/float/float_format.py

@@ -2,14 +2,25 @@
 
 # general rounding
 for val in (116, 1111, 1234, 5010, 11111):
-    print("%.0f" % val)
-    print("%.1f" % val)
-    print("%.3f" % val)
+    print("Test on %d / 1000:" % val)
+    for fmt in ("%.5e", "%.3e", "%.1e", "%.0e", "%.3f", "%.1f", "%.0f", "%.3g", "%.1g", "%.0g"):
+        print(fmt, fmt % (val / 1000))
+
+# make sure round-up to the next unit is handled properly
+for val in range(4, 9):
+    divi = 10**val
+    print("Test on 99994 / (10 ** %d):" % val)
+    for fmt in ("%.5e", "%.3e", "%.1e", "%.0e", "%.3f", "%.1f", "%.0f", "%.3g", "%.1g", "%.0g"):
+        print(fmt, fmt % (99994 / divi))
 
 # make sure rounding is done at the correct precision
 for prec in range(8):
     print(("%%.%df" % prec) % 6e-5)
 
+# make sure trailing zeroes are added properly
+for prec in range(8):
+    print(("%%.%df" % prec) % 1e19)
+
 # check certain cases that had a digit value of 10 render as a ":" character
 print("%.2e" % float("9" * 51 + "e-39"))
 print("%.2e" % float("9" * 40 + "e-21"))

tests/float/float_format_accuracy.py

@@ -0,0 +1,73 @@
+# Test accuracy of `repr` conversions.
+# This test also increases code coverage for corner cases.
+
+try:
+    import array, math, random
+except ImportError:
+    print("SKIP")
+    raise SystemExit
+
+# The largest errors come from seldom used very small numbers, near the
+# limit of the representation. So we keep them out of this test to keep
+# the max relative error display useful.
+if float("1e-100") == 0.0:
+    # single-precision
+    float_type = "f"
+    float_size = 4
+    # testing range
+    min_expo = -96  # i.e. not smaller than 1.0e-29
+    # Expected results (given >=50'000 samples):
+    # - MICROPY_FLTCONV_IMPL_EXACT: 100% exact conversions
+    # - MICROPY_FLTCONV_IMPL_APPROX: >=98.53% exact conversions, max relative error <= 1.01e-7
+    min_success = 0.980  # with only 1200 samples, the success rate is lower
+    max_rel_err = 1.1e-7
+    # REPR_C is typically used with FORMAT_IMPL_BASIC, which has a larger error
+    is_REPR_C = float("1.0000001") == float("1.0")
+    if is_REPR_C:  # REPR_C
+        min_success = 0.83
+        max_rel_err = 5.75e-07
+else:
+    # double-precision
+    float_type = "d"
+    float_size = 8
+    # testing range
+    min_expo = -845  # i.e. not smaller than 1.0e-254
+    # Expected results (given >=200'000 samples):
+    # - MICROPY_FLTCONV_IMPL_EXACT: 100% exact conversions
+    # - MICROPY_FLTCONV_IMPL_APPROX: >=99.83% exact conversions, max relative error <= 2.7e-16
+    min_success = 0.997  # with only 1200 samples, the success rate is lower
+    max_rel_err = 2.7e-16
+
+
+# Deterministic pseudorandom generator. Designed to be uniform
+# on mantissa values and exponents, not on the represented number
+def pseudo_randfloat():
+    rnd_buff = bytearray(float_size)
+    for _ in range(float_size):
+        rnd_buff[_] = random.getrandbits(8)
+    return array.array(float_type, rnd_buff)[0]
+
+
+random.seed(42)
+stats = 0
+N = 1200
+max_err = 0
+for _ in range(N):
+    f = pseudo_randfloat()
+    while type(f) is not float or math.isinf(f) or math.isnan(f) or math.frexp(f)[1] <= min_expo:
+        f = pseudo_randfloat()
+
+    str_f = repr(f)
+    f2 = float(str_f)
+    if f2 == f:
+        stats += 1
+    else:
+        error = abs((f2 - f) / f)
+        if max_err < error:
+            max_err = error
+
+print(N, "values converted")
+if stats / N >= min_success and max_err <= max_rel_err:
+    print("float format accuracy OK")
+else:
+    print("FAILED: repr rate=%.3f%% max_err=%.3e" % (100 * stats / N, max_err))

tests/float/float_format_ints.py

@@ -12,14 +12,42 @@ for b in [13, 123, 457, 23456]:
             print(title, "with format", f_fmt, "gives", f_fmt.format(f))
             print(title, "with format", g_fmt, "gives", g_fmt.format(f))
 
+# The tests below check border cases involving all mantissa bits.
+# In case of REPR_C, where the mantissa is missing two bits, the
+# the string representation for such numbers might not always be exactly
+# the same but nevertheless be correct, so we must allow a few exceptions.
+is_REPR_C = float("1.0000001") == float("1.0")
+
 # 16777215 is 2^24 - 1, the largest integer that can be completely held
 # in a float32.
-print("{:f}".format(16777215))
+val_str = "{:f}".format(16777215)
+
+# When using REPR_C, 16777215.0 is the same as 16777212.0 or 16777214.4
+# (depending on the implementation of pow() function, the result may differ)
+if is_REPR_C and (val_str == "16777212.000000" or val_str == "16777214.400000"):
+    val_str = "16777215.000000"
+
+print(val_str)
+
 # 4294967040 = 16777215 * 128 is the largest integer that is exactly
 # represented by a float32 and that will also fit within a (signed) int32.
 # The upper bound of our integer-handling code is actually double this,
 # but that constant might cause trouble on systems using 32 bit ints.
-print("{:f}".format(2147483520))
+val_str = "{:f}".format(2147483520)
+
+# When using FLOAT_IMPL_FLOAT, 2147483520.0 == 2147483500.0
+# Both representations are valid, the second being "simpler"
+is_float32 = float("1e300") == float("inf")
+if is_float32 and val_str == "2147483500.000000":
+    val_str = "2147483520.000000"
+
+# When using REPR_C, 2147483520.0 is the same as 2147483200.0
+# Both representations are valid, the second being "simpler"
+if is_REPR_C and val_str == "2147483200.000000":
+    val_str = "2147483520.000000"
+
+print(val_str)
+
 # Very large positive integers can be a test for precision and resolution.
 # This is a weird way to represent 1e38 (largest power of 10 for float32).
 print("{:.6e}".format(float("9" * 30 + "e8")))

tests/float/float_struct_e.py

@@ -32,7 +32,7 @@ for j in test_values:
     for i in (j, -j):
         x = struct.pack("<e", i)
         v = struct.unpack("<e", x)[0]
-        print("%.7f %s %.15f %s" % (i, x, v, i == v))
+        print("%.7f %s %.7f %s" % (i, x, v, i == v))
 
 # In CPython, packing a float that doesn't fit into a half-float raises OverflowError.
 # But in MicroPython it does not, but rather stores the value as inf.

tests/float/float_struct_e_doubleprec.py

@@ -0,0 +1,43 @@
+# Test struct pack/unpack with 'e' typecode.
+
+try:
+    import struct
+except ImportError:
+    print("SKIP")
+    raise SystemExit
+
+test_values = (
+    1e-7,
+    2e-7,
+    1e-6,
+    1e-5,
+    1e-4,
+    1e-3,
+    1e-2,
+    0.1,
+    0,
+    1,
+    2,
+    4,
+    8,
+    10,
+    100,
+    1e3,
+    1e4,
+    6e4,
+    float("inf"),
+)
+
+for j in test_values:
+    for i in (j, -j):
+        x = struct.pack("<e", i)
+        v = struct.unpack("<e", x)[0]
+        print("%.7f %s %.15f %s" % (i, x, v, i == v))
+
+# In CPython, packing a float that doesn't fit into a half-float raises OverflowError.
+# But in MicroPython it does not, but rather stores the value as inf.
+# This test is here for coverage.
+try:
+    struct.pack("e", 1e15)
+except OverflowError:
+    pass

tests/float/float_struct_e_fp30.py

@@ -0,0 +1,43 @@
+# Test struct pack/unpack with 'e' typecode.
+
+try:
+    import struct
+except ImportError:
+    print("SKIP")
+    raise SystemExit
+
+test_values = (
+    1e-7,
+    2e-7,
+    1e-6,
+    1e-5,
+    1e-4,
+    1e-3,
+    1e-2,
+    0.1,
+    0,
+    1,
+    2,
+    4,
+    8,
+    10,
+    100,
+    1e3,
+    1e4,
+    6e4,
+    float("inf"),
+)
+
+for j in test_values:
+    for i in (j, -j):
+        x = struct.pack("<e", i)
+        v = struct.unpack("<e", x)[0]
+        print("%.7f %s %.5f %s" % (i, x, v, i == v))
+
+# In CPython, packing a float that doesn't fit into a half-float raises OverflowError.
+# But in MicroPython it does not, but rather stores the value as inf.
+# This test is here for coverage.
+try:
+    struct.pack("e", 1e15)
+except OverflowError:
+    pass

tests/float/string_format_modulo3.py

@@ -1,3 +1,3 @@
-# uPy and CPython outputs differ for the following
+# Test corner cases where MicroPython and CPython outputs used to differ in the past
 print("%.1g" % -9.9)  # round up 'g' with '-' sign
 print("%.2g" % 99.9)  # round up

tests/float/string_format_modulo3.py.exp

@@ -1,2 +0,0 @@
--10
-100

tests/ports/unix/extra_coverage.py.exp

@@ -127,10 +127,6 @@ TypeError: can't convert NoneType to int
 TypeError: can't convert NoneType to int
 ValueError: 
 Warning: test
-# format float
-?
-+1e+00
-+1e+00
 # binary
 123
 456

tests/ports/unix/ffi_float2.py

@@ -29,4 +29,5 @@ except OSError:
 
 for fun in (tgammaf,):
     for val in (0.5, 1, 1.0, 1.5, 4, 4.0):
-        print("%.6f" % fun(val))
+        # limit to 5 decimals in order to pass with REPR_C with FORMAT_IMPL_BASIC
+        print("%.5f" % fun(val))

tests/ports/unix/ffi_float2.py.exp

@@ -1,6 +1,6 @@
-1.772454
-1.000000
-1.000000
-0.886227
-6.000000
-6.000000
+1.77245
+1.00000
+1.00000
+0.88623
+6.00000
+6.00000

tests/run-tests.py

@@ -785,12 +785,16 @@ def run_tests(pyb, tests, args, result_dir, num_threads=1):
         skip_tests.add(
             "float/float2int_intbig.py"
         )  # requires fp32, there's float2int_fp30_intbig.py instead
+        skip_tests.add(
+            "float/float_struct_e.py"
+        )  # requires fp32, there's float_struct_e_fp30.py instead
         skip_tests.add("float/bytes_construct.py")  # requires fp32
         skip_tests.add("float/bytearray_construct.py")  # requires fp32
         skip_tests.add("float/float_format_ints_power10.py")  # requires fp32
     if upy_float_precision < 64:
         skip_tests.add("float/float_divmod.py")  # tested by float/float_divmod_relaxed.py instead
         skip_tests.add("float/float2int_doubleprec_intbig.py")
+        skip_tests.add("float/float_struct_e_doubleprec.py")
         skip_tests.add("float/float_format_ints_doubleprec.py")
         skip_tests.add("float/float_parse_doubleprec.py")