lib: Add halffloat implementation

Probe for and make an API available for tests to use f16c intrinsics to
generate needed fp16 pixel data.

Also import a pure c fp32 <-> fp16 conversion implementation from Mesa
18.3.4, which will act as a fallback when f16c is unavailable.

rfc2:
- Change API to reduce number of function calls (Maarten)

v1:
- Move pragma so AVX code isn't emitted for fallbacks (Ville)
- Change edx to ecx (Ville)

Signed-off-by: Kevin Strasser <kevin.strasser@intel.com>
Reviewed-by: Ville Syrjälä <ville.syrjala@linux.intel.com>
Reviewed-by: Maarten Lankhorst <maarten.lankhorst@linux.intel.com>
Signed-off-by: Ville Syrjälä <ville.syrjala@linux.intel.com>

6 files changed, 274 insertions, 0 deletions

diff --git a/lib/Makefile.sources b/lib/Makefile.sources
index 62219349..97685823 100644
--- a/lib/Makefile.sources
+++ b/lib/Makefile.sources
@@ -33,6 +33,8 @@ lib_source_list =	 	\
 	igt_gt.h		\
 	igt_gvt.c		\
 	igt_gvt.h		\
+	igt_halffloat.c		\
+	igt_halffloat.h		\
 	igt_matrix.c		\
 	igt_matrix.h		\
 	igt_primes.c		\
