add a 'map-reduce' example; fix

Author: kevinushey

examples/boost-simd-map-reduce.cpp

@@ -0,0 +1,51 @@
+// [[Rcpp::depends(RcppParallel)]]
+#define RCPP_PARALLEL_USE_SIMD
+#include <RcppParallel.h>
+using namespace RcppParallel;
+
+#include <Rcpp.h>
+using namespace Rcpp;
+
+class SumOfSquaresReducer
+{
+public:
+   
+   explicit SumOfSquaresReducer(double mean)
+      : mean_(mean)
+   {}
+   
+   template <typename T>
+   void map(const T& self, T* pBuffer)
+   {
+      *pBuffer += boost::simd::sqr(self - mean_);
+   }
+   
+   template <typename T, typename U>
+   void reduce(const T& data, U* pBuffer)
+   {
+      *pBuffer += boost::simd::sum(data);
+   }
+   
+private:
+   double mean_;
+};
+
+// [[Rcpp::export]]
+double simdVar(NumericVector x)
+{
+   double total = simdReduce(x.begin(), x.end(), 0.0, simd_ops::plus());
+   double n = x.size();
+   double mean = total / n;
+   
+   double ssq = simdMapReduce(x.begin(), x.end(), 0.0, SumOfSquaresReducer(mean));
+   
+   return ssq / (n - 1);
+}
+
+/*** R
+x <- rnorm(1024 * 10000)
+var(x)
+simdVar(x)
+library(microbenchmark)
+microbenchmark(var(x), simdVar(x))
+*/

inst/include/RcppParallel/simd/algorithm.h

@@ -60,8 +60,8 @@ U simdMapReduce(const T* it, const T* end, U init, MapReducer mapper)
    for (; it != aligned_end; it += N)
       mapper.map(boost::simd::aligned_load<vT>(it), &buffer);
 
-   // Reduce our buffer, and join it with the initial value
-   mapper.map(mapper.reduce(buffer), &init);
+   // Reduce the buffer, joining it into the scalar vale
+   mapper.reduce(buffer, &init);
 
    // Leftover unaligned region
    for (; it != end; ++it)