diff --git a/src/bit64_conversion.c b/src/bit64_conversion.c deleted file mode 100644 index a6a53fe..0000000 --- a/src/bit64_conversion.c +++ /dev/null @@ -1,70 +0,0 @@ -#include "bit64_conversion.h" - - -void uint32_to_int32(const void* in_buf, size_t n, void* out_buf) { - - R_xlen_t i; - R_xlen_t n_overflow = 0; - - for (i = 0; i < n; i++) { - if (((uint32_t *)in_buf)[i] > INT_MAX) { - ((int32_t *)out_buf)[i] = INT_MIN; - n_overflow++; - } else { - ((int32_t *)out_buf)[i] = ((uint32_t *)in_buf)[i]; - } - } - - if (n_overflow > 0) { - Rf_warning( - "Integer overflow on %zu elements: converting 32bit unsigned integer to 32bit signed integer resulted in NA values", - n_overflow - ); - } - -} - -void int64_to_int32(const void* in_buf, size_t n, void* out_buf, bool is_signed) { - - R_xlen_t i; - R_xlen_t n_overflow = 0; - R_xlen_t n_underflow = 0; - - if (is_signed) { - for (i=0; i INT_MAX) { - ((int32_t *)out_buf)[i] = INT_MIN; - n_overflow++; - } - else if (((int64_t *)in_buf)[i] < INT_MIN) { - ((int32_t *)out_buf)[i] = INT_MIN; - n_underflow++; - } else { - ((int32_t *)out_buf)[i] = ((int64_t *)in_buf)[i]; - } - } - } else { - for (i=0; i INT_MAX) { - ((int32_t *)out_buf)[i] = INT_MIN; - n_overflow++; - } else { - ((int *)out_buf)[i] = ((uint64_t *)in_buf)[i]; - } - } - } - - if (n_overflow > 0) { - Rf_warning( - "Integer overflow on %zu elements: converting 64bit integer to 32bit integer resulted in NA values", - n_overflow - ); - } - if (n_underflow > 0) { - Rf_warning( - "Integer underflow on %zu elements: converting 64bit integer to 32bit integer resulted in NA values", - n_underflow - ); - } - -} diff --git a/src/bit64_conversion.h b/src/bit64_conversion.h index a2b8360..112000b 100644 --- a/src/bit64_conversion.h +++ b/src/bit64_conversion.h @@ -1,4 +1,87 @@ #include "grumpy.h" -void uint32_to_int32(const void* in_buf, size_t n, void* out_buf); -void int64_to_int32(const void* in_buf, size_t n, void* out_buf, bool is_signed); +static inline void uint32_to_int32(const void* restrict in_buf, size_t n, void* restrict out_buf) { + + const uint32_t* restrict in = (const uint32_t *)in_buf; + int32_t* restrict out = (int32_t *)out_buf; + R_xlen_t i; + + // This needs to be done in two passes to allow for auto-vectorization of the conversion loop by + // the compiler. + for (i = 0; i < n; i++) { + out[i] = (in[i] > INT_MAX) ? INT_MIN : (int32_t)in[i]; + } + + R_xlen_t n_overflow = 0; + for (i = 0; i < n; i++) { + n_overflow += (in[i] > INT_MAX); + } + + if (n_overflow > 0) { + Rf_warning( + "Integer overflow on %zu elements: converting 32bit unsigned integer to 32bit signed integer resulted in NA values", + n_overflow + ); + } + +} + +static inline void int64_to_int32(const void* restrict in_buf, size_t n, void* restrict out_buf) { + + const int64_t* restrict in = (const int64_t *)in_buf; + int32_t* restrict out = (int32_t *)out_buf; + R_xlen_t i; + + // This needs to be done in two passes to allow for auto-vectorization of the conversion loop by + // the compiler. + for (i = 0; i < n; i++) { + out[i] = (in[i] > INT_MAX || in[i] < INT_MIN) ? INT_MIN : (int32_t)in[i]; + } + + R_xlen_t n_overflow = 0; + R_xlen_t n_underflow = 0; + for (i = 0; i < n; i++) { + n_overflow += (in[i] > INT_MAX); + n_underflow += (in[i] < INT_MIN); + } + + if (n_overflow > 0) { + Rf_warning( + "Integer overflow on %zu elements: converting 64bit integer to 32bit integer resulted in NA values", + n_overflow + ); + } + if (n_underflow > 0) { + Rf_warning( + "Integer underflow on %zu elements: converting 64bit integer to 32bit integer resulted in NA values", + n_underflow + ); + } + +} + +static inline void uint64_to_int32(const void* restrict in_buf, size_t n, void* restrict out_buf) { + + const uint64_t* restrict in = (const uint64_t *)in_buf; + int32_t* restrict out = (int32_t *)out_buf; + R_xlen_t i; + + // This needs to be done in two passes to allow for auto-vectorization of the conversion loop by + // the compiler. + for (i = 0; i < n; i++) { + out[i] = (in[i] > INT_MAX) ? INT_MIN : (int32_t)in[i]; + } + + R_xlen_t n_overflow = 0; + for (i = 0; i < n; i++) { + n_overflow += (in[i] > INT_MAX); + } + + if (n_overflow > 0) { + Rf_warning( + "Integer overflow on %zu elements: converting 64bit unsigned integer to 32bit signed integer resulted in NA values", + n_overflow + ); + } + +} diff --git a/src/float16_conversion.c b/src/float16_conversion.c deleted file mode 100644 index ce1bc0f..0000000 --- a/src/float16_conversion.c +++ /dev/null @@ -1,41 +0,0 @@ -#include "float16_conversion.h" - -/* this function is based on the float16->float32 implementation found at - * https://gist.github.com/milhidaka/95863906fe828198f47991c813dbe233 - * as well as the process described - * https://fgiesen.wordpress.com/2012/03/28/half-to-float-done-quic/ - */ -double float16_to_float64(uint16_t float16_value) { - // float16=1bit: sign, 5bit: exponent, 10bit: fraction - // float64=1bit: sign, 11bit: exponent, 52bit: fraction - const uint64_t sign = float16_value >> 15; - uint64_t exponent = (float16_value >> 10) & 0x1F; - uint64_t fraction = (float16_value & 0x3FF); - uint64_t float64_value; - double res; - if (exponent == 0) { - if (fraction == 0) { - /* zero */ - float64_value = (sign << 63); - } else { - /* denormalised number */ - exponent = 1023 - 14; - while ((fraction & (1 << 10)) == 0) { - exponent--; - fraction <<= 1; - } - fraction &= 0x3FF; - float64_value = (sign << 63) | (exponent << 52) | (fraction << 42); - } - } else if (exponent == 0x1F) { - /* Inf or NaN */ - float64_value = (sign << 63) | (0x7FFULL << 52) | (fraction << 42); - } else { - /* ordinary number */ - float64_value = (sign << 63) | ((exponent + (1023-15)) << 52) | (fraction << 42); - } - - // we do this to avoid GCC warnings about casting uint64_t to double - memcpy(&res, &float64_value, sizeof(double)); - return res; -} diff --git a/src/float16_conversion.h b/src/float16_conversion.h index 7705052..7ad0b26 100644 --- a/src/float16_conversion.h +++ b/src/float16_conversion.h @@ -1,3 +1,42 @@ #include "grumpy.h" -double float16_to_float64(uint16_t float16_value); +/* this function is based on the float16->float32 implementation found at + * https://gist.github.com/milhidaka/95863906fe828198f47991c813dbe233 + * as well as the process described + * https://fgiesen.wordpress.com/2012/03/28/half-to-float-done-quic/ + */ +static inline double float16_to_float64(uint16_t float16_value) { + // float16=1bit: sign, 5bit: exponent, 10bit: fraction + // float64=1bit: sign, 11bit: exponent, 52bit: fraction + const uint64_t sign = float16_value >> 15; + uint64_t exponent = (float16_value >> 10) & 0x1F; + uint64_t fraction = (float16_value & 0x3FF); + uint64_t float64_value; + double res; + if (exponent == 0) { + if (fraction == 0) { + /* zero */ + float64_value = (sign << 63); + } else { + /* denormalised number */ + exponent = 1023 - 14; + while ((fraction & (1 << 10)) == 0) { + exponent--; + fraction <<= 1; + } + fraction &= 0x3FF; + float64_value = (sign << 63) | (exponent << 52) | (fraction << 42); + } + } else if (exponent == 0x1F) { + /* Inf or NaN */ + float64_value = (sign << 63) | (0x7FFULL << 52) | (fraction << 42); + } else { + /* ordinary number */ + float64_value = (sign << 63) | ((exponent + (1023-15)) << 52) | (fraction << 42); + } + + // we do this to avoid GCC warnings about casting uint64_t to double + memcpy(&res, &float64_value, sizeof(double)); + return res; +} + diff --git a/src/type_conversion.c b/src/type_conversion.c index 051cb67..ab2c300 100644 --- a/src/type_conversion.c +++ b/src/type_conversion.c @@ -50,9 +50,8 @@ SEXP type_convert_int(SEXP input, SEXP _n_bytes) { p_data[i] = ((const int8_t *)raw_buffer)[i]; } } else if(n_bytes == 2) { - const int16_t *mock_buffer = (const int16_t *)raw_buffer; for (i = 0; i < data_length; i++) { - p_data[i] = mock_buffer[i]; + p_data[i] = ((const int16_t *)raw_buffer)[i]; } } else if(n_bytes == 4) { memcpy(p_data, raw_buffer, length); @@ -60,7 +59,7 @@ SEXP type_convert_int(SEXP input, SEXP _n_bytes) { // for now we convert to 32bit int and overflow values are NA_integer int bit64conversion = 0; if (bit64conversion == 0) { - int64_to_int32(raw_buffer, data_length, p_data, true); + int64_to_int32(raw_buffer, data_length, p_data); } } @@ -88,9 +87,8 @@ SEXP type_convert_uint(SEXP input, SEXP _n_bytes) { p_data[i] = ((const uint8_t *)raw_buffer)[i]; } } else if(n_bytes == 2) { - const uint16_t *mock_buffer = (const uint16_t *)raw_buffer; for (i = 0; i < data_length; i++) { - p_data[i] = mock_buffer[i]; + p_data[i] = ((const uint16_t *)raw_buffer)[i]; } } else if(n_bytes == 4) { uint32_to_int32(raw_buffer, data_length, p_data); @@ -98,7 +96,7 @@ SEXP type_convert_uint(SEXP input, SEXP _n_bytes) { // for now we convert to 32bit int and overflow values are NA_integer int bit64conversion = 0; if (bit64conversion == 0) { - int64_to_int32(raw_buffer, data_length, p_data, false); + uint64_to_int32(raw_buffer, data_length, p_data); } } @@ -129,9 +127,8 @@ SEXP type_convert_float(SEXP input, SEXP _n_bytes) { } else if(n_bytes == 4) { - const float *mock_buffer = (const float *)raw_buffer; for (i = 0; i < data_length; i++) { - p_data[i] = (double)mock_buffer[i]; + p_data[i] = (double)((const float *)raw_buffer)[i]; } } else if (n_bytes == 8) { diff --git a/touchstone/script.R b/touchstone/script.R index 58a0a71..07323ea 100644 --- a/touchstone/script.R +++ b/touchstone/script.R @@ -24,6 +24,15 @@ touchstone::benchmark_run( n = 100 ) +# float16 +touchstone::benchmark_run( + { + library(grumpy) + }, + read_float16 = read_npy("inst/extdata/test_float16.npy"), + n = 100 +) + # float32 touchstone::benchmark_run( {