diff --git a/lib/igt_halffloat.c b/lib/igt_halffloat.c
new file mode 100644
index 00000000..263b28c1
--- /dev/null
+++ b/lib/igt_halffloat.c
@@ -0,0 +1,234 @@
+/*
+ * Copyright (C) 1999-2007  Brian Paul   All Rights Reserved.
+ * Copyright 2015 Philip Taylor <philip@zaynar.co.uk>
+ * Copyright 2018 Advanced Micro Devices, Inc.
+ *
+ * Permission is hereby granted, free of charge, to any person obtaining a
+ * copy of this software and associated documentation files (the "Software"),
+ * to deal in the Software without restriction, including without limitation
+ * the rights to use, copy, modify, merge, publish, distribute, sublicense,
+ * and/or sell copies of the Software, and to permit persons to whom the
+ * Software is furnished to do so, subject to the following conditions:
+ *
+ * The above copyright notice and this permission notice shall be included
+ * in all copies or substantial portions of the Software.
+ *
+ * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
+ * OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
+ * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT.  IN NO EVENT SHALL
+ * THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR
+ * OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE,
+ * ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR
+ * OTHER DEALINGS IN THE SOFTWARE.
+ */
+
+#include <assert.h>
+#include <math.h>
+
+#include "igt_halffloat.h"
+#include "igt_x86.h"
+
+typedef union { float f; int32_t i; uint32_t u; } fi_type;
+
+/**
+ * Convert a 4-byte float to a 2-byte half float.
+ *
+ * Not all float32 values can be represented exactly as a float16 value. We
+ * round such intermediate float32 values to the nearest float16. When the
+ * float32 lies exactly between to float16 values, we round to the one with
+ * an even mantissa.
+ *
+ * This rounding behavior has several benefits:
+ *   - It has no sign bias.
+ *
+ *   - It reproduces the behavior of real hardware: opcode F32TO16 in Intel's
+ *     GPU ISA.
+ *
+ *   - By reproducing the behavior of the GPU (at least on Intel hardware),
+ *     compile-time evaluation of constant packHalf2x16 GLSL expressions will
+ *     result in the same value as if the expression were executed on the GPU.
+ */
+static inline uint16_t _float_to_half(float val)
+{
+	const fi_type fi = {val};
+	const int flt_m = fi.i & 0x7fffff;
+	const int flt_e = (fi.i >> 23) & 0xff;
+	const int flt_s = (fi.i >> 31) & 0x1;
+	int s, e, m = 0;
+	uint16_t result;
+
+	/* sign bit */
+	s = flt_s;
+
+	/* handle special cases */
+	if ((flt_e == 0) && (flt_m == 0)) {
+		/* zero */
+		/* m = 0; - already set */
+		e = 0;
+	} else if ((flt_e == 0) && (flt_m != 0)) {
+		/* denorm -- denorm float maps to 0 half */
+		/* m = 0; - already set */
+		e = 0;
+	} else if ((flt_e == 0xff) && (flt_m == 0)) {
+		/* infinity */
+		/* m = 0; - already set */
+		e = 31;
+	} else if ((flt_e == 0xff) && (flt_m != 0)) {
+		/* NaN */
+		m = 1;
+		e = 31;
+	} else {
+		/* regular number */
+		const int new_exp = flt_e - 127;
+		if (new_exp < -14) {
+			/* The float32 lies in the range (0.0, min_normal16) and
+			 * is rounded to a nearby float16 value. The result will
+			 * be either zero, subnormal, or normal.
+			 */
+			e = 0;
+			m = lrintf((1 << 24) * fabsf(fi.f));
+		} else if (new_exp > 15) {
+			/* map this value to infinity */
+			/* m = 0; - already set */
+			e = 31;
+		} else {
+			/* The float32 lies in the range
+			 *   [min_normal16, max_normal16 + max_step16)
+			 * and is rounded to a nearby float16 value. The result
+			 * will be either normal or infinite.
+			 */
+			e = new_exp + 15;
+			m = lrintf(flt_m / (float)(1 << 13));
+		}
+	}
+
+	assert(0 <= m && m <= 1024);
+	if (m == 1024) {
+		/* The float32 was rounded upwards into the range of the next
+		 * exponent, so bump the exponent. This correctly handles the
+		 * case where f32 should be rounded up to float16 infinity.
+		 */
+		++e;
+		m = 0;
+	}
+
+	result = (s << 15) | (e << 10) | m;
+	return result;
+}
+
+/**
+ * Convert a 2-byte half float to a 4-byte float.
+ * Based on code from:
+ * http://www.opengl.org/discussion_boards/ubb/Forum3/HTML/008786.html
+ */
+static inline float _half_to_float(uint16_t val)
+{
+	/* XXX could also use a 64K-entry lookup table */
+	const int m = val & 0x3ff;
+	const int e = (val >> 10) & 0x1f;
+	const int s = (val >> 15) & 0x1;
+	int flt_m, flt_e, flt_s;
+	fi_type fi;
+
+	/* sign bit */
+	flt_s = s;
+
+	/* handle special cases */
+	if ((e == 0) && (m == 0)) {
+		/* zero */
+		flt_m = 0;
+		flt_e = 0;
+	} else if ((e == 0) && (m != 0)) {
+		/* denorm -- denorm half will fit in non-denorm single */
+		const float half_denorm = 1.0f / 16384.0f; /* 2^-14 */
+		float mantissa = ((float) (m)) / 1024.0f;
+		float sign = s ? -1.0f : 1.0f;
+		return sign * mantissa * half_denorm;
+	} else if ((e == 31) && (m == 0)) {
+		/* infinity */
+		flt_e = 0xff;
+		flt_m = 0;
+	} else if ((e == 31) && (m != 0)) {
+		/* NaN */
+		flt_e = 0xff;
+		flt_m = 1;
+	} else {
+		/* regular */
+		flt_e = e + 112;
+		flt_m = m << 13;
+	}
+
+	fi.i = (flt_s << 31) | (flt_e << 23) | flt_m;
+	return fi.f;
+}
+
+#if defined(__x86_64__) && !defined(__clang__)
+#pragma GCC push_options
+#pragma GCC target("f16c")
+
+#include <immintrin.h>
+
+static void float_to_half_f16c(const float *f, uint16_t *h, unsigned int num)
+{
+	for (int i = 0; i < num; i++)
+		h[i] = _cvtss_sh(f[i], 0);
+}
+
+static void half_to_float_f16c(const uint16_t *h, float *f, unsigned int num)
+{
+	for (int i = 0; i < num; i++)
+		f[i] = _cvtsh_ss(h[i]);
+}
+
+#pragma GCC pop_options
+
+static void float_to_half(const float *f, uint16_t *h, unsigned int num)
+{
+	for (int i = 0; i < num; i++)
+		h[i] = _float_to_half(f[i]);
+}
+
+static void half_to_float(const uint16_t *h, float *f, unsigned int num)
+{
+	for (int i = 0; i < num; i++)
+		f[i] = _half_to_float(h[i]);
+}
+
+static void (*resolve_float_to_half(void))(const float *f, uint16_t *h, unsigned int num)
+{
+	if (igt_x86_features() & F16C)
+		return float_to_half_f16c;
+
+	return float_to_half;
+}
+
+void igt_float_to_half(const float *f, uint16_t *h, unsigned int num)
+	__attribute__((ifunc("resolve_float_to_half")));
+
+static void (*resolve_half_to_float(void))(const uint16_t *h, float *f, unsigned int num)
+{
+	if (igt_x86_features() & F16C)
+		return half_to_float_f16c;
+
+	return half_to_float;
+}
+
+void igt_half_to_float(const uint16_t *h, float *f, unsigned int num)
+	__attribute__((ifunc("resolve_half_to_float")));
+
+#else
+
+void igt_float_to_half(const float *f, uint16_t *h, unsigned int num)
+{
+	for (int i = 0; i < num; i++)
+		h[i] = _float_to_half(f[i]);
+}
+
+float igt_half_to_float(const uint16_t *h, float *f, unsigned int num)
+{
+	for (int i = 0; i < num; i++)
+		f[i] = _half_to_float(h[i]);
+}
+
+#endif
+
diff --git a/lib/igt_halffloat.h b/lib/igt_halffloat.h
new file mode 100644
index 00000000..8db448a1
--- /dev/null
+++ b/lib/igt_halffloat.h
@@ -0,0 +1,27 @@
+/*
+ * Copyright (C) 1999-2007  Brian Paul   All Rights Reserved.
+ *
+ * Permission is hereby granted, free of charge, to any person obtaining a
+ * copy of this software and associated documentation files (the "Software"),
+ * to deal in the Software without restriction, including without limitation
+ * the rights to use, copy, modify, merge, publish, distribute, sublicense,
+ * and/or sell copies of the Software, and to permit persons to whom the
+ * Software is furnished to do so, subject to the following conditions:
+ *
+ * The above copyright notice and this permission notice shall be included
+ * in all copies or substantial portions of the Software.
+ *
+ * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
+ * OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
+ * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT.  IN NO EVENT SHALL
+ * THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR
+ * OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE,
+ * ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR
+ * OTHER DEALINGS IN THE SOFTWARE.
+ */
+
+#include <stdint.h>
+
+void igt_float_to_half(const float *f, uint16_t *h, unsigned int num);
+void igt_half_to_float(const uint16_t *h, float *f, unsigned int num);
+
diff --git a/lib/igt_x86.c b/lib/igt_x86.c
index 88e514df..6ac700df 100644
--- a/lib/igt_x86.c
+++ b/lib/igt_x86.c
@@ -88,6 +88,10 @@
 #define bit_AVX		(1 << 28)
 #endif
 
