Compiler optimizations

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
+    );
+  }
+
+}

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;
+}
+

src/type_conversion.c

@@ -50,17 +50,16 @@ 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);
   } else if (n_bytes == 8) {
     // 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,17 +87,16 @@ 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);
   } else if (n_bytes == 8) {
     // 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) {

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(
   {