tests/tcg: add float_madds test to multiarch

This is a generic floating point multiply and accumulate test for
single precision floating point values. I've split of the common float
functions into a helper library so additional tests can use the same
common code.

As I don't have references for all architectures I've allowed some
flexibility for tests to pass without reference files. They can be
added as we get collect them.

Signed-off-by: Alex Bennée <alex.bennee@linaro.org>
Reviewed-by: Richard Henderson <richard.henderson@linaro.org>
Tested-by: Philippe Mathieu-Daudé <philmd@redhat.com>

4 files changed, 370 insertions, 1 deletions

diff --git a/tests/tcg/multiarch/Makefile.target b/tests/tcg/multiarch/Makefile.target
index 657a04f..30ecb8b 100644
--- a/tests/tcg/multiarch/Makefile.target
+++ b/tests/tcg/multiarch/Makefile.target
@@ -10,12 +10,22 @@ MULTIARCH_SRC=$(SRC_PATH)/tests/tcg/multiarch
 # Set search path for all sources
 VPATH 		+= $(MULTIARCH_SRC)
 MULTIARCH_SRCS   =$(notdir $(wildcard $(MULTIARCH_SRC)/*.c))
-MULTIARCH_TESTS  =$(MULTIARCH_SRCS:.c=)
+MULTIARCH_TESTS  =$(filter-out float_helpers, $(MULTIARCH_SRCS:.c=))
 
 #
 # The following are any additional rules needed to build things
 #
 
+
+float_madds: LDFLAGS+=-lm
+float_madds: float_madds.c float_helpers.c
+	$(CC) $(CFLAGS) $(EXTRA_CFLAGS) $< $(MULTIARCH_SRC)/float_helpers.c -o $@ $(LDFLAGS)
+
+run-float_madds: float_madds
+	$(call run-test,$<, $(QEMU) $(QEMU_OPTS) $<,"$< on $(TARGET_NAME)")
+	$(call conditional-diff-out,$<,$(SRC_PATH)/tests/tcg/$(TARGET_NAME)/$<.ref)
+
+
 testthread: LDFLAGS+=-lpthread
 
 # We define the runner for test-mmap after the individual
diff --git a/tests/tcg/multiarch/float_helpers.c b/tests/tcg/multiarch/float_helpers.c
new file mode 100644
index 0000000..8ee7903
--- /dev/null
+++ b/tests/tcg/multiarch/float_helpers.c
@@ -0,0 +1,230 @@
+/*
+ * Common Float Helpers
+ *
+ * This contains a series of useful utility routines and a set of
+ * floating point constants useful for exercising the edge cases in
+ * floating point tests.
+ *
+ * Copyright (c) 2019 Linaro
+ *
+ * SPDX-License-Identifier: GPL-3.0-or-later
+ */
+
+/* we want additional float type definitions */
+#define __STDC_WANT_IEC_60559_BFP_EXT__
+#define __STDC_WANT_IEC_60559_TYPES_EXT__
+
+#define _GNU_SOURCE
+#include <stdio.h>
+#include <stdlib.h>
+#include <inttypes.h>
+#include <math.h>
+#include <float.h>
+#include <fenv.h>
+
+#include "float_helpers.h"
+
+#define ARRAY_SIZE(x) (sizeof(x) / sizeof((x)[0]))
+
+/*
+ * Half Precision Numbers
+ *
+ * Not yet well standardised so we return a plain uint16_t for now.
+ */
+
+/* no handy defines for these numbers */
+static uint16_t f16_numbers[] = {
+    0xffff, /* -NaN / AHP -Max */
+    0xfcff, /* -NaN / AHP */
+    0xfc01, /* -NaN / AHP */
+    0xfc00, /* -Inf */
+    0xfbff, /* -Max */
+    0xc000, /* -2 */
+    0xbc00, /* -1 */
+    0x8001, /* -MIN subnormal */
+    0x8000, /* -0 */
+    0x0000, /* +0 */
+    0x0001, /* MIN subnormal */
+    0x3c00, /* 1 */
+    0x7bff, /* Max */
+    0x7c00, /* Inf */
+    0x7c01, /* NaN / AHP */
+    0x7cff, /* NaN / AHP */
+    0x7fff, /* NaN / AHP +Max*/
+};
+
+static const int num_f16 = ARRAY_SIZE(f16_numbers);
+
+int get_num_f16(void)
+{
+    return num_f16;
+}
+
+uint16_t get_f16(int i)
+{
+    return f16_numbers[i % num_f16];
+}
+
+/* only display as hex */
+char *fmt_16(uint16_t num)
+{
+    char *fmt;
+    asprintf(&fmt, "f16(%#04x)", num);
+    return fmt;
+}
+
+/*
+ * Single Precision Numbers
+ */
+
+#ifndef SNANF
+/* Signaling NaN macros, if supported.  */
+# if __GNUC_PREREQ(3, 3)
+#  define SNANF (__builtin_nansf (""))
+#  define SNAN (__builtin_nans (""))
+#  define SNANL (__builtin_nansl (""))
+# endif
+#endif
+
+static float f32_numbers[] = {
+    -SNANF,
+    -NAN,
+    -INFINITY,
+    -FLT_MAX,
+    -0x1.1874b2p+103,
+    -0x1.c0bab6p+99,
+    -0x1.31f75p-40,
+    -0x1.505444p-66,
+    -FLT_MIN,
+    0.0,
+    FLT_MIN,
+    0x1p-25,
+    0x1.ffffe6p-25, /* min positive FP16 subnormal */
+    0x1.ff801ap-15, /* max subnormal FP16 */
+    0x1.00000cp-14, /* min positive normal FP16 */
+    1.0,
+    0x1.004p+0, /* smallest float after 1.0 FP16 */
+    2.0,
+    M_E, M_PI,
+    0x1.ffbep+15,
+    0x1.ffcp+15, /* max FP16 */
+    0x1.ffc2p+15,
+    0x1.ffbfp+16,
+    0x1.ffcp+16, /* max AFP */
+    0x1.ffc1p+16,
+    0x1.c0bab6p+99,
+    FLT_MAX,
+    INFINITY,
+    NAN,
+    SNANF
+};
+
+static const int num_f32 = ARRAY_SIZE(f32_numbers);
+
+int get_num_f32(void)
+{
+    return num_f32;
+}
+
+float get_f32(int i)
+{
+    return f32_numbers[i % num_f32];
+}
+
+char *fmt_f32(float num)
+{
+    uint32_t single_as_hex = *(uint32_t *) &num;
+    char *fmt;
+    asprintf(&fmt, "f32(%02.20a:%#010x)", num, single_as_hex);
+    return fmt;
+}
+
+
+/* This allows us to initialise some doubles as pure hex */
+typedef union {
+    double d;
+    uint64_t h;
+} test_doubles;
+
+static test_doubles f64_numbers[] = {
+    {SNAN},
+    {-NAN},
+    {-INFINITY},
+    {-DBL_MAX},
+    {-FLT_MAX-1.0},
+    {-FLT_MAX},
+    {-1.111E+31},
+    {-1.111E+30}, /* half prec */
+    {-2.0}, {-1.0},
+    {-DBL_MIN},
+    {-FLT_MIN},
+    {0.0},
+    {FLT_MIN},
+    {2.98023224e-08},
+    {5.96046E-8}, /* min positive FP16 subnormal */
+    {6.09756E-5}, /* max subnormal FP16 */
+    {6.10352E-5}, /* min positive normal FP16 */
+    {1.0},
+    {1.0009765625}, /* smallest float after 1.0 FP16 */
+    {DBL_MIN},
+    {1.3789972848607228e-308},
+    {1.4914738736681624e-308},
+    {1.0}, {2.0},
+    {M_E}, {M_PI},
+    {65503.0},
+    {65504.0}, /* max FP16 */
+    {65505.0},
+    {131007.0},
+    {131008.0}, /* max AFP */
+    {131009.0},
+    {.h = 0x41dfffffffc00000 }, /* to int = 0x7fffffff */
+    {FLT_MAX},
+    {FLT_MAX + 1.0},
+    {DBL_MAX},
+    {INFINITY},
+    {NAN},
+    {.h = 0x7ff0000000000001}, /* SNAN */
+    {SNAN},
+};
+
+static const int num_f64 = ARRAY_SIZE(f64_numbers);
+
+int get_num_f64(void)
+{
+    return num_f64;
+}
+
+double get_f64(int i)
+{
+    return f64_numbers[i % num_f64].d;
+}
+
+char *fmt_f64(double num)
+{
+    uint64_t double_as_hex = *(uint64_t *) &num;
+    char *fmt;
+    asprintf(&fmt, "f64(%02.20a:%#020" PRIx64 ")", num, double_as_hex);
+    return fmt;
+}
+
+/*
+ * Float flags
+ */
+char *fmt_flags(void)
+{
+    int flags = fetestexcept(FE_ALL_EXCEPT);
+    char *fmt;
+
+    if (flags) {
+        asprintf(&fmt, "%s%s%s%s%s",
+                 flags & FE_OVERFLOW ? "OVERFLOW " : "",
+                 flags & FE_UNDERFLOW ? "UNDERFLOW " : "",
+                 flags & FE_DIVBYZERO ? "DIV0 " : "",
+                 flags & FE_INEXACT ? "INEXACT " : "",
+                 flags & FE_INVALID ? "INVALID" : "");
+    } else {
+        asprintf(&fmt, "OK");
+    }
+
+    return fmt;
+}
diff --git a/tests/tcg/multiarch/float_helpers.h b/tests/tcg/multiarch/float_helpers.h
new file mode 100644
index 0000000..6337bc6
--- /dev/null
+++ b/tests/tcg/multiarch/float_helpers.h
@@ -0,0 +1,26 @@
+/*
+ * Common Float Helpers
+ *
+ * Copyright (c) 2019 Linaro
+ *
+ * SPDX-License-Identifier: GPL-3.0-or-later
+ */
+
+#include <inttypes.h>
+
+/* Number of constants in each table */
+int get_num_f16(void);
+int get_num_f32(void);
+int get_num_f64(void);
+
+/* Accessor helpers, overflows will automatically wrap */
+uint16_t get_f16(int i); /* use _Float16 when we can */
+float    get_f32(int i);
+double   get_f64(int i);
+
+/* Return format strings, free after use */
+char * fmt_f16(uint16_t);
+char * fmt_f32(float);
+char * fmt_f64(double);
+/* exception flags */
+char * fmt_flags(void);
diff --git a/tests/tcg/multiarch/float_madds.c b/tests/tcg/multiarch/float_madds.c
new file mode 100644
index 0000000..eceb4ae
--- /dev/null
+++ b/tests/tcg/multiarch/float_madds.c
@@ -0,0 +1,103 @@
+/*
+ * Fused Multiply Add (Single)
+ *
+ * Copyright (c) 2019 Linaro
+ *
+ * SPDX-License-Identifier: GPL-3.0-or-later
+ */
+
+#include <stdio.h>
+#include <stdlib.h>
+#include <math.h>
+#include <float.h>
+#include <fenv.h>
+
+#include "float_helpers.h"
+
+#define ARRAY_SIZE(x) (sizeof(x) / sizeof((x)[0]))
+
+typedef struct {
+    int flag;
+    char *desc;
+} float_mapping;
+
+float_mapping round_flags[] = {
+    { FE_TONEAREST, "to nearest" },
+#ifdef FE_UPWARD
+    { FE_UPWARD, "upwards" },
+#endif
+#ifdef FE_DOWNWARD
+    { FE_DOWNWARD, "downwards" },
+#endif
+    { FE_TOWARDZERO, "to zero" }
+};
+
+
+static void print_inputs(float a, float b, float c)
+{
+    char *a_fmt, *b_fmt, *c_fmt;
+
+    a_fmt = fmt_f32(a);
+    b_fmt = fmt_f32(b);
+    c_fmt = fmt_f32(c);
+
+    printf("op : %s * %s + %s\n", a_fmt, b_fmt, c_fmt);
+
+    free(a_fmt);
+    free(b_fmt);
+    free(c_fmt);
+}
+
+static void print_result(float r, int j, int k)
+{
+    char *r_fmt, *flag_fmt;
+
+    r_fmt = fmt_f32(r);
+    flag_fmt = fmt_flags();
+
+    printf("res: %s flags=%s (%d/%d)\n", r_fmt, flag_fmt, j, k);
+
+    free(r_fmt);
+    free(flag_fmt);
+}
+
+static void do_madds(float a, float b, float c, int j, int k)
+{
+    float r;
+
+    print_inputs(a, b, c);
+
+    feclearexcept(FE_ALL_EXCEPT);
+    r = __builtin_fmaf(a, b, c);
+
+    print_result(r, j, k);
+}
+
+int main(int argc, char *argv[argc])
+{
+    int i, j, k, nums = get_num_f32();
+    float a, b, c;
+
+    for (i = 0; i < ARRAY_SIZE(round_flags); ++i) {
+        if (fesetround(round_flags[i].flag) != 0) {
+            printf("### Rounding %s skipped\n", round_flags[i].desc);
+            continue;
+        }
+        printf("### Rounding %s\n", round_flags[i].desc);
+        for (j = 0; j < nums; j++) {
+            for (k = 0; k < 3; k++) {
+                a = get_f32(j + ((k)%3));
+                b = get_f32(j + ((k+1)%3));
+                c = get_f32(j + ((k+2)%3));
+                do_madds(a, b, c, j, k);
+            }
+        }
+
+        /* From https://bugs.launchpad.net/qemu/+bug/1841491 */
+        printf("# LP184149\n");
+        do_madds(0x1.ffffffffffffcp-1022, 0x1.0000000000001p-1, 0x0.0000000000001p-1022, j, 0);
+        do_madds(0x8p-152, 0x8p-152, 0x8p-152, j+1, 0);
+    }
+
+    return 0;
+}