+#ifndef bit_F16C
+#define bit_F16C	(1 << 29)
+#endif
+
 #ifndef bit_AVX2
 #define bit_AVX2	(1<<5)
 #endif
@@ -138,6 +142,9 @@ unsigned igt_x86_features(void)
 
 		if (edx & bit_SSE2)
 			features |= SSE2;
+
+		if (ecx & bit_F16C)
+			features |= F16C;
 	}
 
 	if (max >= 7) {
@@ -174,6 +181,8 @@ char *igt_x86_features_to_string(unsigned features, char *line)
 		line += sprintf(line, ", avx");
 	if (features & AVX2)
 		line += sprintf(line, ", avx2");
+	if (features & F16C)
+		line += sprintf(line, ", f16c");
 
 	(void)line;
 
diff --git a/lib/igt_x86.h b/lib/igt_x86.h
index d4f8c343..c7b84dec 100644
--- a/lib/igt_x86.h
+++ b/lib/igt_x86.h
@@ -39,6 +39,7 @@
 #define SSE4_2	0x40
 #define AVX	0x80
 #define AVX2	0x100
+#define F16C	0x200
 
 #if defined(__x86_64__) || defined(__i386__)
 unsigned igt_x86_features(void);
diff --git a/lib/meson.build b/lib/meson.build
index 20c0e3e6..97f701c7 100644
--- a/lib/meson.build
+++ b/lib/meson.build
@@ -12,6 +12,7 @@ lib_sources = [
 	'igt_gpu_power.c',
 	'igt_gt.c',
 	'igt_gvt.c',
+	'igt_halffloat.c',
 	'igt_matrix.c',
 	'igt_perf.c',
 	'igt_primes.c',