1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
|
/**
* @title Transforming a Matrix in Parallel using RcppParallel
* @author JJ Allaire
* @license GPL (>= 2)
*/
#include <Rcpp.h>
using namespace Rcpp;
#include <cmath>
#include <algorithm>
double squareRoot(double x) {
return ::sqrt(x);
}
// [[Rcpp::export]]
NumericMatrix matrixSqrt(NumericMatrix orig) {
// allocate the matrix we will return
NumericMatrix mat(orig.nrow(), orig.ncol());
// transform it
std::transform(orig.begin(), orig.end(), mat.begin(), squareRoot);
// return the new matrix
return mat;
}
// [[Rcpp::depends(RcppParallel)]]
#include <RcppParallel.h>
using namespace RcppParallel;
struct SquareRoot : public Worker
{
// source matrix
const RMatrix<double> input;
// destination matrix
RMatrix<double> output;
// initialize with source and destination
SquareRoot(const NumericMatrix input, NumericMatrix output)
: input(input), output(output) {}
// take the square root of the range of elements requested
void operator()(std::size_t begin, std::size_t end) {
std::transform(input.begin() + begin,
input.begin() + end,
output.begin() + begin,
squareRoot);
}
};
// [[Rcpp::export]]
NumericMatrix parallelMatrixSqrt(NumericMatrix x) {
// allocate the output matrix
NumericMatrix output(x.nrow(), x.ncol());
// SquareRoot functor (pass input and output matrixes)
SquareRoot squareRoot(x, output);
// call parallelFor to do the work
parallelFor(0, x.nrow() * x.ncol(), squareRoot);
// return the output matrix
return output;
}
|
